1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Authors: Svein Lindal, Daniel Lohner *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
16 ////////////////////////////////////////////////
17 //---------------------------------------------
18 // Class handling all kinds of selection cuts for
19 // Gamma Conversion analysis
20 //---------------------------------------------
21 ////////////////////////////////////////////////
23 #include "AliConversionCuts.h"
25 #include "AliKFVertex.h"
26 #include "AliAODTrack.h"
27 #include "AliESDtrack.h"
28 #include "AliAnalysisManager.h"
29 #include "AliInputEventHandler.h"
30 #include "AliMCEventHandler.h"
31 #include "AliAODHandler.h"
32 #include "AliPIDResponse.h"
37 #include "AliAODConversionMother.h"
38 #include "TObjString.h"
39 #include "AliAODEvent.h"
40 #include "AliESDEvent.h"
41 #include "AliCentrality.h"
45 #include "AliGenCocktailEventHeader.h"
46 #include "AliGenDPMjetEventHeader.h"
47 #include "AliGenPythiaEventHeader.h"
48 #include "AliGenHijingEventHeader.h"
49 #include "AliTriggerAnalysis.h"
50 #include "AliV0ReaderV1.h"
51 #include "AliAODMCParticle.h"
52 #include "AliAODMCHeader.h"
58 ClassImp(AliConversionCuts)
61 const char* AliConversionCuts::fgkCutNames[AliConversionCuts::kNCuts] = {
66 "MultiplicityMethod",//4
68 "RejectExtraSignals",//6
76 "piMomdedxSigmaCut",//14
77 "piMaxMomdedxSigmaCut",//15
78 "LowPRejectionSigmaCut",//16
83 "DoPhotonAsymmetryCut",//21
84 "CosinePointingAngle", //22
85 "SharedElectronCuts", //23
86 "RejectToCloseV0s", //24
92 //________________________________________________________________________
93 AliConversionCuts::AliConversionCuts(const char *name,const char *title) :
94 AliAnalysisCuts(name,title),
110 fLineCutZRSlopeMin(0.),
111 fLineCutZValueMin(0),
112 fChi2CutConversion(1000),
113 fPIDProbabilityCutNegativeParticle(0),
114 fPIDProbabilityCutPositiveParticle(0),
115 fDodEdxSigmaCut(kTRUE),
116 fDoTOFsigmaCut(kFALSE),
117 fPIDTRDEfficiency(1),
119 fPIDnSigmaAboveElectronLine(100),
120 fPIDnSigmaBelowElectronLine(-100),
121 fTofPIDnSigmaAboveElectronLine(100),
122 fTofPIDnSigmaBelowElectronLine(-100),
123 fPIDnSigmaAbovePionLine(0),
124 fPIDnSigmaAbovePionLineHighPt(-100),
125 fPIDMinPnSigmaAbovePionLine(0),
126 fPIDMaxPnSigmaAbovePionLine(0),
127 fDoKaonRejectionLowP(kFALSE),
128 fDoProtonRejectionLowP(kFALSE),
129 fDoPionRejectionLowP(kFALSE),
130 fPIDnSigmaAtLowPAroundKaonLine(0),
131 fPIDnSigmaAtLowPAroundProtonLine(0),
132 fPIDnSigmaAtLowPAroundPionLine(0),
133 fPIDMinPKaonRejectionLowP(1.5),
134 fPIDMinPProtonRejectionLowP(2),
135 fPIDMinPPionRejectionLowP(0),
136 fDoQtGammaSelection(kTRUE),
137 fDoHighPtQtGammaSelection(kFALSE),
145 fUseEtaMinCut(kFALSE),
146 fUseOnFlyV0Finder(kTRUE),
147 fDoPhotonAsymmetryCut(kTRUE),
148 fMinPPhotonAsymmetryCut(100.),
149 fMinPhotonAsymmetry(0.),
151 fDetectorCentrality(0),
152 fModCentralityClass(0),
156 fUseCorrectedTPCClsInfo(kFALSE),
158 fMultiplicityMethod(0),
160 fRemovePileUp(kFALSE),
161 fOpeningAngle(0.005),
163 fCosPAngleCut(10000),
164 fDoToCloseV0sCut(kFALSE),
165 fRejectExtraSignals(0),
167 fDoSharedElecCut(kFALSE),
168 fOfflineTriggerMask(0),
172 fElectronArraySize(500),
173 fElectronLabelArray(NULL),
174 fDCAZPrimVtxCut(1000),
175 fDCARPrimVtxCut(1000),
176 fConversionPointXArray(0.0),
177 fConversionPointYArray(0.0),
178 fConversionPointZArray(0.0),
180 fNotRejectedStart(NULL),
181 fNotRejectedEnd(NULL),
182 fGeneratorNames(NULL),
187 fDoReweightHistoMCPi0(kFALSE),
188 fDoReweightHistoMCEta(kFALSE),
189 fDoReweightHistoMCK0s(kFALSE),
190 fPathTrFReweighting(""),
191 fNameHistoReweightingPi0(""),
192 fNameHistoReweightingEta(""),
193 fNameHistoReweightingK0s(""),
198 hTPCdEdxbefore(NULL),
200 hTPCdEdxSigbefore(NULL),
201 hTPCdEdxSigafter(NULL),
207 hInvMassbefore(NULL),
208 hArmenterosbefore(NULL),
210 hArmenterosafter(NULL),
211 hAcceptanceCuts(NULL),
215 hCentralityVsNumberOfPrimaryTracks(NULL),
218 hTriggerClassSelected(NULL),
219 hReweightMCHistPi0(NULL),
220 hReweightMCHistEta(NULL),
221 hReweightMCHistK0s(NULL),
226 fTriggerSelectedManually(kFALSE),
227 fSpecialTriggerName("")
231 for(Int_t jj=0;jj<kNCuts;jj++){fCuts[jj]=0;}
232 fCutString=new TObjString((GetCutNumber()).Data());
234 fElectronLabelArray = new Int_t[fElectronArraySize];
235 fUtils = new AliAnalysisUtils();
236 //if you do not want to apply the cut on the distance between the SPD and TRK vertex:
237 //fUtils->SetCutOnZVertexSPD(kFALSE);
242 //________________________________________________________________________
243 AliConversionCuts::AliConversionCuts(const AliConversionCuts &ref) :
244 AliAnalysisCuts(ref),
246 fHeaderList(ref.fHeaderList),
248 fEventQuality(ref.fEventQuality),
251 fEtaCut(ref.fEtaCut),
252 fEtaCutMin(ref.fEtaCutMin),
254 fSinglePtCut(ref.fSinglePtCut),
256 fMinClsTPC(ref.fMinClsTPC),
257 fMinClsTPCToF(ref.fMinClsTPCToF),
258 fLineCutZRSlope(ref.fLineCutZRSlope),
259 fLineCutZValue(ref.fLineCutZValue),
260 fLineCutZRSlopeMin(ref.fLineCutZRSlopeMin),
261 fLineCutZValueMin(ref.fLineCutZValueMin),
262 fChi2CutConversion(ref.fChi2CutConversion),
263 fPIDProbabilityCutNegativeParticle(ref.fPIDProbabilityCutNegativeParticle),
264 fPIDProbabilityCutPositiveParticle(ref.fPIDProbabilityCutPositiveParticle),
265 fDodEdxSigmaCut(ref. fDodEdxSigmaCut),
266 fDoTOFsigmaCut(ref.fDoTOFsigmaCut),
267 fPIDTRDEfficiency(ref.fPIDTRDEfficiency),
268 fDoTRDPID(ref.fDoTRDPID),
269 fPIDnSigmaAboveElectronLine(ref.fPIDnSigmaAboveElectronLine),
270 fPIDnSigmaBelowElectronLine(ref.fPIDnSigmaBelowElectronLine),
271 fTofPIDnSigmaAboveElectronLine(ref.fTofPIDnSigmaAboveElectronLine),
272 fTofPIDnSigmaBelowElectronLine(ref.fTofPIDnSigmaBelowElectronLine),
273 fPIDnSigmaAbovePionLine(ref.fPIDnSigmaAbovePionLine),
274 fPIDnSigmaAbovePionLineHighPt(ref.fPIDnSigmaAbovePionLineHighPt),
275 fPIDMinPnSigmaAbovePionLine(ref.fPIDMinPnSigmaAbovePionLine),
276 fPIDMaxPnSigmaAbovePionLine(ref.fPIDMaxPnSigmaAbovePionLine),
277 fDoKaonRejectionLowP(ref.fDoKaonRejectionLowP),
278 fDoProtonRejectionLowP(ref.fDoProtonRejectionLowP),
279 fDoPionRejectionLowP(ref.fDoPionRejectionLowP),
280 fPIDnSigmaAtLowPAroundKaonLine(ref.fPIDnSigmaAtLowPAroundKaonLine),
281 fPIDnSigmaAtLowPAroundProtonLine(ref.fPIDnSigmaAtLowPAroundProtonLine),
282 fPIDnSigmaAtLowPAroundPionLine(ref.fPIDnSigmaAtLowPAroundPionLine),
283 fPIDMinPKaonRejectionLowP(ref.fPIDMinPKaonRejectionLowP),
284 fPIDMinPProtonRejectionLowP(ref.fPIDMinPProtonRejectionLowP),
285 fPIDMinPPionRejectionLowP(ref.fPIDMinPPionRejectionLowP),
286 fDoQtGammaSelection(ref.fDoQtGammaSelection),
287 fDoHighPtQtGammaSelection(ref.fDoHighPtQtGammaSelection),
289 fHighPtQtMax(ref.fHighPtQtMax),
290 fPtBorderForQt(ref.fPtBorderForQt),
291 fXVertexCut(ref.fXVertexCut),
292 fYVertexCut(ref.fYVertexCut),
293 fZVertexCut(ref.fZVertexCut),
294 fNSigmaMass(ref.fNSigmaMass),
295 fUseEtaMinCut(ref.fUseEtaMinCut),
296 fUseOnFlyV0Finder(ref.fUseOnFlyV0Finder),
297 fDoPhotonAsymmetryCut(ref.fDoPhotonAsymmetryCut),
298 fMinPPhotonAsymmetryCut(ref.fMinPPhotonAsymmetryCut),
299 fMinPhotonAsymmetry(ref.fMinPhotonAsymmetry),
300 fIsHeavyIon(ref.fIsHeavyIon),
301 fDetectorCentrality(ref.fDetectorCentrality),
302 fModCentralityClass(ref.fModCentralityClass),
303 fMaxVertexZ(ref.fMaxVertexZ),
304 fCentralityMin(ref.fCentralityMin),
305 fCentralityMax(ref.fCentralityMax),
306 fUseCorrectedTPCClsInfo(ref.fUseCorrectedTPCClsInfo),
307 fUseTOFpid(ref.fUseTOFpid),
308 fMultiplicityMethod(ref.fMultiplicityMethod),
309 fSpecialTrigger(ref.fSpecialTrigger),
310 fRemovePileUp(ref.fRemovePileUp),
311 fOpeningAngle(ref.fOpeningAngle),
312 fPsiPairCut(ref.fPsiPairCut),
313 fCosPAngleCut(ref.fCosPAngleCut),
314 fDoToCloseV0sCut(ref.fDoToCloseV0sCut),
315 fRejectExtraSignals(ref.fRejectExtraSignals),
316 fminV0Dist(ref.fminV0Dist),
317 fDoSharedElecCut(ref.fDoSharedElecCut),
318 fOfflineTriggerMask(ref.fOfflineTriggerMask),
319 fHasV0AND(ref.fHasV0AND),
320 fIsSDDFired(ref.fIsSDDFired),
321 fRandom(ref.fRandom),
322 fElectronArraySize(ref.fElectronArraySize),
323 fElectronLabelArray(NULL),
324 fDCAZPrimVtxCut(ref.fDCAZPrimVtxCut),
325 fDCARPrimVtxCut(ref.fDCAZPrimVtxCut),
326 fConversionPointXArray(ref.fConversionPointXArray),
327 fConversionPointYArray(ref.fConversionPointYArray),
328 fConversionPointZArray(ref.fConversionPointZArray),
329 fnHeaders(ref.fnHeaders),
330 fNotRejectedStart(NULL),
331 fNotRejectedEnd(NULL),
332 fGeneratorNames(ref.fGeneratorNames),
335 fEtaShift(ref.fEtaShift),
336 fDoEtaShift(ref.fDoEtaShift),
337 fDoReweightHistoMCPi0(ref.fDoReweightHistoMCPi0),
338 fDoReweightHistoMCEta(ref.fDoReweightHistoMCEta),
339 fDoReweightHistoMCK0s(ref.fDoReweightHistoMCK0s),
340 fPathTrFReweighting(ref.fPathTrFReweighting),
341 fNameHistoReweightingPi0(ref.fNameHistoReweightingPi0),
342 fNameHistoReweightingEta(ref.fNameHistoReweightingEta),
343 fNameHistoReweightingK0s(ref.fNameHistoReweightingK0s),
344 fNameFitDataPi0(ref.fNameFitDataPi0),
345 fNameFitDataEta(ref.fNameFitDataEta),
346 fNameFitDataK0s(ref.fNameFitDataK0s),
348 hTPCdEdxbefore(NULL),
350 hTPCdEdxSigbefore(NULL),
351 hTPCdEdxSigafter(NULL),
357 hInvMassbefore(NULL),
358 hArmenterosbefore(NULL),
360 hArmenterosafter(NULL),
361 hAcceptanceCuts(NULL),
365 hCentralityVsNumberOfPrimaryTracks(NULL),
368 hTriggerClassSelected(NULL),
369 hReweightMCHistPi0(ref.hReweightMCHistPi0),
370 hReweightMCHistEta(ref.hReweightMCHistEta),
371 hReweightMCHistK0s(ref.hReweightMCHistK0s),
372 fFitDataPi0(ref.fFitDataPi0),
373 fFitDataEta(ref.fFitDataEta),
374 fFitDataK0s(ref.fFitDataK0s),
375 fPreSelCut(ref.fPreSelCut),
376 fTriggerSelectedManually(ref.fTriggerSelectedManually),
377 fSpecialTriggerName(ref.fSpecialTriggerName)
380 for(Int_t jj=0;jj<kNCuts;jj++){fCuts[jj]=ref.fCuts[jj];}
381 fCutString=new TObjString((GetCutNumber()).Data());
382 fElectronLabelArray = new Int_t[fElectronArraySize];
383 fUtils = new AliAnalysisUtils();
384 // dont copy histograms (if you like histograms, call InitCutHistograms())
389 //________________________________________________________________________
390 AliConversionCuts::~AliConversionCuts() {
392 //Deleting fHistograms leads to seg fault it it's added to output collection of a task
394 // delete fHistograms;
395 // fHistograms = NULL;
396 if(fCutString != NULL){
400 if(fElectronLabelArray){
401 delete fElectronLabelArray;
402 fElectronLabelArray = NULL;
404 if(fNotRejectedStart){
405 delete[] fNotRejectedStart;
406 fNotRejectedStart = NULL;
409 delete[] fNotRejectedEnd;
410 fNotRejectedEnd = NULL;
413 delete[] fGeneratorNames;
414 fGeneratorNames = NULL;
423 //________________________________________________________________________
424 void AliConversionCuts::InitCutHistograms(TString name, Bool_t preCut){
426 // Initialize Cut Histograms for QA (only initialized and filled if function is called)
427 TH1::AddDirectory(kFALSE);
429 if(fHistograms != NULL){
433 if(fHistograms==NULL){
434 fHistograms=new TList();
435 fHistograms->SetOwner(kTRUE);
436 if(name=="")fHistograms->SetName(Form("ConvCuts_%s",GetCutNumber().Data()));
437 else fHistograms->SetName(Form("%s_%s",name.Data(),GetCutNumber().Data()));
440 if (hReweightMCHistPi0){
441 hReweightMCHistPi0->SetName("MCInputForWeightingPi0");
442 fHistograms->Add(hReweightMCHistPi0);
444 if (hReweightMCHistEta){
445 hReweightMCHistEta->SetName("MCInputForWeightingEta");
446 fHistograms->Add(hReweightMCHistEta);
448 if (hReweightMCHistK0s){
449 hReweightMCHistK0s->SetName("MCInputForWeightingK0s");
450 fHistograms->Add(hReweightMCHistK0s);
453 fFitDataPi0->SetName("DataFitForWeightingPi0");
454 fHistograms->Add(fFitDataPi0);
457 fFitDataEta->SetName("DataFitForWeightingEta");
458 fHistograms->Add(fFitDataEta);
461 fFitDataK0s->SetName("DataFitForWeightingK0s");
462 fHistograms->Add(fFitDataK0s);
465 hCutIndex=new TH1F(Form("IsPhotonSelected %s",GetCutNumber().Data()),"IsPhotonSelected",10,-0.5,9.5);
466 hCutIndex->GetXaxis()->SetBinLabel(kPhotonIn+1,"in");
467 hCutIndex->GetXaxis()->SetBinLabel(kOnFly+1,"onfly");
468 hCutIndex->GetXaxis()->SetBinLabel(kNoTracks+1,"no tracks");
469 hCutIndex->GetXaxis()->SetBinLabel(kdEdxCuts+1,"dEdx");
470 hCutIndex->GetXaxis()->SetBinLabel(kTrackCuts+1,"Track cuts");
471 hCutIndex->GetXaxis()->SetBinLabel(kConvPointFail+1,"ConvPoint fail");
472 hCutIndex->GetXaxis()->SetBinLabel(kPhotonCuts+1,"PhotonCuts");
473 hCutIndex->GetXaxis()->SetBinLabel(kPhotonOut+1,"out");
474 fHistograms->Add(hCutIndex);
477 hTrackCuts=new TH1F(Form("TrackCuts %s",GetCutNumber().Data()),"TrackCuts",9,-0.5,8.5);
478 hTrackCuts->GetXaxis()->SetBinLabel(1,"in");
479 hTrackCuts->GetXaxis()->SetBinLabel(2,"likesign");
480 hTrackCuts->GetXaxis()->SetBinLabel(3,"ntpccl");
481 hTrackCuts->GetXaxis()->SetBinLabel(4,"acceptance");
482 hTrackCuts->GetXaxis()->SetBinLabel(5,"singlept");
483 hTrackCuts->GetXaxis()->SetBinLabel(6,"TPCrefit");
484 hTrackCuts->GetXaxis()->SetBinLabel(7,"kink");
485 hTrackCuts->GetXaxis()->SetBinLabel(8,"out");
486 fHistograms->Add(hTrackCuts);
489 hPhotonCuts=new TH1F(Form("PhotonCuts %s",GetCutNumber().Data()),"PhotonCuts",14,-0.5,13.5);
490 hPhotonCuts->GetXaxis()->SetBinLabel(1,"in");
491 hPhotonCuts->GetXaxis()->SetBinLabel(2,"qtcut");
492 hPhotonCuts->GetXaxis()->SetBinLabel(3,"chi2");
493 hPhotonCuts->GetXaxis()->SetBinLabel(4,"acceptance");
494 hPhotonCuts->GetXaxis()->SetBinLabel(5,"asymmetry");
495 hPhotonCuts->GetXaxis()->SetBinLabel(6,"pidprob");
496 hPhotonCuts->GetXaxis()->SetBinLabel(7,"cortpcclinfo");
497 hPhotonCuts->GetXaxis()->SetBinLabel(8,"PsiPair");
498 hPhotonCuts->GetXaxis()->SetBinLabel(9,"CosPAngle");
499 hPhotonCuts->GetXaxis()->SetBinLabel(10,"DCA R");
500 hPhotonCuts->GetXaxis()->SetBinLabel(11,"DCA Z");
501 hPhotonCuts->GetXaxis()->SetBinLabel(12,"out");
502 fHistograms->Add(hPhotonCuts);
505 hInvMassbefore=new TH1F(Form("InvMass_before %s",GetCutNumber().Data()),"InvMass_before",1000,0,0.3);
506 fHistograms->Add(hInvMassbefore);
507 hArmenterosbefore=new TH2F(Form("Armenteros_before %s",GetCutNumber().Data()),"Armenteros_before",200,-1,1,1000,0,1.);
508 fHistograms->Add(hArmenterosbefore);
510 hInvMassafter=new TH1F(Form("InvMass_after %s",GetCutNumber().Data()),"InvMass_after",1000,0,0.3);
511 fHistograms->Add(hInvMassafter);
512 hArmenterosafter=new TH2F(Form("Armenteros_after %s",GetCutNumber().Data()),"Armenteros_after",200,-1,1,250,0,0.25);
513 fHistograms->Add(hArmenterosafter);
515 hAcceptanceCuts=new TH1F(Form("PhotonAcceptanceCuts %s",GetCutNumber().Data()),"PhotonAcceptanceCuts",10,-0.5,9.5);
516 hAcceptanceCuts->GetXaxis()->SetBinLabel(1,"in");
517 hAcceptanceCuts->GetXaxis()->SetBinLabel(2,"maxR");
518 hAcceptanceCuts->GetXaxis()->SetBinLabel(3,"minR");
519 hAcceptanceCuts->GetXaxis()->SetBinLabel(4,"line");
520 hAcceptanceCuts->GetXaxis()->SetBinLabel(5,"maxZ");
521 hAcceptanceCuts->GetXaxis()->SetBinLabel(6,"eta");
522 hAcceptanceCuts->GetXaxis()->SetBinLabel(7,"minpt");
523 hAcceptanceCuts->GetXaxis()->SetBinLabel(8,"out");
524 fHistograms->Add(hAcceptanceCuts);
527 hdEdxCuts=new TH1F(Form("dEdxCuts %s",GetCutNumber().Data()),"dEdxCuts",10,-0.5,9.5);
528 hdEdxCuts->GetXaxis()->SetBinLabel(1,"in");
529 hdEdxCuts->GetXaxis()->SetBinLabel(2,"TPCelectron");
530 hdEdxCuts->GetXaxis()->SetBinLabel(3,"TPCpion");
531 hdEdxCuts->GetXaxis()->SetBinLabel(4,"TPCpionhighp");
532 hdEdxCuts->GetXaxis()->SetBinLabel(5,"TPCkaonlowprej");
533 hdEdxCuts->GetXaxis()->SetBinLabel(6,"TPCprotonlowprej");
534 hdEdxCuts->GetXaxis()->SetBinLabel(7,"TPCpionlowprej");
535 hdEdxCuts->GetXaxis()->SetBinLabel(8,"TOFelectron");
536 hdEdxCuts->GetXaxis()->SetBinLabel(9,"TRDelectron");
537 hdEdxCuts->GetXaxis()->SetBinLabel(10,"out");
538 fHistograms->Add(hdEdxCuts);
540 TAxis *AxisBeforedEdx = NULL;
541 TAxis *AxisBeforedEdxSig = NULL;
542 TAxis *AxisBeforeTOF = NULL;
543 TAxis *AxisBeforeTOFSig = NULL;
545 hTPCdEdxbefore=new TH2F(Form("Gamma_dEdx_before %s",GetCutNumber().Data()),"dEdx Gamma before" ,150,0.03,20,800,0,200);
546 fHistograms->Add(hTPCdEdxbefore);
547 AxisBeforedEdx = hTPCdEdxbefore->GetXaxis();
548 hTPCdEdxSigbefore=new TH2F(Form("Gamma_dEdxSig_before %s",GetCutNumber().Data()),"dEdx Sigma Gamma before" ,150,0.03,20,400,-10,10);
549 fHistograms->Add(hTPCdEdxSigbefore);
550 AxisBeforedEdxSig = hTPCdEdxSigbefore->GetXaxis();
552 hTOFbefore=new TH2F(Form("Gamma_TOF_before %s",GetCutNumber().Data()),"TOF Gamma before" ,150,0.03,20,11000,-1000,10000);
553 fHistograms->Add(hTOFbefore);
554 AxisBeforeTOF = hTOFbefore->GetXaxis();
555 hTOFSigbefore=new TH2F(Form("Gamma_TOFSig_before %s",GetCutNumber().Data()),"TOF Sigma Gamma before" ,150,0.03,20,400,-6,10);
556 fHistograms->Add(hTOFSigbefore);
557 AxisBeforeTOFSig = hTOFSigbefore->GetXaxis();
560 hTPCdEdxSigafter=new TH2F(Form("Gamma_dEdxSig_after %s",GetCutNumber().Data()),"dEdx Sigma Gamma after" ,150,0.03,20,400, -10,10);
561 fHistograms->Add(hTPCdEdxSigafter);
563 hTPCdEdxafter=new TH2F(Form("Gamma_dEdx_after %s",GetCutNumber().Data()),"dEdx Gamma after" ,150,0.03,20,800,0,200);
564 fHistograms->Add(hTPCdEdxafter);
566 hTOFSigafter=new TH2F(Form("Gamma_TOFSig_after %s",GetCutNumber().Data()),"TOF Sigma Gamma after" ,150,0.03,20,400,-6,10);
567 fHistograms->Add(hTOFSigafter);
569 TAxis *AxisAfter = hTPCdEdxSigafter->GetXaxis();
570 Int_t bins = AxisAfter->GetNbins();
571 Double_t from = AxisAfter->GetXmin();
572 Double_t to = AxisAfter->GetXmax();
573 Double_t *newBins = new Double_t[bins+1];
575 Double_t factor = TMath::Power(to/from, 1./bins);
576 for(Int_t i=1; i<=bins; ++i) newBins[i] = factor * newBins[i-1];
577 AxisAfter->Set(bins, newBins);
578 AxisAfter = hTOFSigafter->GetXaxis();
579 AxisAfter->Set(bins, newBins);
580 AxisAfter = hTPCdEdxafter->GetXaxis();
581 AxisAfter->Set(bins, newBins);
583 AxisBeforedEdx->Set(bins, newBins);
584 AxisBeforeTOF->Set(bins, newBins);
585 AxisBeforedEdxSig->Set(bins, newBins);
586 AxisBeforeTOFSig->Set(bins, newBins);
590 // Event Cuts and Info
592 hV0EventCuts=new TH1F(Form("ESD_EventCuts %s",GetCutNumber().Data()),"Event Cuts",7,-0.5,6.5);
593 hV0EventCuts->GetXaxis()->SetBinLabel(1,"in");
594 hV0EventCuts->GetXaxis()->SetBinLabel(2,"OfflineTrigger");
595 hV0EventCuts->GetXaxis()->SetBinLabel(3,"nvtxcontr");
596 hV0EventCuts->GetXaxis()->SetBinLabel(4,"VertexZ");
597 hV0EventCuts->GetXaxis()->SetBinLabel(5,"pileup");
598 hV0EventCuts->GetXaxis()->SetBinLabel(6,"centrsel");
599 hV0EventCuts->GetXaxis()->SetBinLabel(7,"out");
600 fHistograms->Add(hV0EventCuts);
602 hCentrality=new TH1F(Form("Centrality %s",GetCutNumber().Data()),"Centrality",1000,0,100);
603 fHistograms->Add(hCentrality);
604 hCentralityVsNumberOfPrimaryTracks=new TH2F(Form("Centrality vs Primary Tracks %s",GetCutNumber().Data()),"Centrality vs Primary Tracks ",100,0,100,4000,0,4000);
605 fHistograms->Add(hCentralityVsNumberOfPrimaryTracks);
606 hVertexZ=new TH1F(Form("VertexZ %s",GetCutNumber().Data()),"VertexZ",1000,-50,50);
607 fHistograms->Add(hVertexZ);
609 hTriggerClass= new TH1F(Form("OfflineTrigger %s",GetCutNumber().Data()),"OfflineTrigger",35,-0.5,34.5);
610 hTriggerClass->GetXaxis()->SetBinLabel( 1,"kMB");
611 hTriggerClass->GetXaxis()->SetBinLabel( 2,"kINT7");
612 hTriggerClass->GetXaxis()->SetBinLabel( 3,"kMUON");
613 hTriggerClass->GetXaxis()->SetBinLabel( 4,"kHighMult");
614 hTriggerClass->GetXaxis()->SetBinLabel( 5,"kKEMC1");
615 hTriggerClass->GetXaxis()->SetBinLabel( 6,"kCINT5");
616 hTriggerClass->GetXaxis()->SetBinLabel( 7,"kCMUS5/kMUSPB");
617 hTriggerClass->GetXaxis()->SetBinLabel( 8,"kMUSH7/kMUSHPB");
618 hTriggerClass->GetXaxis()->SetBinLabel( 9,"kMUL7/kMuonLikePB");
619 hTriggerClass->GetXaxis()->SetBinLabel(10,"kMUU7/kMuonUnlikePB");
620 hTriggerClass->GetXaxis()->SetBinLabel(11,"kEMC7/kEMC8");
621 hTriggerClass->GetXaxis()->SetBinLabel(12,"kMUS7");
622 hTriggerClass->GetXaxis()->SetBinLabel(13,"kPHI1");
623 hTriggerClass->GetXaxis()->SetBinLabel(14,"kPHI7/kPHI8/kPHOSPb");
624 hTriggerClass->GetXaxis()->SetBinLabel(15,"kEMCEJE");
625 hTriggerClass->GetXaxis()->SetBinLabel(16,"kEMCEGA");
626 hTriggerClass->GetXaxis()->SetBinLabel(17,"kCentral");
627 hTriggerClass->GetXaxis()->SetBinLabel(18,"kSemiCentral");
628 hTriggerClass->GetXaxis()->SetBinLabel(19,"kDG5");
629 hTriggerClass->GetXaxis()->SetBinLabel(20,"kZED");
630 hTriggerClass->GetXaxis()->SetBinLabel(21,"kSPI7/kSPI");
631 hTriggerClass->GetXaxis()->SetBinLabel(22,"kINT8");
632 hTriggerClass->GetXaxis()->SetBinLabel(23,"kMuonSingleLowPt8");
633 hTriggerClass->GetXaxis()->SetBinLabel(24,"kMuonSingleHighPt8");
634 hTriggerClass->GetXaxis()->SetBinLabel(25,"kMuonLikeLowPt8");
635 hTriggerClass->GetXaxis()->SetBinLabel(26,"kMuonUnlikeLowPt8");
636 hTriggerClass->GetXaxis()->SetBinLabel(27,"kMuonUnlikeLowPt0");
637 hTriggerClass->GetXaxis()->SetBinLabel(28,"kUserDefined");
638 hTriggerClass->GetXaxis()->SetBinLabel(29,"kTRD");
639 hTriggerClass->GetXaxis()->SetBinLabel(30,"kFastOnly");
640 hTriggerClass->GetXaxis()->SetBinLabel(31,"kAnyINT");
641 hTriggerClass->GetXaxis()->SetBinLabel(32,"kAny");
642 hTriggerClass->GetXaxis()->SetBinLabel(33,"V0AND");
643 hTriggerClass->GetXaxis()->SetBinLabel(34,"NOT kFastOnly");
644 hTriggerClass->GetXaxis()->SetBinLabel(35,"failed Physics Selection");
645 fHistograms->Add(hTriggerClass);
648 hTriggerClassSelected= new TH1F(Form("OfflineTriggerSelected %s",GetCutNumber().Data()),"OfflineTriggerSelected",34,-0.5,33.5);
649 hTriggerClassSelected->GetXaxis()->SetBinLabel( 1,"kMB");
650 hTriggerClassSelected->GetXaxis()->SetBinLabel( 2,"kINT7");
651 hTriggerClassSelected->GetXaxis()->SetBinLabel( 3,"kMUON");
652 hTriggerClassSelected->GetXaxis()->SetBinLabel( 4,"kHighMult");
653 hTriggerClassSelected->GetXaxis()->SetBinLabel( 5,"kKEMC1");
654 hTriggerClassSelected->GetXaxis()->SetBinLabel( 6,"kCINT5");
655 hTriggerClassSelected->GetXaxis()->SetBinLabel( 7,"kCMUS5/kMUSPB");
656 hTriggerClassSelected->GetXaxis()->SetBinLabel( 8,"kMUSH7/kMUSHPB");
657 hTriggerClassSelected->GetXaxis()->SetBinLabel( 9,"kMUL7/kMuonLikePB");
658 hTriggerClassSelected->GetXaxis()->SetBinLabel(10,"kMUU7/kMuonUnlikePB");
659 hTriggerClassSelected->GetXaxis()->SetBinLabel(11,"kEMC7/kEMC8");
660 hTriggerClassSelected->GetXaxis()->SetBinLabel(12,"kMUS7");
661 hTriggerClassSelected->GetXaxis()->SetBinLabel(13,"kPHI1");
662 hTriggerClassSelected->GetXaxis()->SetBinLabel(14,"kPHI7/kPHI8/kPHOSPb");
663 hTriggerClassSelected->GetXaxis()->SetBinLabel(15,"kEMCEJE");
664 hTriggerClassSelected->GetXaxis()->SetBinLabel(16,"kEMCEGA");
665 hTriggerClassSelected->GetXaxis()->SetBinLabel(17,"kCentral");
666 hTriggerClassSelected->GetXaxis()->SetBinLabel(18,"kSemiCentral");
667 hTriggerClassSelected->GetXaxis()->SetBinLabel(19,"kDG5");
668 hTriggerClassSelected->GetXaxis()->SetBinLabel(20,"kZED");
669 hTriggerClassSelected->GetXaxis()->SetBinLabel(21,"kSPI7/kSPI");
670 hTriggerClassSelected->GetXaxis()->SetBinLabel(22,"kINT8");
671 hTriggerClassSelected->GetXaxis()->SetBinLabel(23,"kMuonSingleLowPt8");
672 hTriggerClassSelected->GetXaxis()->SetBinLabel(24,"kMuonSingleHighPt8");
673 hTriggerClassSelected->GetXaxis()->SetBinLabel(25,"kMuonLikeLowPt8");
674 hTriggerClassSelected->GetXaxis()->SetBinLabel(26,"kMuonUnlikeLowPt8");
675 hTriggerClassSelected->GetXaxis()->SetBinLabel(27,"kMuonUnlikeLowPt0");
676 hTriggerClassSelected->GetXaxis()->SetBinLabel(28,"kUserDefined");
677 hTriggerClassSelected->GetXaxis()->SetBinLabel(29,"kTRD");
678 hTriggerClassSelected->GetXaxis()->SetBinLabel(30,"kFastOnly");
679 hTriggerClassSelected->GetXaxis()->SetBinLabel(31,"kAnyINT");
680 hTriggerClassSelected->GetXaxis()->SetBinLabel(32,"kAny");
681 hTriggerClassSelected->GetXaxis()->SetBinLabel(33,"V0AND");
682 hTriggerClassSelected->GetXaxis()->SetBinLabel(34,"NOT kFastOnly");
683 fHistograms->Add(hTriggerClassSelected);
685 TH1::AddDirectory(kTRUE);
688 //________________________________________________________________________
689 Bool_t AliConversionCuts::InitPIDResponse(){
690 // Set Pointer to AliPIDResponse
692 AliAnalysisManager *man=AliAnalysisManager::GetAnalysisManager();
694 AliInputEventHandler* inputHandler = (AliInputEventHandler*) (man->GetInputEventHandler());
695 fPIDResponse = (AliPIDResponse*)inputHandler->GetPIDResponse();
696 if(fPIDResponse)return kTRUE;
703 ///________________________________________________________________________
704 Bool_t AliConversionCuts::EventIsSelected(AliVEvent *fInputEvent, AliVEvent *fMCEvent){
705 // Process Event Selection
708 if(hV0EventCuts)hV0EventCuts->Fill(cutindex);
711 // Check for MC event
712 if(fMCEvent && fInputEvent->IsA()==AliESDEvent::Class()){
713 // Check if MC event is correctly loaded
714 AliMCEventHandler* mcHandler = dynamic_cast<AliMCEventHandler*> (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
719 if (!mcHandler->InitOk() ){
723 if (!mcHandler->TreeK() ){
727 if (!mcHandler->TreeTR() ) {
734 // cout << "before event trigger" << endl;
735 if(!IsTriggerSelected(fInputEvent)){
736 if(hV0EventCuts)hV0EventCuts->Fill(cutindex);
742 if(fInputEvent->IsA()==AliESDEvent::Class()){
743 AliTriggerAnalysis fTriggerAnalysis;// = new AliTriggerAnalysis;
744 fHasV0AND = fTriggerAnalysis.IsOfflineTriggerFired((AliESDEvent*)fInputEvent, AliTriggerAnalysis::kV0AND);
745 if(fHasV0AND&&hTriggerClass)hTriggerClass->Fill(32);
747 // cout << "event number " << ((AliESDEvent*)fInputEvent)->GetEventNumberInFile() << " entered"<< endl;
750 // Number of Contributors Cut
751 if(GetNumberOfContributorsVtx(fInputEvent)<=0) {
752 if(hV0EventCuts)hV0EventCuts->Fill(cutindex);
758 // Z Vertex Position Cut
759 if(!VertexZCut(fInputEvent)){
760 if(hV0EventCuts)hV0EventCuts->Fill(cutindex);
769 if(fInputEvent->IsPileupFromSPD(3,0.8,3.,2.,5.)){
770 if(hV0EventCuts)hV0EventCuts->Fill(cutindex);
777 // Centrality Selection
778 if(!IsCentralitySelected(fInputEvent,fMCEvent)){
779 if(hV0EventCuts)hV0EventCuts->Fill(cutindex);
785 // Fill Event Histograms
786 if(hV0EventCuts)hV0EventCuts->Fill(cutindex);
787 if(hVertexZ)hVertexZ->Fill(fInputEvent->GetPrimaryVertex()->GetZ());
788 if(hCentrality)hCentrality->Fill(GetCentrality(fInputEvent));
789 if(hCentralityVsNumberOfPrimaryTracks)
790 hCentralityVsNumberOfPrimaryTracks->Fill(GetCentrality(fInputEvent),
791 ((AliV0ReaderV1*)AliAnalysisManager::GetAnalysisManager()
792 ->GetTask("V0ReaderV1"))->GetNumberOfPrimaryTracks());
797 ///________________________________________________________________________
798 Bool_t AliConversionCuts::PhotonIsSelectedMC(TParticle *particle,AliStack *fMCStack,Bool_t checkForConvertedGamma){
799 // MonteCarlo Photon Selection
801 if(!fMCStack)return kFALSE;
803 if (particle->GetPdgCode() == 22){
806 if( particle->Eta() > (fEtaCut + fEtaShift) || particle->Eta() < (-fEtaCut + fEtaShift) )
809 if( particle->Eta() < (fEtaCutMin + fEtaShift) && particle->Eta() > (-fEtaCutMin + fEtaShift) )
813 if(particle->GetMother(0) >-1 && fMCStack->Particle(particle->GetMother(0))->GetPdgCode() == 22){
814 return kFALSE; // no photon as mothers!
817 if(particle->GetMother(0) >= fMCStack->GetNprimary()){
818 return kFALSE; // the gamma has a mother, and it is not a primary particle
821 if(!checkForConvertedGamma) return kTRUE; // return in case of accepted gamma
823 // looking for conversion gammas (electron + positron from pairbuilding (= 5) )
824 TParticle* ePos = NULL;
825 TParticle* eNeg = NULL;
827 if(particle->GetNDaughters() >= 2){
828 for(Int_t daughterIndex=particle->GetFirstDaughter();daughterIndex<=particle->GetLastDaughter();daughterIndex++){
829 TParticle *tmpDaughter = fMCStack->Particle(daughterIndex);
830 if(tmpDaughter->GetUniqueID() == 5){
831 if(tmpDaughter->GetPdgCode() == 11){
833 } else if(tmpDaughter->GetPdgCode() == -11){
840 if(ePos == NULL || eNeg == NULL){ // means we do not have two daughters from pair production
844 if(ePos->Pt()<fSinglePtCut || eNeg->Pt()<fSinglePtCut){
845 return kFALSE; // no reconstruction below the Pt cut
848 if( ePos->Eta() > (fEtaCut + fEtaShift) || ePos->Eta() < (-fEtaCut + fEtaShift) ||
849 eNeg->Eta() > (fEtaCut + fEtaShift) || eNeg->Eta() < (-fEtaCut + fEtaShift) )
852 if(fEtaCutMin > -0.1){
853 if( (ePos->Eta() < (fEtaCutMin + fEtaShift) && ePos->Eta() > (-fEtaCutMin + fEtaShift)) ||
854 (eNeg->Eta() < (fEtaCutMin + fEtaShift) && eNeg->Eta() > (-fEtaCutMin + fEtaShift)) )
859 return kFALSE; // cuts on distance from collision point
862 if(abs(ePos->Vz()) > fMaxZ){
863 return kFALSE; // outside material
865 if(abs(eNeg->Vz()) > fMaxZ){
866 return kFALSE; // outside material
869 if( ePos->R() <= ((abs(ePos->Vz()) * fLineCutZRSlope) - fLineCutZValue)){
870 return kFALSE; // line cut to exclude regions where we do not reconstruct
871 } else if ( fEtaCutMin != -0.1 && ePos->R() >= ((abs(ePos->Vz()) * fLineCutZRSlopeMin) - fLineCutZValueMin)){
875 if( eNeg->R() <= ((abs(eNeg->Vz()) * fLineCutZRSlope) - fLineCutZValue)){
876 return kFALSE; // line cut to exclude regions where we do not reconstruct
877 } else if ( fEtaCutMin != -0.1 && eNeg->R() >= ((abs(eNeg->Vz()) * fLineCutZRSlopeMin) - fLineCutZValueMin)){
882 //if(AcceptanceCut(particle,ePos,eNeg))return kTRUE;
886 ///________________________________________________________________________
887 Bool_t AliConversionCuts::PhotonIsSelectedAODMC(AliAODMCParticle *particle,TClonesArray *aodmcArray,Bool_t checkForConvertedGamma){
888 // MonteCarlo Photon Selection
890 if(!aodmcArray)return kFALSE;
892 if (particle->GetPdgCode() == 22){
893 if( particle->Eta() > (fEtaCut + fEtaShift) || particle->Eta() < (-fEtaCut + fEtaShift) )
896 if( particle->Eta() < (fEtaCutMin + fEtaShift) && particle->Eta() > (-fEtaCutMin + fEtaShift) )
900 if(particle->GetMother() > -1){
901 if((static_cast<AliAODMCParticle*>(aodmcArray->At(particle->GetMother())))->GetPdgCode() == 22){
902 return kFALSE; // no photon as mothers!
904 if(!(static_cast<AliAODMCParticle*>(aodmcArray->At(particle->GetMother()))->IsPrimary())){
905 return kFALSE; // the gamma has a mother, and it is not a primary particle
909 if(!checkForConvertedGamma) return kTRUE; // return in case of accepted gamma
911 // looking for conversion gammas (electron + positron from pairbuilding (= 5) )
912 AliAODMCParticle* ePos = NULL;
913 AliAODMCParticle* eNeg = NULL;
915 if(particle->GetNDaughters() >= 2){
916 for(Int_t daughterIndex=particle->GetDaughter(0);daughterIndex<=particle->GetDaughter(1);daughterIndex++){
917 AliAODMCParticle *tmpDaughter = static_cast<AliAODMCParticle*>(aodmcArray->At(daughterIndex));
918 if(!tmpDaughter) continue;
919 if(((tmpDaughter->GetMCProcessCode())) == 5){ // STILL A BUG IN ALIROOT >>8 HAS TPO BE REMOVED AFTER FIX
920 if(tmpDaughter->GetPdgCode() == 11){
922 } else if(tmpDaughter->GetPdgCode() == -11){
929 if(ePos == NULL || eNeg == NULL){ // means we do not have two daughters from pair production
933 if(ePos->Pt()<fSinglePtCut || eNeg->Pt()<fSinglePtCut){
934 return kFALSE; // no reconstruction below the Pt cut
937 if( ePos->Eta() > (fEtaCut + fEtaShift) || ePos->Eta() < (-fEtaCut + fEtaShift) ||
938 eNeg->Eta() > (fEtaCut + fEtaShift) || eNeg->Eta() < (-fEtaCut + fEtaShift) )
941 if(fEtaCutMin > -0.1){
942 if( (ePos->Eta() < (fEtaCutMin + fEtaShift) && ePos->Eta() > (-fEtaCutMin + fEtaShift)) ||
943 (eNeg->Eta() < (fEtaCutMin + fEtaShift) && eNeg->Eta() > (-fEtaCutMin + fEtaShift)) )
947 Double_t rPos = sqrt( (ePos->Xv()*ePos->Xv()) + (ePos->Yv()*ePos->Yv()) );
948 Double_t rNeg = sqrt( (eNeg->Xv()*eNeg->Xv()) + (eNeg->Yv()*eNeg->Yv()) );
951 return kFALSE; // cuts on distance from collision point
953 if(abs(ePos->Zv()) > fMaxZ){
954 return kFALSE; // outside material
956 if(abs(eNeg->Zv()) > fMaxZ){
957 return kFALSE; // outside material
960 if( rPos <= ((abs(ePos->Zv()) * fLineCutZRSlope) - fLineCutZValue)){
961 return kFALSE; // line cut to exclude regions where we do not reconstruct
962 } else if ( fEtaCutMin != -0.1 && rPos >= ((abs(ePos->Zv()) * fLineCutZRSlopeMin) - fLineCutZValueMin)){
966 if( rNeg <= ((abs(eNeg->Zv()) * fLineCutZRSlope) - fLineCutZValue)){
967 return kFALSE; // line cut to exclude regions where we do not reconstruct
968 } else if ( fEtaCutMin != -0.1 && rNeg >= ((abs(eNeg->Zv()) * fLineCutZRSlopeMin) - fLineCutZValueMin)){
973 //if(AcceptanceCut(particle,ePos,eNeg))return kTRUE;
977 ///________________________________________________________________________
978 Bool_t AliConversionCuts::PhotonCuts(AliConversionPhotonBase *photon,AliVEvent *event)
979 { // Specific Photon Cuts
982 if(hPhotonCuts)hPhotonCuts->Fill(cutIndex);
985 // Fill Histos before Cuts
986 if(hInvMassbefore)hInvMassbefore->Fill(photon->GetMass());
988 if(hArmenterosbefore)hArmenterosbefore->Fill(photon->GetArmenterosAlpha(),photon->GetArmenterosQt());
990 // Gamma selection based on QT from Armenteros
991 if(fDoQtGammaSelection == kTRUE){
992 if(!ArmenterosQtCut(photon)){
993 if(hPhotonCuts)hPhotonCuts->Fill(cutIndex); //1
1000 if(photon->GetChi2perNDF() > fChi2CutConversion || photon->GetChi2perNDF() <=0){
1002 if(hPhotonCuts)hPhotonCuts->Fill(cutIndex); //2
1008 // Reconstruction Acceptance Cuts
1009 if(!AcceptanceCuts(photon)){
1010 if(hPhotonCuts)hPhotonCuts->Fill(cutIndex); //3
1016 if(fDoPhotonAsymmetryCut == kTRUE){
1017 if(!AsymmetryCut(photon,event)){
1018 if(hPhotonCuts)hPhotonCuts->Fill(cutIndex); //4
1023 //Check the pid probability
1025 if(!PIDProbabilityCut(photon, event)) {
1026 if(hPhotonCuts)hPhotonCuts->Fill(cutIndex); //5
1031 if(!CorrectedTPCClusterCut(photon, event)) {
1032 if(hPhotonCuts)hPhotonCuts->Fill(cutIndex); //6
1038 if(!PsiPairCut(photon)) {
1039 if(hPhotonCuts)hPhotonCuts->Fill(cutIndex); //7
1044 if(!CosinePAngleCut(photon, event)) {
1045 if(hPhotonCuts)hPhotonCuts->Fill(cutIndex); //8
1049 AliAODConversionPhoton* photonAOD = dynamic_cast<AliAODConversionPhoton*>(photon);
1051 photonAOD->CalculateDistanceOfClossetApproachToPrimVtx(event->GetPrimaryVertex());
1054 if(photonAOD->GetDCArToPrimVtx() > fDCARPrimVtxCut) { //DCA R cut of photon to primary vertex
1055 if(hPhotonCuts)hPhotonCuts->Fill(cutIndex); //9
1060 if(abs(photonAOD->GetDCAzToPrimVtx()) > fDCAZPrimVtxCut) { //DCA Z cut of photon to primary vertex
1061 if(hPhotonCuts)hPhotonCuts->Fill(cutIndex); //10
1069 if(hPhotonCuts)hPhotonCuts->Fill(cutIndex); //11
1071 // Histos after Cuts
1072 if(hInvMassafter)hInvMassafter->Fill(photon->GetMass());
1073 if(hArmenterosafter)hArmenterosafter->Fill(photon->GetArmenterosAlpha(),photon->GetArmenterosQt());
1079 ///________________________________________________________________________
1080 Bool_t AliConversionCuts::CorrectedTPCClusterCut(AliConversionPhotonBase *photon, AliVEvent * event)
1081 { //Cut on corrected TPC Cluster Info
1083 AliVTrack * negTrack = GetTrack(event, photon->GetTrackLabelNegative());
1084 AliVTrack * posTrack = GetTrack(event, photon->GetTrackLabelPositive());
1086 if(!negTrack||!posTrack)return kFALSE;
1088 Double_t negclsToF=0;
1090 if (!fUseCorrectedTPCClsInfo ){
1091 if(negTrack->GetTPCNclsF()!=0){
1092 negclsToF = (Double_t)negTrack->GetNcls(1)/(Double_t)negTrack->GetTPCNclsF();}// Ncluster/Nfindablecluster
1095 negclsToF = negTrack->GetTPCClusterInfo(2,0,GetFirstTPCRow(photon->GetConversionRadius()));
1098 Double_t posclsToF = 0.;
1099 if (!fUseCorrectedTPCClsInfo ){
1100 if(posTrack->GetTPCNclsF()!=0){
1101 posclsToF = (Double_t)posTrack->GetNcls(1)/(Double_t)posTrack->GetTPCNclsF();
1104 posclsToF = posTrack->GetTPCClusterInfo(2,0,GetFirstTPCRow(photon->GetConversionRadius()));
1107 if( negclsToF < fMinClsTPCToF || posclsToF < fMinClsTPCToF ){
1114 ///________________________________________________________________________
1115 Bool_t AliConversionCuts::PhotonIsSelected(AliConversionPhotonBase *photon, AliVEvent * event)
1117 //Selection of Reconstructed Photons
1119 FillPhotonCutIndex(kPhotonIn);
1121 if(event->IsA()==AliESDEvent::Class()) {
1122 if(!SelectV0Finder( ( ((AliESDEvent*)event)->GetV0(photon->GetV0Index()))->GetOnFlyStatus() ) ){
1123 FillPhotonCutIndex(kOnFly);
1127 // else if(event->IsA()==AliAODEvent::Class()) {
1128 // if(!SelectV0Finder( ( ((AliAODEvent*)event)->GetV0(photon->GetV0Index())) ) ){
1129 // FillPhotonCutIndex(kOnFly);
1135 AliVTrack * negTrack = GetTrack(event, photon->GetTrackLabelNegative());
1136 AliVTrack * posTrack = GetTrack(event, photon->GetTrackLabelPositive());
1138 if(!negTrack || !posTrack) {
1139 FillPhotonCutIndex(kNoTracks);
1142 photon->DeterminePhotonQuality(negTrack,posTrack);
1144 if(!TracksAreSelected(negTrack, posTrack)){
1145 FillPhotonCutIndex(kTrackCuts);
1150 if(!dEdxCuts(negTrack) || !dEdxCuts(posTrack)) {
1151 FillPhotonCutIndex(kdEdxCuts);
1156 if(!PhotonCuts(photon,event)){
1157 FillPhotonCutIndex(kPhotonCuts);
1161 // Photon passed cuts
1162 FillPhotonCutIndex(kPhotonOut);
1166 ///________________________________________________________________________
1167 Bool_t AliConversionCuts::ArmenterosQtCut(AliConversionPhotonBase *photon)
1168 { // Armenteros Qt Cut
1170 if(fDoHighPtQtGammaSelection){
1171 if(photon->GetPhotonPt() < fPtBorderForQt){
1172 if(photon->GetArmenterosQt()>fQtMax){
1176 if(photon->GetArmenterosQt()>fHighPtQtMax){
1182 if(photon->GetArmenterosQt()>fQtMax){
1190 ///________________________________________________________________________
1191 Bool_t AliConversionCuts::AcceptanceCuts(AliConversionPhotonBase *photon) {
1192 // Exclude certain areas for photon reconstruction
1195 if(hAcceptanceCuts)hAcceptanceCuts->Fill(cutIndex);
1198 if(photon->GetConversionRadius()>fMaxR){ // cuts on distance from collision point
1199 if(hAcceptanceCuts)hAcceptanceCuts->Fill(cutIndex);
1204 if(photon->GetConversionRadius()<fMinR){ // cuts on distance from collision point
1205 if(hAcceptanceCuts)hAcceptanceCuts->Fill(cutIndex);
1210 if(photon->GetConversionRadius() <= ((abs(photon->GetConversionZ())*fLineCutZRSlope)-fLineCutZValue)){
1211 if(hAcceptanceCuts)hAcceptanceCuts->Fill(cutIndex);
1214 else if (fUseEtaMinCut && photon->GetConversionRadius() >= ((abs(photon->GetConversionZ())*fLineCutZRSlopeMin)-fLineCutZValueMin )){
1215 if(hAcceptanceCuts)hAcceptanceCuts->Fill(cutIndex);
1220 if(abs(photon->GetConversionZ()) > fMaxZ ){ // cuts out regions where we do not reconstruct
1221 if(hAcceptanceCuts)hAcceptanceCuts->Fill(cutIndex);
1227 if( photon->GetPhotonEta() > (fEtaCut + fEtaShift) || photon->GetPhotonEta() < (-fEtaCut + fEtaShift) ){
1228 if(hAcceptanceCuts)hAcceptanceCuts->Fill(cutIndex);
1231 if(fEtaCutMin>-0.1){
1232 if( photon->GetPhotonEta() < (fEtaCutMin + fEtaShift) && photon->GetPhotonEta() > (-fEtaCutMin + fEtaShift) ){
1233 if(hAcceptanceCuts)hAcceptanceCuts->Fill(cutIndex);
1239 if(photon->GetPhotonPt()<fPtCut){
1240 if(hAcceptanceCuts)hAcceptanceCuts->Fill(cutIndex);
1245 if(hAcceptanceCuts)hAcceptanceCuts->Fill(cutIndex);
1251 ///________________________________________________________________________
1252 Bool_t AliConversionCuts::SpecificTrackCuts(AliAODTrack * negTrack, AliAODTrack * posTrack,Int_t &cutIndex) {
1253 // Track Cuts which require AOD/ESD specific implementation
1255 if( !negTrack->IsOn(AliESDtrack::kTPCrefit) || !posTrack->IsOn(AliESDtrack::kTPCrefit) ) {
1256 if(hTrackCuts)hTrackCuts->Fill(cutIndex);
1261 AliAODVertex * NegVtxType=negTrack->GetProdVertex();
1262 AliAODVertex * PosVtxType=posTrack->GetProdVertex();
1263 if( (NegVtxType->GetType())==AliAODVertex::kKink || (PosVtxType->GetType())==AliAODVertex::kKink) {
1264 if(hTrackCuts)hTrackCuts->Fill(cutIndex);
1272 ///________________________________________________________________________
1273 Bool_t AliConversionCuts::SpecificTrackCuts(AliESDtrack * negTrack, AliESDtrack * posTrack,Int_t &cutIndex) {
1274 // Track Cuts which require AOD/ESD specific implementation
1276 if( !negTrack->IsOn(AliESDtrack::kTPCrefit) || !posTrack->IsOn(AliESDtrack::kTPCrefit) ) {
1277 if(hTrackCuts)hTrackCuts->Fill(cutIndex);
1282 if(negTrack->GetKinkIndex(0) > 0 || posTrack->GetKinkIndex(0) > 0 ) {
1283 if(hTrackCuts)hTrackCuts->Fill(cutIndex);
1291 ///________________________________________________________________________
1292 Bool_t AliConversionCuts::TracksAreSelected(AliVTrack * negTrack, AliVTrack * posTrack) {
1293 // Track Selection for Photon Reconstruction
1296 if(hTrackCuts)hTrackCuts->Fill(cutIndex);
1300 if(negTrack->Charge() == posTrack->Charge()) {
1301 if(hTrackCuts)hTrackCuts->Fill(cutIndex);
1306 // Number of TPC Clusters
1309 if( negTrack->GetNcls(1) < fMinClsTPC || posTrack->GetNcls(1) < fMinClsTPC ) {
1310 if(hTrackCuts)hTrackCuts->Fill(cutIndex);
1316 if( posTrack->Eta() > (fEtaCut + fEtaShift) || posTrack->Eta() < (-fEtaCut + fEtaShift) ||
1317 negTrack->Eta() > (fEtaCut + fEtaShift) || negTrack->Eta() < (-fEtaCut + fEtaShift) ){
1318 if(hTrackCuts)hTrackCuts->Fill(cutIndex);
1321 if(fEtaCutMin>-0.1){
1322 if( (posTrack->Eta() < (fEtaCutMin + fEtaShift) && posTrack->Eta() > (-fEtaCutMin + fEtaShift)) ||
1323 (negTrack->Eta() < (fEtaCutMin + fEtaShift) && negTrack->Eta() > (-fEtaCutMin + fEtaShift)) ){
1324 if(hTrackCuts)hTrackCuts->Fill(cutIndex);
1331 if( negTrack->Pt()< fSinglePtCut || posTrack->Pt()< fSinglePtCut){
1332 if(hTrackCuts)hTrackCuts->Fill(cutIndex);
1337 // AOD ESD specific cuts
1338 Bool_t passCuts = kTRUE;
1340 if(negTrack->IsA()==AliAODTrack::Class()) {
1341 passCuts = passCuts * SpecificTrackCuts(static_cast<AliAODTrack*>(negTrack), static_cast<AliAODTrack*>(posTrack),cutIndex);
1343 passCuts = passCuts * SpecificTrackCuts(static_cast<AliESDtrack*>(negTrack), static_cast<AliESDtrack*>(posTrack),cutIndex);
1347 if(hTrackCuts)hTrackCuts->Fill(cutIndex);
1352 if(hTrackCuts)hTrackCuts->Fill(cutIndex);
1358 ///________________________________________________________________________
1359 Bool_t AliConversionCuts::dEdxCuts(AliVTrack *fCurrentTrack){
1360 // Electron Identification Cuts for Photon reconstruction
1362 if(!fPIDResponse){InitPIDResponse();}// Try to reinitialize PID Response
1363 if(!fPIDResponse){AliError("No PID Response"); return kTRUE;}// if still missing fatal error
1366 if(hdEdxCuts)hdEdxCuts->Fill(cutIndex);
1367 if(hTPCdEdxSigbefore)hTPCdEdxSigbefore->Fill(fCurrentTrack->P(),fPIDResponse->NumberOfSigmasTPC(fCurrentTrack, AliPID::kElectron));
1368 if(hTPCdEdxbefore)hTPCdEdxbefore->Fill(fCurrentTrack->P(),fCurrentTrack->GetTPCsignal());
1372 if(fDodEdxSigmaCut == kTRUE){
1373 // TPC Electron Line
1374 if( fPIDResponse->NumberOfSigmasTPC(fCurrentTrack,AliPID::kElectron)<fPIDnSigmaBelowElectronLine ||
1375 fPIDResponse->NumberOfSigmasTPC(fCurrentTrack,AliPID::kElectron)>fPIDnSigmaAboveElectronLine){
1377 if(hdEdxCuts)hdEdxCuts->Fill(cutIndex);
1383 if( fCurrentTrack->P()>fPIDMinPnSigmaAbovePionLine && fCurrentTrack->P()<fPIDMaxPnSigmaAbovePionLine ){
1384 if(fPIDResponse->NumberOfSigmasTPC(fCurrentTrack,AliPID::kElectron)>fPIDnSigmaBelowElectronLine &&
1385 fPIDResponse->NumberOfSigmasTPC(fCurrentTrack,AliPID::kElectron)<fPIDnSigmaAboveElectronLine&&
1386 fPIDResponse->NumberOfSigmasTPC(fCurrentTrack,AliPID::kPion)<fPIDnSigmaAbovePionLine){
1388 if(hdEdxCuts)hdEdxCuts->Fill(cutIndex);
1395 if( fCurrentTrack->P()>fPIDMaxPnSigmaAbovePionLine ){
1396 if(fPIDResponse->NumberOfSigmasTPC(fCurrentTrack,AliPID::kElectron)>fPIDnSigmaBelowElectronLine &&
1397 fPIDResponse->NumberOfSigmasTPC(fCurrentTrack,AliPID::kElectron)<fPIDnSigmaAboveElectronLine &&
1398 fPIDResponse->NumberOfSigmasTPC(fCurrentTrack,AliPID::kPion)<fPIDnSigmaAbovePionLineHighPt){
1400 if(hdEdxCuts)hdEdxCuts->Fill(cutIndex);
1408 if(fDoKaonRejectionLowP == kTRUE){
1409 if(fCurrentTrack->P()<fPIDMinPKaonRejectionLowP ){
1410 if( abs(fPIDResponse->NumberOfSigmasTPC(fCurrentTrack,AliPID::kKaon))<fPIDnSigmaAtLowPAroundKaonLine){
1412 if(hdEdxCuts)hdEdxCuts->Fill(cutIndex);
1419 if(fDoProtonRejectionLowP == kTRUE){
1420 if( fCurrentTrack->P()<fPIDMinPProtonRejectionLowP ){
1421 if( abs(fPIDResponse->NumberOfSigmasTPC(fCurrentTrack,AliPID::kProton))<fPIDnSigmaAtLowPAroundProtonLine){
1423 if(hdEdxCuts)hdEdxCuts->Fill(cutIndex);
1430 if(fDoPionRejectionLowP == kTRUE){
1431 if( fCurrentTrack->P()<fPIDMinPPionRejectionLowP ){
1432 if( abs(fPIDResponse->NumberOfSigmasTPC(fCurrentTrack,AliPID::kPion))<fPIDnSigmaAtLowPAroundPionLine){
1434 if(hdEdxCuts)hdEdxCuts->Fill(cutIndex);
1442 // cout<<"Start"<<endl;
1443 // AliPIDResponse::EDetPidStatus status=fPIDResponse->CheckPIDStatus(AliPIDResponse::kTOF,fCurrentTrack);
1445 // if( ( (status & AliVTrack::kTOFout) == AliVTrack::kTOFout ) && ( (status & AliVTrack::kTIME) == AliVTrack::kTIME ))
1446 // {cout<<"TOF DA"<<endl;}
1447 // if(status == AliPIDResponse::kDetPidOk){
1448 // Float_t probMis = fPIDResponse->GetTOFMismatchProbability(fCurrentTrack);
1449 // cout<<"--> "<<probMis<<endl;
1450 // if(probMis > 0.01){
1455 if((fCurrentTrack->GetStatus() & AliESDtrack::kTOFpid) && !(fCurrentTrack->GetStatus() & AliESDtrack::kTOFmismatch)){
1457 Double_t t0 = fPIDResponse->GetTOFResponse().GetStartTime(fCurrentTrack->P());
1459 fCurrentTrack->GetIntegratedTimes(times);
1460 Double_t TOFsignal = fCurrentTrack->GetTOFsignal();
1461 Double_t dT = TOFsignal - t0 - times[0];
1462 hTOFbefore->Fill(fCurrentTrack->P(),dT);
1464 if(hTOFSigbefore) hTOFSigbefore->Fill(fCurrentTrack->P(),fPIDResponse->NumberOfSigmasTOF(fCurrentTrack, AliPID::kElectron));
1466 if(fPIDResponse->NumberOfSigmasTOF(fCurrentTrack, AliPID::kElectron)>fTofPIDnSigmaAboveElectronLine ||
1467 fPIDResponse->NumberOfSigmasTOF(fCurrentTrack, AliPID::kElectron)<fTofPIDnSigmaBelowElectronLine ){
1468 if(hdEdxCuts)hdEdxCuts->Fill(cutIndex);
1472 if(hTOFSigafter)hTOFSigafter->Fill(fCurrentTrack->P(),fPIDResponse->NumberOfSigmasTOF(fCurrentTrack, AliPID::kElectron));
1478 if(!fPIDResponse->IdentifiedAsElectronTRD(fCurrentTrack,fPIDTRDEfficiency)){
1479 if(hdEdxCuts)hdEdxCuts->Fill(cutIndex);
1485 if(hdEdxCuts)hdEdxCuts->Fill(cutIndex);
1486 if(hTPCdEdxSigafter)hTPCdEdxSigafter->Fill(fCurrentTrack->P(),fPIDResponse->NumberOfSigmasTPC(fCurrentTrack, AliPID::kElectron));
1487 if(hTPCdEdxafter)hTPCdEdxafter->Fill(fCurrentTrack->P(),fCurrentTrack->GetTPCsignal());
1491 ///________________________________________________________________________
1492 Bool_t AliConversionCuts::AsymmetryCut(AliConversionPhotonBase * photon,AliVEvent *event) {
1493 // Cut on Energy Assymetry
1495 for(Int_t ii=0;ii<2;ii++){
1497 AliVTrack *track=GetTrack(event,photon->GetTrackLabel(ii));
1499 if( track->P() > fMinPPhotonAsymmetryCut ){
1500 Double_t trackNegAsy=0;
1501 if (photon->GetPhotonP()!=0.){
1502 trackNegAsy= track->P()/photon->GetPhotonP();
1505 if( trackNegAsy<fMinPhotonAsymmetry ||trackNegAsy>(1.- fMinPhotonAsymmetry)){
1514 ///________________________________________________________________________
1515 AliVTrack *AliConversionCuts::GetTrack(AliVEvent * event, Int_t label){
1516 //Returns pointer to the track with given ESD label
1517 //(Important for AOD implementation, since Track array in AOD data is different
1518 //from ESD array, but ESD tracklabels are stored in AOD Tracks)
1520 AliESDEvent * esdEvent = dynamic_cast<AliESDEvent*>(event);
1522 if(label > event->GetNumberOfTracks() ) return NULL;
1523 AliESDtrack * track = esdEvent->GetTrack(label);
1527 AliVTrack * track = 0x0;
1528 if(((AliV0ReaderV1*)AliAnalysisManager::GetAnalysisManager()->GetTask("V0ReaderV1"))->AreAODsRelabeled()){
1529 track = dynamic_cast<AliVTrack*>(event->GetTrack(label));
1533 for(Int_t ii=0; ii<event->GetNumberOfTracks(); ii++) {
1534 track = dynamic_cast<AliVTrack*>(event->GetTrack(ii));
1536 if(track->GetID() == label) {
1543 //AliDebug(5,(Form("track not found %d %d",label,event->GetNumberOfTracks()));
1547 ///________________________________________________________________________
1548 AliESDtrack *AliConversionCuts::GetESDTrack(AliESDEvent * event, Int_t label){
1549 //Returns pointer to the track with given ESD label
1550 //(Important for AOD implementation, since Track array in AOD data is different
1551 //from ESD array, but ESD tracklabels are stored in AOD Tracks)
1554 if(label > event->GetNumberOfTracks() ) return NULL;
1555 AliESDtrack * track = event->GetTrack(label);
1558 //AliDebug(5,(Form("track not found %d %d",label,event->GetNumberOfTracks()));
1564 ///________________________________________________________________________
1565 Bool_t AliConversionCuts::PIDProbabilityCut(AliConversionPhotonBase *photon, AliVEvent * event){
1566 // Cut on Electron Probability for Photon Reconstruction
1568 AliESDEvent * esdEvent = dynamic_cast<AliESDEvent*>(event);
1572 Bool_t iResult=kFALSE;
1574 Double_t *posProbArray = new Double_t[AliPID::kSPECIES];
1575 Double_t *negProbArray = new Double_t[AliPID::kSPECIES];
1577 AliESDtrack* negTrack = esdEvent->GetTrack(photon->GetTrackLabelNegative());
1578 AliESDtrack* posTrack = esdEvent->GetTrack(photon->GetTrackLabelPositive());
1580 if(negProbArray && posProbArray){
1582 negTrack->GetTPCpid(negProbArray);
1583 posTrack->GetTPCpid(posProbArray);
1585 if(negProbArray[AliPID::kElectron]>=fPIDProbabilityCutNegativeParticle && posProbArray[AliPID::kElectron]>=fPIDProbabilityCutPositiveParticle){
1590 delete [] posProbArray;
1591 delete [] negProbArray;
1595 ///Not possible for AODs
1605 ///________________________________________________________________________
1606 Bool_t AliConversionCuts::AcceptanceCut(TParticle *particle, TParticle * ePos,TParticle* eNeg){
1607 // MC Acceptance Cuts
1608 //(Certain areas were excluded for photon reconstruction)
1610 if(particle->R()>fMaxR){
1613 if(ePos->R()>fMaxR){
1617 if(ePos->R()<fMinR){
1621 if( ePos->R() <= ((abs(ePos->Vz())*fLineCutZRSlope)-fLineCutZValue)){
1624 else if (fUseEtaMinCut && ePos->R() >= ((abs(ePos->Vz())*fLineCutZRSlopeMin)-fLineCutZValueMin )){
1628 if(abs(eNeg->Vz()) > fMaxZ ){ // cuts out regions where we do not reconstruct
1632 if(eNeg->Vz()!=ePos->Vz()||eNeg->R()!=ePos->R()){
1636 if(abs(ePos->Vz()) > fMaxZ ){ // cuts out regions where we do not reconstruct
1641 if( particle->Eta() > (fEtaCut + fEtaShift) || particle->Eta() < (-fEtaCut + fEtaShift) ){
1644 if( ePos->Eta() > (fEtaCut + fEtaShift) || ePos->Eta() < (-fEtaCut + fEtaShift) ){
1647 if( eNeg->Eta() > (fEtaCut + fEtaShift) || eNeg->Eta() < (-fEtaCut + fEtaShift) ){
1650 if(fEtaCutMin>-0.1){
1651 if( particle->Eta() < (fEtaCutMin + fEtaShift) && particle->Eta() > (-fEtaCutMin + fEtaShift) ){
1654 if( ePos->Eta() < (fEtaCutMin + fEtaShift) && ePos->Eta() > (-fEtaCutMin + fEtaShift) ){
1657 if( eNeg->Eta() < (fEtaCutMin + fEtaShift) && eNeg->Eta() > (-fEtaCutMin + fEtaShift) ){
1662 if( ePos->Pt()< fSinglePtCut || eNeg->Pt()< fSinglePtCut){
1666 if(particle->Pt()<fPtCut){
1672 ///________________________________________________________________________
1673 Bool_t AliConversionCuts::UpdateCutString() {
1674 ///Update the cut string (if it has been created yet)
1676 if(fCutString && fCutString->GetString().Length() == kNCuts) {
1677 fCutString->SetString(GetCutNumber());
1683 ///________________________________________________________________________
1684 void AliConversionCuts::LoadReweightingHistosMCFromFile() {
1686 AliInfo("Entering loading of histograms for weighting");
1687 TFile *f = TFile::Open(fPathTrFReweighting.Data());
1689 AliError(Form("file for weighting %s not found",fPathTrFReweighting.Data()));
1692 if (fNameHistoReweightingPi0.CompareTo("") != 0 && fDoReweightHistoMCPi0 ){
1693 TH1D *hReweightMCHistPi0temp = (TH1D*)f->Get(fNameHistoReweightingPi0.Data());
1694 hReweightMCHistPi0 = new TH1D(*hReweightMCHistPi0temp);
1695 if (hReweightMCHistPi0) AliInfo(Form("%s has been loaded from %s", fNameHistoReweightingPi0.Data(),fPathTrFReweighting.Data() ));
1696 else AliWarning(Form("%s not found in %s", fNameHistoReweightingPi0.Data() ,fPathTrFReweighting.Data()));
1698 if (fNameFitDataPi0.CompareTo("") != 0 && fDoReweightHistoMCPi0 ){
1699 TF1 *fFitDataPi0temp = (TF1*)f->Get(fNameFitDataPi0.Data());
1700 fFitDataPi0 = new TF1(*fFitDataPi0temp);
1701 if (fFitDataPi0) AliInfo(Form("%s has been loaded from %s", fNameFitDataPi0.Data(),fPathTrFReweighting.Data() ));
1702 else AliWarning(Form("%s not found in %s",fPathTrFReweighting.Data(), fNameFitDataPi0.Data() ));
1705 if (fNameHistoReweightingEta.CompareTo("") != 0 && fDoReweightHistoMCEta){
1706 TH1D *hReweightMCHistEtatemp = (TH1D*)f->Get(fNameHistoReweightingEta.Data());
1707 hReweightMCHistEta = new TH1D(*hReweightMCHistEtatemp);
1708 if (hReweightMCHistEta) AliInfo(Form("%s has been loaded from %s", fNameHistoReweightingEta.Data(),fPathTrFReweighting.Data() ));
1709 else AliWarning(Form("%s not found in %s", fNameHistoReweightingEta.Data(),fPathTrFReweighting.Data() ));
1712 if (fNameFitDataEta.CompareTo("") != 0 && fDoReweightHistoMCEta){
1713 TF1 *fFitDataEtatemp = (TF1*)f->Get(fNameFitDataEta.Data());
1714 fFitDataEta = new TF1(*fFitDataEtatemp);
1715 if (fFitDataEta) AliInfo(Form("%s has been loaded from %s", fNameFitDataEta.Data(),fPathTrFReweighting.Data() ));
1716 else AliWarning(Form("%s not found in %s", fNameFitDataEta.Data(),fPathTrFReweighting.Data() ));
1719 if (fNameHistoReweightingK0s.CompareTo("") != 0 && fDoReweightHistoMCK0s){
1720 TH1D *hReweightMCHistK0stemp = (TH1D*)f->Get(fNameHistoReweightingK0s.Data());
1721 hReweightMCHistK0s = new TH1D(*hReweightMCHistK0stemp);
1722 if (hReweightMCHistK0s) AliInfo(Form("%s has been loaded from %s", fNameHistoReweightingK0s.Data(),fPathTrFReweighting.Data() ));
1723 else AliWarning(Form("%s not found in %s", fNameHistoReweightingK0s.Data(),fPathTrFReweighting.Data() ));
1726 if (fNameFitDataK0s.CompareTo("") != 0 && fDoReweightHistoMCK0s){
1727 TF1 *fFitDataK0stemp = (TF1*)f->Get(fNameFitDataK0s.Data());
1728 fFitDataK0s = new TF1(*fFitDataK0stemp);
1729 if (fFitDataK0s) AliInfo(Form("%s has been loaded from %s", fNameFitDataK0s.Data(),fPathTrFReweighting.Data() ));
1730 else AliWarning(Form("%s not found in %s", fNameFitDataK0s.Data(),fPathTrFReweighting.Data() ));
1736 ///________________________________________________________________________
1737 Bool_t AliConversionCuts::InitializeCutsFromCutString(const TString analysisCutSelection ) {
1738 // Initialize Cuts from a given Cut string
1739 if(fDoReweightHistoMCPi0 || fDoReweightHistoMCEta || fDoReweightHistoMCK0s) {
1740 AliInfo("Weighting was enabled");
1741 LoadReweightingHistosMCFromFile();
1744 AliInfo(Form("Set Photoncut Number: %s",analysisCutSelection.Data()));
1745 if(analysisCutSelection.Length()!=kNCuts) {
1746 AliError(Form("Cut selection has the wrong length! size is %d, number of cuts is %d", analysisCutSelection.Length(), kNCuts));
1749 if(!analysisCutSelection.IsDigit()){
1750 AliError("Cut selection contains characters");
1754 const char *cutSelection = analysisCutSelection.Data();
1755 #define ASSIGNARRAY(i) fCuts[i] = cutSelection[i] - '0'
1756 for(Int_t ii=0;ii<kNCuts;ii++){
1760 // Set Individual Cuts
1761 for(Int_t ii=0;ii<kNCuts;ii++){
1762 if(!SetCut(cutIds(ii),fCuts[ii]))return kFALSE;
1769 ///________________________________________________________________________
1770 Bool_t AliConversionCuts::SetCut(cutIds cutID, const Int_t value) {
1771 ///Set individual cut ID
1776 if( SetV0Finder(value)) {
1777 fCuts[kv0FinderType] = value;
1780 } else return kFALSE;
1782 case kededxSigmaCut:
1783 if( SetTPCdEdxCutElectronLine(value)) {
1784 fCuts[kededxSigmaCut] = value;
1787 } else return kFALSE;
1789 case kpidedxSigmaCut:
1790 if( SetTPCdEdxCutPionLine(value)) {
1791 fCuts[kpidedxSigmaCut] = value;
1794 } else return kFALSE;
1796 case kpiMomdedxSigmaCut:
1797 if( SetMinMomPiondEdxCut(value)) {
1798 fCuts[kpiMomdedxSigmaCut] = value;
1801 } else return kFALSE;
1804 if( SetChi2GammaCut(value)) {
1805 fCuts[kchi2GammaCut] = value;
1808 } else return kFALSE;
1811 if( SetSinglePtCut(value)) {
1812 fCuts[ksinglePtCut] = value;
1815 } else return kFALSE;
1818 if( SetTPCClusterCut(value)) {
1819 fCuts[kclsTPCCut] = value;
1822 } else return kFALSE;
1825 if( SetEtaCut(value)) {
1826 fCuts[ketaCut] = value;
1829 } else return kFALSE;
1831 case kLowPRejectionSigmaCut:
1832 if( SetLowPRejectionCuts(value)) {
1833 fCuts[kLowPRejectionSigmaCut] = value;
1836 } else return kFALSE;
1839 if( SetQtMaxCut(value)) {
1840 fCuts[kQtMaxCut] = value;
1843 } else return kFALSE;
1845 case kpiMaxMomdedxSigmaCut:
1846 if( SetMaxMomPiondEdxCut(value)) {
1847 fCuts[kpiMaxMomdedxSigmaCut] = value;
1850 } else return kFALSE;
1853 if( SetRCut(value)) {
1854 fCuts[kRCut] = value;
1857 } else return kFALSE;
1860 if( SetRemovePileUp(value)) {
1861 fCuts[kremovePileUp] = value;
1864 } else return kFALSE;
1867 if( SetSelectSpecialTrigger(value)) {
1868 fCuts[kselectV0AND] = value;
1871 } else return kFALSE;
1873 case kmultiplicityMethod:
1874 if( SetMultiplicityMethod(value)) {
1875 fCuts[kmultiplicityMethod] = value;
1878 } else return kFALSE;
1881 if( SetIsHeavyIon(value)) {
1882 fCuts[kisHeavyIon] = value;
1885 } else return kFALSE;
1887 case kCentralityMin:
1888 if( SetCentralityMin(value)) {
1889 fCuts[kCentralityMin] = value;
1892 } else return kFALSE;
1894 case kCentralityMax:
1895 if( SetCentralityMax(value)) {
1896 fCuts[kCentralityMax] = value;
1899 } else return kFALSE;
1901 case kTOFelectronPID:
1902 if( SetTOFElectronPIDCut(value)) {
1903 fCuts[kTOFelectronPID] = value;
1906 } else return kFALSE;
1908 case kdoPhotonAsymmetryCut:
1909 if( SetPhotonAsymmetryCut(value)) {
1910 fCuts[kdoPhotonAsymmetryCut] = value;
1913 } else return kFALSE;
1916 if( SetPsiPairCut(value)) {
1917 fCuts[kPsiPair] = value;
1920 } else return kFALSE;
1923 if( SetCosPAngleCut(value)) {
1924 fCuts[kCosPAngle] = value;
1927 } else return kFALSE;
1931 if( SetSharedElectronCut(value)) {
1932 fCuts[kElecShare] = value;
1935 } else return kFALSE;
1938 if( SetToCloseV0sCut(value)) {
1939 fCuts[kToCloseV0s] = value;
1942 } else return kFALSE;
1945 if( SetRejectExtraSignalsCut(value)) {
1946 fCuts[kExtraSignals] = value;
1949 } else return kFALSE;
1952 if( SetDCARPhotonPrimVtxCut(value)) {
1953 fCuts[kDcaRPrimVtx] = value;
1956 } else return kFALSE;
1959 if( SetDCAZPhotonPrimVtxCut(value)) {
1960 fCuts[kDcaZPrimVtx] = value;
1963 } else return kFALSE;
1967 AliError("Cut id out of range");
1971 AliError("Cut id %d not recognized");
1976 ///________________________________________________________________________
1977 void AliConversionCuts::PrintCuts() {
1978 // Print out current Cut Selection
1979 for(Int_t ic = 0; ic < kNCuts; ic++) {
1980 printf("%-30s : %d \n", fgkCutNames[ic], fCuts[ic]);
1983 ///________________________________________________________________________
1984 Bool_t AliConversionCuts::SetIsHeavyIon(Int_t isHeavyIon)
1992 fDetectorCentrality=0;
1996 fDetectorCentrality=1;
1998 case 3: //allows to select centrality 0-45% in steps of 5% for V0 Multiplicity
2000 fDetectorCentrality=0;
2001 fModCentralityClass=1;
2003 case 4: //allows to select centrality 45-90% in steps of 5% for V0 Multiplicity
2005 fDetectorCentrality=0;
2006 fModCentralityClass=2;
2008 case 5: //strict cut on v0 tracks for MC
2010 fDetectorCentrality=0;
2011 fModCentralityClass=3;
2013 case 6: //allows to select centrality 0-45% in steps of 5% for track mult
2014 //strict cut on v0 tracks for MC
2016 fDetectorCentrality=0;
2017 fModCentralityClass=4;
2019 case 7: //allows to select centrality 45-90% in steps of 5% for V0 Multiplicity
2020 //strict cut on v0 tracks for MC
2022 fDetectorCentrality=0;
2023 fModCentralityClass=5;
2027 fDetectorCentrality=0;
2031 fDetectorCentrality=1;
2034 AliError(Form("SetHeavyIon not defined %d",isHeavyIon));
2040 //___________________________________________________________________
2041 Bool_t AliConversionCuts::SetCentralityMin(Int_t minCentrality)
2044 if(minCentrality<0||minCentrality>9){
2045 AliError(Form("minCentrality not defined %d",minCentrality));
2049 fCentralityMin=minCentrality;
2052 //___________________________________________________________________
2053 Bool_t AliConversionCuts::SetCentralityMax(Int_t maxCentrality)
2056 if(maxCentrality<0||maxCentrality>9){
2057 AliError(Form("maxCentrality not defined %d",maxCentrality));
2060 fCentralityMax=maxCentrality;
2063 ///________________________________________________________________________
2064 Int_t AliConversionCuts::SetSelectSpecialTrigger(Int_t selectSpecialTrigger)
2067 switch(selectSpecialTrigger){
2069 fSpecialTrigger=0; // dont care
2072 fSpecialTrigger=1; // V0AND
2075 fSpecialTrigger=2; // with SDD requested
2078 fSpecialTrigger=3; // V0AND plus with SDD requested
2080 // allows to run MB & 6 other different trigger classes in parallel with the same photon cut
2082 fSpecialTrigger=4; // different trigger class as MB
2083 fTriggerSelectedManually = kTRUE;
2086 fSpecialTrigger=4; // different trigger class as MB
2087 fTriggerSelectedManually = kTRUE;
2090 fSpecialTrigger=4; // different trigger class as MB
2091 fTriggerSelectedManually = kTRUE;
2094 fSpecialTrigger=4; // different trigger class as MB
2095 fTriggerSelectedManually = kTRUE;
2098 fSpecialTrigger=4; // different trigger class as MB
2099 fTriggerSelectedManually = kTRUE;
2102 fSpecialTrigger=4; // different trigger class as MB
2103 fTriggerSelectedManually = kTRUE;
2106 AliError("Warning: Special Trigger Not known");
2111 ///________________________________________________________________________
2112 Bool_t AliConversionCuts::SetMultiplicityMethod(Int_t multiplicityMethod)
2115 fMultiplicityMethod=multiplicityMethod;
2117 // 0 Photon Multiplicity
2118 // 1 TPC Track multiplicity
2124 ///________________________________________________________________________
2125 Bool_t AliConversionCuts::SetRemovePileUp(Int_t removePileUp)
2127 switch(removePileUp){
2129 fRemovePileUp=kFALSE;
2132 fRemovePileUp=kTRUE;
2135 AliError("RemovePileUpCut not defined");
2140 ///________________________________________________________________________
2141 Bool_t AliConversionCuts::SetRejectExtraSignalsCut(Int_t extraSignal) {
2143 switch(extraSignal){
2145 fRejectExtraSignals = 0;
2146 break; // No Rejection
2148 fRejectExtraSignals = 1;
2149 break; // MinBias Header
2151 fRejectExtraSignals = 2;
2152 break; // User String Array
2154 fRejectExtraSignals = 3;
2155 break; // Rejection for Gamma Correction only
2157 AliError(Form("Extra Signal Rejection not defined %d",extraSignal));
2162 ///________________________________________________________________________
2163 Bool_t AliConversionCuts::SetV0Finder(Int_t v0FinderType)
2165 switch (v0FinderType){
2166 case 0: // on fly V0 finder
2167 cout << "have chosen onfly V0" << endl;
2168 fUseOnFlyV0Finder=kTRUE;
2170 case 1: // offline V0 finder
2171 cout << "have chosen offline V0" << endl;
2172 fUseOnFlyV0Finder=kFALSE;
2175 AliError(Form(" v0FinderType not defined %d",v0FinderType));
2180 ///________________________________________________________________________
2181 Bool_t AliConversionCuts::SetEtaCut(Int_t etaCut)
2184 //Set Standard LineCutZValues
2185 fLineCutZValueMin = -2;
2186 fLineCutZValue = 7.;
2191 fLineCutZRSlope = tan(2*atan(exp(-fEtaCut)));
2193 fLineCutZRSlopeMin = 0.;
2195 case 1: // 1.2 // changed from 1.2 to 0.6 on 2013.06.10
2197 fLineCutZRSlope = tan(2*atan(exp(-fEtaCut)));
2199 fLineCutZRSlopeMin = 0.;
2203 fLineCutZRSlope = tan(2*atan(exp(-fEtaCut)));
2205 fLineCutZRSlopeMin = 0.;
2209 fLineCutZRSlope = tan(2*atan(exp(-fEtaCut)));
2211 fLineCutZRSlopeMin = 0.;
2215 fLineCutZRSlope = tan(2*atan(exp(-fEtaCut)));
2217 fLineCutZRSlopeMin = 0.;
2221 fLineCutZRSlope = tan(2*atan(exp(-fEtaCut)));
2223 fLineCutZRSlopeMin = 0.;
2227 fLineCutZRSlope = tan(2*atan(exp(-fEtaCut)));
2229 fLineCutZRSlopeMin = 0.;
2233 fLineCutZRSlope = tan(2*atan(exp(-fEtaCut)));
2235 fLineCutZRSlopeMin = 0.;
2237 // case 8: // 0.1 - 0.8
2239 // fLineCutZRSlope = tan(2*atan(exp(-fEtaCut)));
2240 // fEtaCutMin = 0.1;
2241 // fLineCutZRSlopeMin = tan(2*atan(exp(-fEtaCutMin)));
2245 fLineCutZRSlope = tan(2*atan(exp(-fEtaCut)));
2247 fLineCutZRSlopeMin = 0.;
2251 fLineCutZRSlope = tan(2*atan(exp(-fEtaCut)));
2253 fLineCutZRSlopeMin = 0.;
2256 AliError(Form(" EtaCut not defined %d",etaCut));
2261 ///________________________________________________________________________
2262 Bool_t AliConversionCuts::SetRCut(Int_t RCut){
2285 // High purity cuts for PbPb (remove first layers of material)
2308 AliError("RCut not defined");
2313 ///________________________________________________________________________
2314 Bool_t AliConversionCuts::SetSinglePtCut(Int_t singlePtCut)
2316 switch(singlePtCut){
2317 case 0: // 0.050 GeV
2318 fSinglePtCut = 0.050;
2320 case 1: // 0.100 GeV
2321 fSinglePtCut = 0.100;
2323 case 2: // 0.150 GeV
2324 fSinglePtCut = 0.150;
2326 case 3: // 0.200 GeV
2327 fSinglePtCut = 0.200;
2329 case 4: // 0.075 GeV
2330 fSinglePtCut = 0.075;
2332 case 5: // 0.125 GeV
2333 fSinglePtCut = 0.125;
2336 fSinglePtCut = 0.040;
2342 AliError(Form("singlePtCut not defined %d",singlePtCut));
2347 ///________________________________________________________________________
2348 Bool_t AliConversionCuts::SetTPCClusterCut(Int_t clsTPCCut)
2363 case 4: // 95% of findable clusters
2364 fMinClsTPCToF= 0.95;
2365 fUseCorrectedTPCClsInfo=1;
2367 case 5: // 0% of findable clusters
2369 fUseCorrectedTPCClsInfo=1;
2371 case 6: // 70% of findable clusters
2373 fUseCorrectedTPCClsInfo=1;
2375 case 7: // 0% of findable clusters
2376 fMinClsTPCToF= 0.35;
2377 fUseCorrectedTPCClsInfo=0;
2380 fMinClsTPCToF= 0.35;
2381 fUseCorrectedTPCClsInfo=1;
2385 fUseCorrectedTPCClsInfo=1;
2388 AliError(Form("Warning: clsTPCCut not defined %d",clsTPCCut));
2393 ///________________________________________________________________________
2394 Bool_t AliConversionCuts::SetTPCdEdxCutElectronLine(Int_t ededxSigmaCut)
2396 switch(ededxSigmaCut){
2398 fPIDnSigmaBelowElectronLine=-10;
2399 fPIDnSigmaAboveElectronLine=10;
2402 fPIDnSigmaBelowElectronLine=-5;
2403 fPIDnSigmaAboveElectronLine=5;
2406 fPIDnSigmaBelowElectronLine=-3;
2407 fPIDnSigmaAboveElectronLine=5;
2410 fPIDnSigmaBelowElectronLine=-4;
2411 fPIDnSigmaAboveElectronLine=5;
2414 fPIDnSigmaBelowElectronLine=-6;
2415 fPIDnSigmaAboveElectronLine=7;
2418 fPIDnSigmaBelowElectronLine=-4;
2419 fPIDnSigmaAboveElectronLine=4;
2422 fPIDnSigmaBelowElectronLine=-2.5;
2423 fPIDnSigmaAboveElectronLine=4;
2426 fPIDnSigmaBelowElectronLine=-2;
2427 fPIDnSigmaAboveElectronLine=3.5;
2430 AliError("TPCdEdxCutElectronLine not defined");
2436 ///________________________________________________________________________
2437 Bool_t AliConversionCuts::SetTPCdEdxCutPionLine(Int_t pidedxSigmaCut)
2440 switch(pidedxSigmaCut){
2442 fPIDnSigmaAbovePionLine=-10;
2443 fPIDnSigmaAbovePionLineHighPt=-10;
2446 fPIDnSigmaAbovePionLine=0;
2447 fPIDnSigmaAbovePionLineHighPt=-10;
2450 fPIDnSigmaAbovePionLine=1;
2451 fPIDnSigmaAbovePionLineHighPt=-10;
2454 fPIDnSigmaAbovePionLine=2.5;
2455 fPIDnSigmaAbovePionLineHighPt=-10;
2458 fPIDnSigmaAbovePionLine=0.5;
2459 fPIDnSigmaAbovePionLineHighPt=-10;
2462 fPIDnSigmaAbovePionLine=2.;
2463 fPIDnSigmaAbovePionLineHighPt=-10;
2466 fPIDnSigmaAbovePionLine=2.;
2467 fPIDnSigmaAbovePionLineHighPt=0.5;
2470 fPIDnSigmaAbovePionLine=3.5;
2471 fPIDnSigmaAbovePionLineHighPt=-10;
2474 fPIDnSigmaAbovePionLine=2.;
2475 fPIDnSigmaAbovePionLineHighPt=1.;
2478 fPIDnSigmaAbovePionLine=3.0; // We need a bit less tight cut on dE/dx
2479 fPIDnSigmaAbovePionLineHighPt=-10;
2482 AliError(Form("Warning: pidedxSigmaCut not defined %d",pidedxSigmaCut));
2487 ///________________________________________________________________________
2488 Bool_t AliConversionCuts::SetMinMomPiondEdxCut(Int_t piMomdedxSigmaCut)
2490 switch(piMomdedxSigmaCut){
2492 fPIDMinPnSigmaAbovePionLine=0.5;
2495 fPIDMinPnSigmaAbovePionLine=1.;
2498 fPIDMinPnSigmaAbovePionLine=1.5;
2501 fPIDMinPnSigmaAbovePionLine=20.;
2504 fPIDMinPnSigmaAbovePionLine=50.;
2507 fPIDMinPnSigmaAbovePionLine=0.3;
2510 fPIDMinPnSigmaAbovePionLine=0.25;
2513 fPIDMinPnSigmaAbovePionLine=0.4;
2516 fPIDMinPnSigmaAbovePionLine=0.2;
2519 AliError(Form("piMomdedxSigmaCut not defined %d",piMomdedxSigmaCut));
2524 ///________________________________________________________________________
2525 Bool_t AliConversionCuts::SetMaxMomPiondEdxCut(Int_t piMaxMomdedxSigmaCut)
2527 switch(piMaxMomdedxSigmaCut){
2529 fPIDMaxPnSigmaAbovePionLine=100.;
2532 fPIDMaxPnSigmaAbovePionLine=5.;
2535 fPIDMaxPnSigmaAbovePionLine=4.;
2538 fPIDMaxPnSigmaAbovePionLine=3.5;
2541 fPIDMaxPnSigmaAbovePionLine=3.;
2544 fPIDMaxPnSigmaAbovePionLine=7.;
2547 AliError(Form("piMaxMomdedxSigmaCut not defined %d",piMaxMomdedxSigmaCut));
2552 ///________________________________________________________________________
2553 Bool_t AliConversionCuts::SetLowPRejectionCuts(Int_t LowPRejectionSigmaCut)
2555 switch(LowPRejectionSigmaCut){
2557 fPIDnSigmaAtLowPAroundKaonLine=0;
2558 fPIDnSigmaAtLowPAroundProtonLine=0;
2559 fPIDnSigmaAtLowPAroundPionLine=0;
2560 fDoKaonRejectionLowP = kFALSE;
2561 fDoProtonRejectionLowP = kFALSE;
2562 fDoPionRejectionLowP = kFALSE;
2563 fPIDMinPPionRejectionLowP = fPIDMinPnSigmaAbovePionLine;
2566 fPIDnSigmaAtLowPAroundKaonLine=0.5;
2567 fPIDnSigmaAtLowPAroundProtonLine=0.5;
2568 fPIDnSigmaAtLowPAroundPionLine=0.5;
2569 fDoKaonRejectionLowP = kTRUE;
2570 fDoProtonRejectionLowP = kTRUE;
2571 fDoPionRejectionLowP = kTRUE;
2572 fPIDMinPPionRejectionLowP = fPIDMinPnSigmaAbovePionLine;
2575 fPIDnSigmaAtLowPAroundKaonLine=1;
2576 fPIDnSigmaAtLowPAroundProtonLine=1;
2577 fPIDnSigmaAtLowPAroundPionLine=1;
2578 fDoKaonRejectionLowP = kTRUE;
2579 fDoProtonRejectionLowP = kTRUE;
2580 fDoPionRejectionLowP = kTRUE;
2581 fPIDMinPPionRejectionLowP = fPIDMinPnSigmaAbovePionLine;
2584 fPIDnSigmaAtLowPAroundKaonLine=2.;
2585 fPIDnSigmaAtLowPAroundProtonLine=2.;
2586 fPIDnSigmaAtLowPAroundPionLine=2.;
2587 fDoKaonRejectionLowP = kTRUE;
2588 fDoProtonRejectionLowP = kTRUE;
2589 fDoPionRejectionLowP = kTRUE;
2590 fPIDMinPPionRejectionLowP = fPIDMinPnSigmaAbovePionLine;
2593 fPIDnSigmaAtLowPAroundKaonLine=0.;
2594 fPIDnSigmaAtLowPAroundProtonLine=0.;
2595 fPIDnSigmaAtLowPAroundPionLine=1;
2596 fDoKaonRejectionLowP = kFALSE;
2597 fDoProtonRejectionLowP = kFALSE;
2598 fDoPionRejectionLowP = kTRUE;
2599 fPIDMinPPionRejectionLowP = fPIDMinPnSigmaAbovePionLine;
2602 fPIDnSigmaAtLowPAroundKaonLine=0.;
2603 fPIDnSigmaAtLowPAroundProtonLine=0.;
2604 fPIDnSigmaAtLowPAroundPionLine=1.5;
2605 fDoKaonRejectionLowP = kFALSE;
2606 fDoProtonRejectionLowP = kFALSE;
2607 fDoPionRejectionLowP = kTRUE;
2608 fPIDMinPPionRejectionLowP = fPIDMinPnSigmaAbovePionLine;
2611 fPIDnSigmaAtLowPAroundKaonLine=0.;
2612 fPIDnSigmaAtLowPAroundProtonLine=0.;
2613 fPIDnSigmaAtLowPAroundPionLine=2.;
2614 fDoKaonRejectionLowP = kFALSE;
2615 fDoProtonRejectionLowP = kFALSE;
2616 fDoPionRejectionLowP = kTRUE;
2617 fPIDMinPPionRejectionLowP = fPIDMinPnSigmaAbovePionLine;
2620 fPIDnSigmaAtLowPAroundKaonLine=0.;
2621 fPIDnSigmaAtLowPAroundProtonLine=0.;
2622 fPIDnSigmaAtLowPAroundPionLine=0.5;
2623 fDoKaonRejectionLowP = kFALSE;
2624 fDoProtonRejectionLowP = kFALSE;
2625 fDoPionRejectionLowP = kTRUE;
2626 fPIDMinPPionRejectionLowP = fPIDMinPnSigmaAbovePionLine;
2629 AliError(Form("LowPRejectionSigmaCut not defined %d",LowPRejectionSigmaCut));
2634 ///________________________________________________________________________
2635 Bool_t AliConversionCuts::SetTOFElectronPIDCut(Int_t TOFelectronPID){
2637 switch(TOFelectronPID){
2639 fUseTOFpid = kFALSE;
2640 fTofPIDnSigmaBelowElectronLine=-100;
2641 fTofPIDnSigmaAboveElectronLine=100;
2645 fTofPIDnSigmaBelowElectronLine=-7;
2646 fTofPIDnSigmaAboveElectronLine=7;
2650 fTofPIDnSigmaBelowElectronLine=-5;
2651 fTofPIDnSigmaAboveElectronLine=5;
2655 fTofPIDnSigmaBelowElectronLine=-3;
2656 fTofPIDnSigmaAboveElectronLine=5;
2660 fTofPIDnSigmaBelowElectronLine=-2;
2661 fTofPIDnSigmaAboveElectronLine=3;
2665 fTofPIDnSigmaBelowElectronLine=-3;
2666 fTofPIDnSigmaAboveElectronLine=3;
2669 AliError(Form("TOFElectronCut not defined %d",TOFelectronPID));
2674 ///________________________________________________________________________
2675 Bool_t AliConversionCuts::SetQtMaxCut(Int_t QtMaxCut)
2680 fDoQtGammaSelection=kFALSE;
2681 fDoHighPtQtGammaSelection=kFALSE;
2683 fPtBorderForQt=100.;
2687 fDoHighPtQtGammaSelection=kFALSE;
2689 fPtBorderForQt=100.;
2693 fDoHighPtQtGammaSelection=kFALSE;
2695 fPtBorderForQt=100.;
2699 fDoHighPtQtGammaSelection=kFALSE;
2701 fPtBorderForQt=100.;
2705 fDoHighPtQtGammaSelection=kFALSE;
2707 fPtBorderForQt=100.;
2711 fDoHighPtQtGammaSelection=kFALSE;
2713 fPtBorderForQt=100.;
2717 fDoHighPtQtGammaSelection=kTRUE;
2723 fDoHighPtQtGammaSelection=kFALSE;
2725 fPtBorderForQt=100.;
2728 AliError(Form("Warning: QtMaxCut not defined %d",QtMaxCut));
2733 ///________________________________________________________________________
2734 Bool_t AliConversionCuts::SetChi2GammaCut(Int_t chi2GammaCut)
2737 switch(chi2GammaCut){
2739 fChi2CutConversion = 100.;
2742 fChi2CutConversion = 50.;
2745 fChi2CutConversion = 30.;
2748 fChi2CutConversion = 200.;
2751 fChi2CutConversion = 500.;
2754 fChi2CutConversion = 100000.;
2757 fChi2CutConversion = 5.;
2760 fChi2CutConversion = 10.;
2763 fChi2CutConversion = 20.;
2766 fChi2CutConversion = 15.;
2769 AliError(Form("Warning: Chi2GammaCut not defined %d",chi2GammaCut));
2774 ///________________________________________________________________________
2775 Bool_t AliConversionCuts::SetPsiPairCut(Int_t psiCut) {
2779 fPsiPairCut = 10000; //
2782 fPsiPairCut = 0.1; //
2785 fPsiPairCut = 0.05; // Standard
2788 fPsiPairCut = 0.035; //
2791 fPsiPairCut = 0.15; //
2794 fPsiPairCut = 0.2; //
2797 fPsiPairCut = 0.03; //
2800 fPsiPairCut = 0.025; //
2803 fPsiPairCut = 0.01; //
2806 fPsiPairCut = 0.5; //
2809 AliError(Form("PsiPairCut not defined %d",psiCut));
2815 ///________________________________________________________________________
2816 Bool_t AliConversionCuts::SetPhotonAsymmetryCut(Int_t doPhotonAsymmetryCut){
2818 switch(doPhotonAsymmetryCut){
2820 fDoPhotonAsymmetryCut=0;
2821 fMinPPhotonAsymmetryCut=100.;
2822 fMinPhotonAsymmetry=0.;
2825 fDoPhotonAsymmetryCut=1;
2826 fMinPPhotonAsymmetryCut=3.5;
2827 fMinPhotonAsymmetry=0.04;
2830 fDoPhotonAsymmetryCut=1;
2831 fMinPPhotonAsymmetryCut=3.5;
2832 fMinPhotonAsymmetry=0.06;
2835 AliError(Form("PhotonAsymmetryCut not defined %d",doPhotonAsymmetryCut));
2838 fCuts[kdoPhotonAsymmetryCut]=doPhotonAsymmetryCut;
2841 ///________________________________________________________________________
2842 Bool_t AliConversionCuts::SetCosPAngleCut(Int_t cosCut) {
2846 fCosPAngleCut = TMath::Pi(); // -1
2849 fCosPAngleCut = 0.1; // 0.99500
2852 fCosPAngleCut = 0.05; // 0.99875
2855 fCosPAngleCut = 0.025; // 0.99969
2858 fCosPAngleCut = 0.01; // 0.99995
2861 fCosPAngleCut = 0.2; // 0.98007
2864 fCosPAngleCut = 0.5; // 0.87758
2867 fCosPAngleCut = 0.075; // 0.73169
2870 AliError(Form("Cosine Pointing Angle cut not defined %d",cosCut));
2876 ///________________________________________________________________________
2877 Bool_t AliConversionCuts::SetSharedElectronCut(Int_t sharedElec) {
2881 fDoSharedElecCut = kFALSE;
2884 fDoSharedElecCut = kTRUE;
2887 AliError(Form("Shared Electron Cut not defined %d",sharedElec));
2893 ///________________________________________________________________________
2894 Bool_t AliConversionCuts::SetToCloseV0sCut(Int_t toClose) {
2898 fDoToCloseV0sCut = kFALSE;
2902 fDoToCloseV0sCut = kTRUE;
2906 fDoToCloseV0sCut = kTRUE;
2910 fDoToCloseV0sCut = kTRUE;
2914 AliError(Form("Shared Electron Cut not defined %d",toClose));
2919 ///________________________________________________________________________
2920 Bool_t AliConversionCuts::SetTRDElectronCut(Int_t TRDElectronCut)
2922 switch(TRDElectronCut){
2928 fPIDTRDEfficiency=0.1;
2932 fPIDTRDEfficiency=0.8;
2936 fPIDTRDEfficiency=0.9;
2939 AliError(Form("TRDElectronCut not defined %d",TRDElectronCut));
2946 ///________________________________________________________________________
2947 Bool_t AliConversionCuts::SetDCAZPhotonPrimVtxCut(Int_t DCAZPhotonPrimVtx){
2949 switch(DCAZPhotonPrimVtx){
2951 fDCAZPrimVtxCut = 1000;
2954 fDCAZPrimVtxCut = 10;
2957 fDCAZPrimVtxCut = 5;
2960 fDCAZPrimVtxCut = 4;
2963 fDCAZPrimVtxCut = 3;
2966 fDCAZPrimVtxCut = 2.5;
2969 fDCAZPrimVtxCut = 2;
2972 fDCAZPrimVtxCut = 1.5;
2975 fDCAZPrimVtxCut = 1;
2978 fDCAZPrimVtxCut = 0.5;
2981 cout<<"Warning: DCAZPhotonPrimVtx not defined "<<DCAZPhotonPrimVtx<<endl;
2987 ///________________________________________________________________________
2988 Bool_t AliConversionCuts::SetDCARPhotonPrimVtxCut(Int_t DCARPhotonPrimVtx){
2990 switch(DCARPhotonPrimVtx){
2992 fDCARPrimVtxCut = 1000;
2995 fDCARPrimVtxCut = 10;
2998 fDCARPrimVtxCut = 5;
3001 fDCARPrimVtxCut = 4;
3004 fDCARPrimVtxCut = 3;
3007 fDCARPrimVtxCut = 2.5;
3010 fDCARPrimVtxCut = 2;
3013 fDCARPrimVtxCut = 1.5;
3016 fDCARPrimVtxCut = 1;
3019 fDCARPrimVtxCut = 0.5;
3022 cout<<"Warning: DCARPhotonPrimVtx not defined "<<DCARPhotonPrimVtx<<endl;
3029 //-------------------------------------------------------------
3030 Double_t AliConversionCuts::GetCentrality(AliVEvent *event)
3031 { // Get Event Centrality
3033 AliESDEvent *esdEvent=dynamic_cast<AliESDEvent*>(event);
3035 AliCentrality *fESDCentrality=(AliCentrality*)esdEvent->GetCentrality();
3037 if(fDetectorCentrality==0){
3038 if (fIsHeavyIon==2){
3039 return fESDCentrality->GetCentralityPercentile("V0A"); // default for pPb
3041 return fESDCentrality->GetCentralityPercentile("V0M"); // default
3044 if(fDetectorCentrality==1){
3045 return fESDCentrality->GetCentralityPercentile("CL1");
3049 AliAODEvent *aodEvent=dynamic_cast<AliAODEvent*>(event);
3051 if(aodEvent->GetHeader()){return aodEvent->GetHeader()->GetCentrality();}
3056 //-------------------------------------------------------------
3057 Bool_t AliConversionCuts::IsCentralitySelected(AliVEvent *event, AliVEvent *fMCEvent)
3058 { // Centrality Selection
3059 if(!fIsHeavyIon)return kTRUE;
3061 if(fCentralityMin == fCentralityMax ) return kTRUE;//0-100%
3062 else if(fCentralityMax==0) fCentralityMax=10; //CentralityRange = fCentralityMin-100%
3064 Double_t centrality=GetCentrality(event);
3065 if(centrality<0)return kFALSE;
3067 Int_t centralityC=0;
3068 if (fModCentralityClass == 0){
3069 centralityC= Int_t(centrality/10);
3070 if(centralityC >= fCentralityMin && centralityC < fCentralityMax)
3074 else if (fModCentralityClass ==1){
3075 centralityC= Int_t(centrality);
3076 if(centralityC >= fCentralityMin*5 && centralityC < fCentralityMax*5){
3078 } else return kFALSE;
3080 else if (fModCentralityClass ==2){
3081 centralityC= Int_t(centrality);
3082 if(centralityC >= ((fCentralityMin*5)+45) && centralityC < ((fCentralityMax*5)+45))
3087 Int_t nprimaryTracks = ((AliV0ReaderV1*)AliAnalysisManager::GetAnalysisManager()->GetTask("V0ReaderV1"))->GetNumberOfPrimaryTracks();
3088 Int_t PrimaryTracks10[10][2] =
3101 Int_t PrimaryTracks5a[10][2] =
3114 Int_t PrimaryTracks5b[10][2] =
3129 if(event->IsA()==AliESDEvent::Class()) column = 0;
3130 if(event->IsA()==AliAODEvent::Class()) column = 1;
3132 if (fModCentralityClass == 3){
3134 if(nprimaryTracks > PrimaryTracks10[fCentralityMax][column] && nprimaryTracks <= PrimaryTracks10[fCentralityMin][column])
3139 centralityC= Int_t(centrality/10);
3140 if(centralityC >= fCentralityMin && centralityC < fCentralityMax)
3145 else if (fModCentralityClass ==4){
3147 if(nprimaryTracks > PrimaryTracks5a[fCentralityMax][column] && nprimaryTracks <= PrimaryTracks5a[fCentralityMin][column])
3152 centralityC= Int_t(centrality);
3153 if(centralityC >= fCentralityMin*5 && centralityC < fCentralityMax*5){
3155 } else return kFALSE;
3158 else if (fModCentralityClass ==5){
3160 if(nprimaryTracks > PrimaryTracks5b[fCentralityMax][column] && nprimaryTracks <= PrimaryTracks5b[fCentralityMin][column])
3165 centralityC= Int_t(centrality);
3166 if(centralityC >= ((fCentralityMin*5)+45) && centralityC < ((fCentralityMax*5)+45))
3174 ///________________________________________________________________________
3175 Bool_t AliConversionCuts::VertexZCut(AliVEvent *event){
3176 // Cut on z position of primary vertex
3177 Double_t fVertexZ=event->GetPrimaryVertex()->GetZ();
3179 if(abs(fVertexZ)>fMaxVertexZ)return kFALSE;
3181 if (fIsHeavyIon == 2){
3182 if(fUtils->IsFirstEventInChunk(event)) return kFALSE;
3183 if(!fUtils->IsVertexSelected2013pA(event)) return kFALSE;
3189 ///________________________________________________________________________
3191 Int_t AliConversionCuts::GetNumberOfContributorsVtx(AliVEvent *event){
3192 // returns number of contributors to the vertex
3194 AliESDEvent *fESDEvent=dynamic_cast<AliESDEvent*>(event);
3196 if (fESDEvent->GetPrimaryVertex() != NULL){
3197 if(fESDEvent->GetPrimaryVertex()->GetNContributors()>0) {
3198 // cout << "accepted global" << fESDEvent->GetEventNumberInFile() << " with NCont: " << fESDEvent->GetPrimaryVertex()->GetNContributors() << endl;
3199 return fESDEvent->GetPrimaryVertex()->GetNContributors();
3203 if(fESDEvent->GetPrimaryVertexSPD() !=NULL){
3204 if(fESDEvent->GetPrimaryVertexSPD()->GetNContributors()>0) {
3205 // cout << "accepted SPD" << fESDEvent->GetEventNumberInFile() << " with NCont: " << fESDEvent->GetPrimaryVertexSPD()->GetNContributors() << endl;
3206 return fESDEvent->GetPrimaryVertexSPD()->GetNContributors();
3208 AliWarning(Form("Number of contributors from bad vertex type:: %s",fESDEvent->GetPrimaryVertex()->GetName()));
3209 // cout << "rejected " << fESDEvent->GetEventNumberInFile() << endl;
3215 AliAODEvent *fAODEvent=dynamic_cast<AliAODEvent*>(event);
3217 if (fAODEvent->GetPrimaryVertex() != NULL){
3218 if(fAODEvent->GetPrimaryVertex()->GetNContributors()>0) {
3219 return fAODEvent->GetPrimaryVertex()->GetNContributors();
3222 if(fAODEvent->GetPrimaryVertexSPD() !=NULL){
3223 if(fAODEvent->GetPrimaryVertexSPD()->GetNContributors()>0) {
3224 return fAODEvent->GetPrimaryVertexSPD()->GetNContributors();
3226 AliWarning(Form("Number of contributors from bad vertex type:: %s",fAODEvent->GetPrimaryVertex()->GetName()));
3231 // cout << "rejected " << fESDEvent->GetEventNumberInFile() << endl;
3235 ///________________________________________________________________________
3237 Bool_t AliConversionCuts::IsTriggerSelected(AliVEvent *fInputEvent)
3240 AliInputEventHandler *fInputHandler=(AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
3242 UInt_t isSelected = AliVEvent::kAny;
3243 if (fInputHandler==NULL) return kFALSE;
3244 if( fInputHandler->GetEventSelection() || fInputEvent->IsA()==AliAODEvent::Class()) {
3245 if (!fTriggerSelectedManually){
3246 if (fPreSelCut) fOfflineTriggerMask = AliVEvent::kAny;
3248 if (fIsHeavyIon == 1) fOfflineTriggerMask = AliVEvent::kMB | AliVEvent::kCentral | AliVEvent::kSemiCentral;
3249 else if (fIsHeavyIon == 2) fOfflineTriggerMask = AliVEvent::kINT7;
3250 else fOfflineTriggerMask = AliVEvent::kMB;
3253 // Get the actual offline trigger mask for the event and AND it with the
3254 // requested mask. If no mask requested select by default the event.
3255 // if (fPreSelCut) cout << "Trigger selected from outside: "<< fTriggerSelectedManually <<"\t Offline Trigger mask for Precut: " << fOfflineTriggerMask << endl;
3256 // else cout << "Trigger selected from outside: "<< fTriggerSelectedManually <<"\t Offline Trigger mask: " << fOfflineTriggerMask << endl;
3258 if (fOfflineTriggerMask)
3259 isSelected = fOfflineTriggerMask & fInputHandler->IsEventSelected();
3261 fIsSDDFired = !(fInputHandler->IsEventSelected() & AliVEvent::kFastOnly);
3265 if (fIsSDDFired) hTriggerClass->Fill(33);
3266 if (fInputHandler->IsEventSelected() & AliVEvent::kMB)hTriggerClass->Fill(0);
3267 if (fInputHandler->IsEventSelected() & AliVEvent::kINT7)hTriggerClass->Fill(1);
3268 if (fInputHandler->IsEventSelected() & AliVEvent::kMUON)hTriggerClass->Fill(2);
3269 if (fInputHandler->IsEventSelected() & AliVEvent::kHighMult)hTriggerClass->Fill(3);
3270 if (fInputHandler->IsEventSelected() & AliVEvent::kEMC1)hTriggerClass->Fill(4);
3271 if (fInputHandler->IsEventSelected() & AliVEvent::kCINT5)hTriggerClass->Fill(5);
3272 if (fInputHandler->IsEventSelected() & AliVEvent::kCMUS5)hTriggerClass->Fill(6);
3273 if (fInputHandler->IsEventSelected() & AliVEvent::kMUSPB)hTriggerClass->Fill(6);
3274 if (fInputHandler->IsEventSelected() & AliVEvent::kMUSH7)hTriggerClass->Fill(7);
3275 if (fInputHandler->IsEventSelected() & AliVEvent::kMUSHPB)hTriggerClass->Fill(7);
3276 if (fInputHandler->IsEventSelected() & AliVEvent::kMUL7)hTriggerClass->Fill(8);
3277 if (fInputHandler->IsEventSelected() & AliVEvent::kMuonLikePB)hTriggerClass->Fill(8);
3278 if (fInputHandler->IsEventSelected() & AliVEvent::kMUU7)hTriggerClass->Fill(9);
3279 if (fInputHandler->IsEventSelected() & AliVEvent::kMuonUnlikePB)hTriggerClass->Fill(9);
3280 if (fInputHandler->IsEventSelected() & AliVEvent::kEMC7)hTriggerClass->Fill(10);
3281 if (fInputHandler->IsEventSelected() & AliVEvent::kEMC8)hTriggerClass->Fill(10);
3282 if (fInputHandler->IsEventSelected() & AliVEvent::kMUS7)hTriggerClass->Fill(11);
3283 if (fInputHandler->IsEventSelected() & AliVEvent::kPHI1)hTriggerClass->Fill(12);
3284 if (fInputHandler->IsEventSelected() & AliVEvent::kPHI7)hTriggerClass->Fill(13);
3285 if (fInputHandler->IsEventSelected() & AliVEvent::kPHI8)hTriggerClass->Fill(13);
3286 if (fInputHandler->IsEventSelected() & AliVEvent::kPHOSPb)hTriggerClass->Fill(13);
3287 if (fInputHandler->IsEventSelected() & AliVEvent::kEMCEJE)hTriggerClass->Fill(14);
3288 if (fInputHandler->IsEventSelected() & AliVEvent::kEMCEGA)hTriggerClass->Fill(15);
3289 if (fInputHandler->IsEventSelected() & AliVEvent::kCentral)hTriggerClass->Fill(16);
3290 if (fInputHandler->IsEventSelected() & AliVEvent::kSemiCentral)hTriggerClass->Fill(17);
3291 if (fInputHandler->IsEventSelected() & AliVEvent::kDG5)hTriggerClass->Fill(18);
3292 if (fInputHandler->IsEventSelected() & AliVEvent::kZED)hTriggerClass->Fill(19);
3293 if (fInputHandler->IsEventSelected() & AliVEvent::kSPI7)hTriggerClass->Fill(20);
3294 if (fInputHandler->IsEventSelected() & AliVEvent::kSPI)hTriggerClass->Fill(20);
3295 if (fInputHandler->IsEventSelected() & AliVEvent::kINT8)hTriggerClass->Fill(21);
3296 if (fInputHandler->IsEventSelected() & AliVEvent::kMuonSingleLowPt8)hTriggerClass->Fill(22);
3297 if (fInputHandler->IsEventSelected() & AliVEvent::kMuonSingleHighPt8)hTriggerClass->Fill(23);
3298 if (fInputHandler->IsEventSelected() & AliVEvent::kMuonLikeLowPt8)hTriggerClass->Fill(24);
3299 if (fInputHandler->IsEventSelected() & AliVEvent::kMuonUnlikeLowPt8)hTriggerClass->Fill(25);
3300 if (fInputHandler->IsEventSelected() & AliVEvent::kMuonUnlikeLowPt0)hTriggerClass->Fill(26);
3301 if (fInputHandler->IsEventSelected() & AliVEvent::kUserDefined)hTriggerClass->Fill(27);
3302 if (fInputHandler->IsEventSelected() & AliVEvent::kTRD)hTriggerClass->Fill(28);
3303 if (fInputHandler->IsEventSelected() & AliVEvent::kFastOnly)hTriggerClass->Fill(29);
3304 if (fInputHandler->IsEventSelected() & AliVEvent::kAnyINT)hTriggerClass->Fill(30);
3305 if (fInputHandler->IsEventSelected() & AliVEvent::kAny)hTriggerClass->Fill(31);
3306 if (!fInputHandler->IsEventSelected()) hTriggerClass->Fill(34);
3309 if(hTriggerClassSelected && isSelected){
3310 if (!fIsSDDFired) hTriggerClassSelected->Fill(33);
3311 if (fInputHandler->IsEventSelected() & AliVEvent::kMB)hTriggerClassSelected->Fill(0);
3312 if (fInputHandler->IsEventSelected() & AliVEvent::kINT7)hTriggerClassSelected->Fill(1);
3313 if (fInputHandler->IsEventSelected() & AliVEvent::kMUON)hTriggerClassSelected->Fill(2);
3314 if (fInputHandler->IsEventSelected() & AliVEvent::kHighMult)hTriggerClassSelected->Fill(3);
3315 if (fInputHandler->IsEventSelected() & AliVEvent::kEMC1)hTriggerClassSelected->Fill(4);
3316 if (fInputHandler->IsEventSelected() & AliVEvent::kCINT5)hTriggerClassSelected->Fill(5);
3317 if (fInputHandler->IsEventSelected() & AliVEvent::kCMUS5)hTriggerClassSelected->Fill(6);
3318 if (fInputHandler->IsEventSelected() & AliVEvent::kMUSPB)hTriggerClassSelected->Fill(6);
3319 if (fInputHandler->IsEventSelected() & AliVEvent::kMUSH7)hTriggerClassSelected->Fill(7);
3320 if (fInputHandler->IsEventSelected() & AliVEvent::kMUSHPB)hTriggerClassSelected->Fill(7);
3321 if (fInputHandler->IsEventSelected() & AliVEvent::kMUL7)hTriggerClassSelected->Fill(8);
3322 if (fInputHandler->IsEventSelected() & AliVEvent::kMuonLikePB)hTriggerClassSelected->Fill(8);
3323 if (fInputHandler->IsEventSelected() & AliVEvent::kMUU7)hTriggerClassSelected->Fill(9);
3324 if (fInputHandler->IsEventSelected() & AliVEvent::kMuonUnlikePB)hTriggerClassSelected->Fill(9);
3325 if (fInputHandler->IsEventSelected() & AliVEvent::kEMC7)hTriggerClassSelected->Fill(10);
3326 if (fInputHandler->IsEventSelected() & AliVEvent::kEMC8)hTriggerClassSelected->Fill(10);
3327 if (fInputHandler->IsEventSelected() & AliVEvent::kMUS7)hTriggerClassSelected->Fill(11);
3328 if (fInputHandler->IsEventSelected() & AliVEvent::kPHI1)hTriggerClassSelected->Fill(12);
3329 if (fInputHandler->IsEventSelected() & AliVEvent::kPHI7)hTriggerClassSelected->Fill(13);
3330 if (fInputHandler->IsEventSelected() & AliVEvent::kPHI8)hTriggerClassSelected->Fill(13);
3331 if (fInputHandler->IsEventSelected() & AliVEvent::kPHOSPb)hTriggerClassSelected->Fill(13);
3332 if (fInputHandler->IsEventSelected() & AliVEvent::kEMCEJE)hTriggerClassSelected->Fill(14);
3333 if (fInputHandler->IsEventSelected() & AliVEvent::kEMCEGA)hTriggerClassSelected->Fill(15);
3334 if (fInputHandler->IsEventSelected() & AliVEvent::kCentral)hTriggerClassSelected->Fill(16);
3335 if (fInputHandler->IsEventSelected() & AliVEvent::kSemiCentral)hTriggerClassSelected->Fill(17);
3336 if (fInputHandler->IsEventSelected() & AliVEvent::kDG5)hTriggerClassSelected->Fill(18);
3337 if (fInputHandler->IsEventSelected() & AliVEvent::kZED)hTriggerClassSelected->Fill(19);
3338 if (fInputHandler->IsEventSelected() & AliVEvent::kSPI7)hTriggerClassSelected->Fill(20);
3339 if (fInputHandler->IsEventSelected() & AliVEvent::kSPI)hTriggerClassSelected->Fill(20);
3340 if (fInputHandler->IsEventSelected() & AliVEvent::kINT8)hTriggerClassSelected->Fill(21);
3341 if (fInputHandler->IsEventSelected() & AliVEvent::kMuonSingleLowPt8)hTriggerClassSelected->Fill(22);
3342 if (fInputHandler->IsEventSelected() & AliVEvent::kMuonSingleHighPt8)hTriggerClassSelected->Fill(23);
3343 if (fInputHandler->IsEventSelected() & AliVEvent::kMuonLikeLowPt8)hTriggerClassSelected->Fill(24);
3344 if (fInputHandler->IsEventSelected() & AliVEvent::kMuonUnlikeLowPt8)hTriggerClassSelected->Fill(25);
3345 if (fInputHandler->IsEventSelected() & AliVEvent::kMuonUnlikeLowPt0)hTriggerClassSelected->Fill(26);
3346 if (fInputHandler->IsEventSelected() & AliVEvent::kUserDefined)hTriggerClassSelected->Fill(27);
3347 if (fInputHandler->IsEventSelected() & AliVEvent::kTRD)hTriggerClassSelected->Fill(28);
3348 if (fInputHandler->IsEventSelected() & AliVEvent::kFastOnly)hTriggerClassSelected->Fill(29);
3349 if (fInputHandler->IsEventSelected() & AliVEvent::kAnyINT)hTriggerClassSelected->Fill(30);
3350 if (fInputHandler->IsEventSelected() & AliVEvent::kAny)hTriggerClassSelected->Fill(31);
3353 if(!isSelected)return kFALSE;
3359 ///________________________________________________________________________
3360 Int_t AliConversionCuts::GetFirstTPCRow(Double_t radius){
3361 // Get first TPC row
3362 Int_t firstTPCRow = 0;
3363 Double_t radiusI = 84.8;
3364 Double_t radiusO = 134.6;
3365 Double_t radiusOB = 198.;
3366 Double_t rSizeI = 0.75;
3367 Double_t rSizeO = 1.;
3368 Double_t rSizeOB = 1.5;
3372 if(radius <= radiusI){
3375 if(radius>radiusI && radius<=radiusO){
3376 firstTPCRow = (Int_t)((radius-radiusI)/rSizeI);
3378 if(radius>radiusO && radius<=radiusOB){
3379 firstTPCRow = (Int_t)(nClsI+(radius-radiusO)/rSizeO);
3382 if(radius>radiusOB){
3383 firstTPCRow =(Int_t)(nClsIO+(radius-radiusOB)/rSizeOB);
3389 Bool_t AliConversionCuts::CosinePAngleCut(const AliConversionPhotonBase * photon, AliVEvent * event) const {
3390 ///Check if passes cosine of pointing angle cut
3391 if(GetCosineOfPointingAngle(photon, event) < (TMath::Cos(fCosPAngleCut))){
3397 Double_t AliConversionCuts::GetCosineOfPointingAngle( const AliConversionPhotonBase * photon, AliVEvent * event) const{
3398 // calculates the pointing angle of the recalculated V0
3400 Double_t momV0[3] = {0,0,0};
3401 if(event->IsA()==AliESDEvent::Class()){
3402 AliESDEvent *esdEvent = dynamic_cast<AliESDEvent*>(event);
3403 if(!esdEvent) return -999;
3404 AliESDv0 *v0 = esdEvent->GetV0(photon->GetV0Index());
3405 if(!v0) return -999;
3406 v0->GetPxPyPz(momV0[0],momV0[1],momV0[2]);
3408 if(event->IsA()==AliAODEvent::Class()){
3409 momV0[0] = photon->GetPx();
3410 momV0[1] = photon->GetPy();
3411 momV0[2] = photon->GetPz();
3414 //Double_t momV0[3] = { photon->GetPx(), photon->GetPy(), photon->GetPz() }; //momentum of the V0
3415 Double_t PosV0[3] = { photon->GetConversionX() - event->GetPrimaryVertex()->GetX(),
3416 photon->GetConversionY() - event->GetPrimaryVertex()->GetY(),
3417 photon->GetConversionZ() - event->GetPrimaryVertex()->GetZ() }; //Recalculated V0 Position vector
3419 Double_t momV02 = momV0[0]*momV0[0] + momV0[1]*momV0[1] + momV0[2]*momV0[2];
3420 Double_t PosV02 = PosV0[0]*PosV0[0] + PosV0[1]*PosV0[1] + PosV0[2]*PosV0[2];
3423 Double_t cosinePointingAngle = -999;
3424 if(momV02*PosV02 > 0.0)
3425 cosinePointingAngle = (PosV0[0]*momV0[0] + PosV0[1]*momV0[1] + PosV0[2]*momV0[2] ) / TMath::Sqrt(momV02 * PosV02);
3427 return cosinePointingAngle;
3430 ///________________________________________________________________________
3431 Bool_t AliConversionCuts::PsiPairCut(const AliConversionPhotonBase * photon) const {
3433 if(photon->GetPsiPair() > fPsiPairCut){
3439 ///________________________________________________________________________
3440 TString AliConversionCuts::GetCutNumber(){
3441 // returns TString with current cut number
3443 for(Int_t ii=0;ii<kNCuts;ii++){
3444 a.Append(Form("%d",fCuts[ii]));
3449 ///________________________________________________________________________
3450 void AliConversionCuts::FillElectonLabelArray(AliAODConversionPhoton* photon, Int_t nV0){
3452 Int_t posLabel = photon->GetTrackLabelPositive();
3453 Int_t negLabel = photon->GetTrackLabelNegative();
3455 fElectronLabelArray[nV0*2] = posLabel;
3456 fElectronLabelArray[(nV0*2)+1] = negLabel;
3458 ///________________________________________________________________________
3459 Bool_t AliConversionCuts::RejectSharedElectronV0s(AliAODConversionPhoton* photon, Int_t nV0, Int_t nV0s){
3461 Int_t posLabel = photon->GetTrackLabelPositive();
3462 Int_t negLabel = photon->GetTrackLabelNegative();
3464 for(Int_t i = 0; i<nV0s*2;i++){
3465 if(i==nV0*2) continue;
3466 if(i==(nV0*2)+1) continue;
3467 if(fElectronLabelArray[i] == posLabel){
3469 if(fElectronLabelArray[i] == negLabel){
3475 ///________________________________________________________________________
3476 Bool_t AliConversionCuts::RejectToCloseV0s(AliAODConversionPhoton* photon, TList *photons, Int_t nV0){
3479 Double_t posX = photon->GetConversionX();
3480 Double_t posY = photon->GetConversionY();
3481 Double_t posZ = photon->GetConversionZ();
3483 for(Int_t i = 0;i<photons->GetEntries();i++){
3484 if(nV0 == i) continue;
3485 AliAODConversionPhoton *photonComp = (AliAODConversionPhoton*) photons->At(i);
3486 Double_t posCompX = photonComp->GetConversionX();
3487 Double_t posCompY = photonComp->GetConversionY();
3488 Double_t posCompZ = photonComp->GetConversionZ();
3490 Double_t dist = pow((posX - posCompX),2)+pow((posY - posCompY),2)+pow((posZ - posCompZ),2);
3492 if(dist < fminV0Dist*fminV0Dist){
3493 if(photon->GetChi2perNDF() < photonComp->GetChi2perNDF()) return kTRUE;
3501 ///________________________________________________________________________
3502 void AliConversionCuts::GetNotRejectedParticles(Int_t rejection, TList *HeaderList, AliVEvent *MCEvent){
3506 if(fNotRejectedStart){
3507 delete[] fNotRejectedStart;
3508 fNotRejectedStart = NULL;
3510 if(fNotRejectedEnd){
3511 delete[] fNotRejectedEnd;
3512 fNotRejectedEnd = NULL;
3514 if(fGeneratorNames){
3515 delete[] fGeneratorNames;
3516 fGeneratorNames = NULL;
3519 if(rejection == 0) return; // No Rejection
3521 AliGenCocktailEventHeader *cHeader = 0x0;
3522 AliAODMCHeader *cHeaderAOD = 0x0;
3523 Bool_t headerFound = kFALSE;
3525 if(MCEvent->IsA()==AliMCEvent::Class()){
3526 cHeader = dynamic_cast<AliGenCocktailEventHeader*>(dynamic_cast<AliMCEvent*>(MCEvent)->GenEventHeader());
3527 if(cHeader) headerFound = kTRUE;
3529 if(MCEvent->IsA()==AliAODEvent::Class()){ // MCEvent is a AODEvent in case of AOD
3530 cHeaderAOD = dynamic_cast<AliAODMCHeader*>(MCEvent->FindListObject(AliAODMCHeader::StdBranchName()));
3531 if(cHeaderAOD) headerFound = kTRUE;
3535 TList *genHeaders = 0x0;
3536 if(cHeader) genHeaders = cHeader->GetHeaders();
3538 genHeaders = cHeaderAOD->GetCocktailHeaders();
3539 if(genHeaders->GetEntries()==1){
3540 SetRejectExtraSignalsCut(0);
3544 AliGenEventHeader* gh = 0;
3546 if(rejection == 1 || rejection == 3) fnHeaders = 1; // MinBiasHeader
3547 if(rejection == 2){ // TList of Headers Names
3548 for(Int_t i = 0; i<genHeaders->GetEntries();i++){
3549 gh = (AliGenEventHeader*)genHeaders->At(i);
3550 TString GeneratorName = gh->GetName();
3551 for(Int_t j = 0; j<HeaderList->GetEntries();j++){
3552 TString GeneratorInList = ((TObjString*)HeaderList->At(j))->GetString();
3553 if(GeneratorName.CompareTo(GeneratorInList) == 0){
3561 fNotRejectedStart = new Int_t[fnHeaders];
3562 fNotRejectedEnd = new Int_t[fnHeaders];
3563 fGeneratorNames = new TString[fnHeaders];
3565 if(rejection == 1 || rejection == 3){
3566 fNotRejectedStart[0] = 0;
3567 fNotRejectedEnd[0] = ((AliGenEventHeader*)genHeaders->At(0))->NProduced()-1;
3568 fGeneratorNames[0] = ((AliGenEventHeader*)genHeaders->At(0))->GetName();
3572 Int_t firstindex = 0;
3573 Int_t lastindex = -1;
3575 for(Int_t i = 0; i<genHeaders->GetEntries();i++){
3576 gh = (AliGenEventHeader*)genHeaders->At(i);
3577 TString GeneratorName = gh->GetName();
3578 lastindex = lastindex + gh->NProduced();
3579 for(Int_t j = 0; j<HeaderList->GetEntries();j++){
3580 TString GeneratorInList = ((TObjString*)HeaderList->At(j))->GetString();
3581 if(GeneratorName.CompareTo(GeneratorInList) == 0){
3582 fNotRejectedStart[nummer] = firstindex;
3583 fNotRejectedEnd[nummer] = lastindex;
3584 fGeneratorNames[nummer] = GeneratorName;
3585 //cout << "Number of particles produced for: " << i << "\t" << GeneratorName.Data() << "\t" << lastindex-firstindex+1 << endl;
3590 firstindex = firstindex + gh->NProduced();
3592 } else { // No Cocktail Header Found
3593 fNotRejectedStart = new Int_t[1];
3594 fNotRejectedEnd = new Int_t[1];
3597 fNotRejectedStart[0] = 0;
3598 fNotRejectedEnd[0] = static_cast<AliMCEvent*>(MCEvent)->Stack()->GetNprimary()-1;
3599 fGeneratorNames = new TString[1];
3600 fGeneratorNames[0] = "NoCocktailGeneratorFound";
3602 AliGenPythiaEventHeader *mcHeaderPythia = dynamic_cast<AliGenPythiaEventHeader*>(static_cast<AliMCEvent*>(MCEvent)->GenEventHeader());
3603 if (mcHeaderPythia) fGeneratorNames[0] = "NoCocktailGeneratorFound_Pythia";
3604 AliGenDPMjetEventHeader *mcHeaderPhojet = dynamic_cast<AliGenDPMjetEventHeader*>(static_cast<AliMCEvent*>(MCEvent)->GenEventHeader());
3605 if (mcHeaderPhojet) fGeneratorNames[0] = "NoCocktailGeneratorFound_Phojet";
3606 AliGenHijingEventHeader *mcHeaderHijing = dynamic_cast<AliGenHijingEventHeader*>(static_cast<AliMCEvent*>(MCEvent)->GenEventHeader());
3607 if (mcHeaderHijing) fGeneratorNames[0] = "NoCocktailGeneratorFound_Hijing";
3609 SetRejectExtraSignalsCut(0);
3613 //_________________________________________________________________________
3614 Int_t AliConversionCuts::IsParticleFromBGEvent(Int_t index, AliStack *MCStack, AliVEvent *InputEvent){
3616 // Not Accepted == kFALSE == 0
3617 // Accepted == kTRUE == 1
3618 // FirstHeader == kTRUE == 3
3619 if(index < 0) return 0; // No Particle
3622 if(!InputEvent || InputEvent->IsA()==AliESDEvent::Class()){
3623 if( index >= MCStack->GetNprimary()){ // Secondary Particle
3624 if( ((TParticle*)MCStack->Particle(index))->GetMother(0) < 0) return 1; // Secondary Particle without Mother??
3625 return IsParticleFromBGEvent(((TParticle*)MCStack->Particle(index))->GetMother(0),MCStack,InputEvent);
3627 for(Int_t i = 0;i<fnHeaders;i++){
3628 if(index >= fNotRejectedStart[i] && index <= fNotRejectedEnd[i]){
3630 if(i == 0) accepted = 2; // MB Header
3634 else if(InputEvent->IsA()==AliAODEvent::Class()){
3635 TClonesArray *AODMCTrackArray = dynamic_cast<TClonesArray*>(InputEvent->FindListObject(AliAODMCParticle::StdBranchName()));
3636 AliAODMCParticle *aodMCParticle = static_cast<AliAODMCParticle*>(AODMCTrackArray->At(index));
3637 if(!aodMCParticle) return 1; // Photon Without a Mother ? --> Accepted
3638 if(!aodMCParticle->IsPrimary()){
3639 if( aodMCParticle->GetMother() < 0) return 1;// Secondary Particle without Mother??
3640 return IsParticleFromBGEvent(aodMCParticle->GetMother(),MCStack,InputEvent);
3642 index = abs(static_cast<AliAODMCParticle*>(AODMCTrackArray->At(index))->GetLabel());
3643 for(Int_t i = 0;i<fnHeaders;i++){
3644 if(index >= fNotRejectedStart[i] && index <= fNotRejectedEnd[i]){
3646 if(i == 0) accepted = 2; // MB Header
3653 //_________________________________________________________________________
3654 Int_t AliConversionCuts::IsEventAcceptedByConversionCut(AliConversionCuts *ReaderCuts, AliVEvent *InputEvent, AliMCEvent *MCEvent, Bool_t isHeavyIon){
3656 if ( !IsTriggerSelected(InputEvent) )
3659 if(isHeavyIon && !(IsCentralitySelected(InputEvent,MCEvent)))
3660 return 1; // Check Centrality --> Not Accepted => eventQuality = 1
3662 if(!isHeavyIon && GetIsFromPileup()){
3663 if(InputEvent->IsPileupFromSPD(3,0.8,3.,2.,5.)){
3665 return 6; // Check Pileup --> Not Accepted => eventQuality = 6
3669 Bool_t hasV0And = ReaderCuts->HasV0AND();
3670 Bool_t isSDDFired = ReaderCuts->IsSDDFired();
3671 if( (IsSpecialTrigger() == 2 || IsSpecialTrigger() == 3) && !isSDDFired && !MCEvent)
3672 return 7; // With SDD requested but no fired
3674 if( (IsSpecialTrigger() == 1 || IsSpecialTrigger() == 3) && !hasV0And)
3675 return 8; // V0AND requested but no fired
3680 //_________________________________________________________________________
3681 Float_t AliConversionCuts::GetWeightForMeson(TString period, Int_t index, AliStack *MCStack, AliVEvent *InputEvent){
3682 if (!(period.CompareTo("LHC12f1a") == 0 || period.CompareTo("LHC12f1b") == 0 || period.CompareTo("LHC12i3") == 0 || period.CompareTo("LHC11a10a") == 0 || period.CompareTo("LHC11a10b") == 0 || period.CompareTo("LHC11a10b_bis") == 0 || period.CompareTo("LHC11a10a_bis") == 0 || period.CompareTo("LHC11a10b_plus") == 0 || period.CompareTo("LHC13d2") == 0)) return 1.;
3685 for (Int_t i = 0; i < fnHeaders; i++){
3686 if (index >= fNotRejectedStart[i] && index < fNotRejectedEnd[i]+1){
3687 if (fGeneratorNames[i].CompareTo("Pythia") == 0){
3689 } else if (fGeneratorNames[i].CompareTo("DPMJET") == 0){
3691 } else if (fGeneratorNames[i].CompareTo("HIJING") == 0 ||
3692 fGeneratorNames[i].CompareTo("Hijing") == 0 ||
3693 fGeneratorNames[i].Contains("hijing")){
3695 } else if (fGeneratorNames[i].CompareTo("BOX") == 0){
3697 } else if (fGeneratorNames[i].CompareTo("PARAM") == 0){
3699 } else if (fGeneratorNames[i].CompareTo("NoCocktailGeneratorFound") == 0){
3701 } else if (fGeneratorNames[i].CompareTo("NoCocktailGeneratorFound_Pythia") == 0){
3703 } else if (fGeneratorNames[i].CompareTo("NoCocktailGeneratorFound_Phojet") == 0){
3705 } else if (fGeneratorNames[i].CompareTo("NoCocktailGeneratorFound_Hijing") == 0){
3708 TString periodName = ((AliV0ReaderV1*)AliAnalysisManager::GetAnalysisManager()
3709 ->GetTask("V0ReaderV1"))->GetPeriodName();
3711 if (periodName.Contains("LHC13d2")){
3716 if (kCaseGen == 0) return 1;
3719 Double_t mesonPt = 0;
3720 Double_t mesonMass = 0;
3722 if(!InputEvent || InputEvent->IsA()==AliESDEvent::Class()){
3723 mesonPt = ((TParticle*)MCStack->Particle(index))->Pt();
3724 mesonMass = ((TParticle*)MCStack->Particle(index))->GetCalcMass();
3725 PDGCode = ((TParticle*)MCStack->Particle(index))->GetPdgCode();
3727 else if(InputEvent->IsA()==AliAODEvent::Class()){
3728 TClonesArray *AODMCTrackArray = dynamic_cast<TClonesArray*>(InputEvent->FindListObject(AliAODMCParticle::StdBranchName()));
3729 AliAODMCParticle *aodMCParticle = static_cast<AliAODMCParticle*>(AODMCTrackArray->At(index));
3730 mesonPt = aodMCParticle->Pt();
3731 mesonMass = aodMCParticle->GetCalcMass();
3732 PDGCode = aodMCParticle->GetPdgCode();
3735 Float_t functionResultMC = 1.;
3736 if (kCaseGen == 1){ // Pythia 6
3737 Float_t dNdyMC = 2.1462;
3738 Float_t nMC = 7.06055;
3739 Float_t tMC = 0.12533;
3740 if ( PDGCode == 111){
3744 } else if ( PDGCode == 221){
3749 functionResultMC = dNdyMC / ( 2 * TMath::Pi())*(nMC-1.)*(nMC-2.) / (nMC*tMC*(nMC*tMC+mesonMass*(nMC-2.))) * TMath::Power(1.+(TMath::Sqrt(mesonPt*mesonPt+mesonMass*mesonMass)-mesonMass)/(nMC*tMC), -nMC);
3750 } else if (kCaseGen == 2){ // Phojet
3751 Float_t dNdyMC = 2.35978;
3752 Float_t nMC = 6.81795;
3753 Float_t tMC = 0.11492;
3754 if ( PDGCode == 111){
3758 } else if ( PDGCode == 221){
3763 functionResultMC = dNdyMC / ( 2 * TMath::Pi())*(nMC-1.)*(nMC-2.) / (nMC*tMC*(nMC*tMC+mesonMass*(nMC-2.))) * TMath::Power(1.+(TMath::Sqrt(mesonPt*mesonPt+mesonMass*mesonMass)-mesonMass)/(nMC*tMC), -nMC);
3764 } else if (kCaseGen == 4){ // BOX generators pp
3765 // functionResultMC = 1./sqrt(1.-mesonMass*mesonMass/((mesonMass*mesonMass+mesonPt*mesonPt)*cosh(mesonY)*cosh(mesonY)));
3766 Float_t a = 0.23437;
3768 Float_t c = -1430.5863;
3769 Float_t d = -0.6966624;
3770 Float_t e = 252.3742;
3771 if ( PDGCode == 111){
3777 } else if ( PDGCode == 221){
3784 functionResultMC = a*TMath::Power(mesonPt,-1.*(b+c/(TMath::Power(mesonPt,d)+e)))*1./mesonPt *1./1.6 *1./(2.* TMath::Pi());
3785 // cout << functionResultMC << endl;
3786 } else if (kCaseGen == 3 ){ // HIJING
3787 if ( PDGCode == 111 && fDoReweightHistoMCPi0 && hReweightMCHistPi0!= 0x0){
3788 functionResultMC = hReweightMCHistPi0->Interpolate(mesonPt);
3790 if ( PDGCode == 221 && fDoReweightHistoMCEta && hReweightMCHistEta!= 0x0){
3791 functionResultMC = hReweightMCHistEta->Interpolate(mesonPt);
3793 if ( PDGCode == 310 && fDoReweightHistoMCK0s && hReweightMCHistK0s!= 0x0){
3794 functionResultMC = hReweightMCHistK0s->Interpolate(mesonPt);
3798 Float_t functionResultData = 1;
3799 if (kCaseGen == 1 || kCaseGen == 2 || kCaseGen == 4 ){
3800 Float_t dNdyData = 2.2328;
3801 Float_t nData = 7.1473;
3802 Float_t tData = 0.1346;
3803 if ( PDGCode == 111){
3807 } else if ( PDGCode == 221){
3808 dNdyData = 0.38992; //be careful this fit is not optimal, eta in data still has problems
3812 functionResultData = dNdyData / ( 2 * TMath::Pi())*(nData-1.)*(nData-2.) / (nData*tData*(nData*tData+mesonMass*(nData-2.))) * TMath::Power(1.+(TMath::Sqrt(mesonPt*mesonPt+mesonMass*mesonMass)-mesonMass)/(nData*tData), -nData);
3813 // cout << functionResultData << endl;
3815 if ( PDGCode == 111 && fDoReweightHistoMCPi0 && fFitDataPi0!= 0x0){
3816 functionResultData = fFitDataPi0->Eval(mesonPt);
3818 if ( PDGCode == 221 && fDoReweightHistoMCEta && fFitDataEta!= 0x0){
3819 functionResultData = fFitDataEta->Eval(mesonPt);
3821 if ( PDGCode == 310 && fDoReweightHistoMCK0s && fFitDataK0s!= 0x0){
3822 functionResultData = fFitDataK0s->Eval(mesonPt);
3830 // if ( PDGCode == 111 ){
3831 // if (fModCentralityClass == 1 && fCentralityMin == 0 && fCentralityMax == 1 ){ // 0-5 % PbPb
3832 // a = 25.8747458223;
3833 // b = 5.8761820045;
3834 // c = -33.9928191673;
3835 // d = 3.0731850142;
3836 // e = 13.2500447620;
3837 // } else if (fModCentralityClass == 1 && fCentralityMin == 1 && fCentralityMax == 2){ // 5-10% PbPb
3838 // a = 21.7518148922;
3839 // b = 5.8441200081;
3840 // c = -17.1497051691;
3841 // d = 2.3799090842;
3842 // e = 5.4346404718;
3843 // } else if (fModCentralityClass == 0 && fCentralityMin == 0 && fCentralityMax == 1){ // 0-10% PbPb
3844 // a = 22.9852133622;
3845 // b = 5.8602063916;
3846 // c = -17.0992478654;
3847 // d = 2.4426218039;
3848 // e = 5.1194526345;
3849 // } else if (fModCentralityClass == 0 && fCentralityMin == 1 && fCentralityMax == 2){ // 10-20% PbPb
3850 // a = 19.3237333776;
3851 // b = 5.8145906958;
3852 // c = -13.8316665424;
3853 // d = 2.3737630637;
3854 // e = 4.7690300693;
3855 // } else if (fModCentralityClass == 0 && fCentralityMin == 2 && fCentralityMax == 4){ // 20-40% PbPb
3856 // a = 11.2656032751;
3857 // b = 5.8003194354;
3858 // c = -13.3936105929;
3859 // d = 2.3371452334;
3860 // e = 4.4726244958;
3861 // } else if (fModCentralityClass == 0 && fCentralityMin == 4 && fCentralityMax == 6){ // 40-60% PbPb
3862 // a = 4.1578154081;
3863 // b = 5.6450005163;
3864 // c = -8.4309375240;
3865 // d = 1.8918308704;
3866 // e = 2.9429194709;
3867 // } else if (fModCentralityClass == 0 && fCentralityMin == 6 && fCentralityMax == 8){ // 60-80% PbPb
3868 // a = 1.0635443810;
3869 // b = 5.1337469970;
3870 // c = -8.5906997238;
3871 // d = 2.9794995997;
3872 // e = 3.9294980048;
3873 // } else if (fModCentralityClass == 0 && fCentralityMin == 0 && fCentralityMax == 2){ // 0-20% PbPb
3874 // a = 21.7018745556;
3875 // b = 5.9019352094;
3876 // c = -14.2295510326;
3877 // d = 2.2104490688;
3878 // e = 4.2969671500;
3879 // } else if (fModCentralityClass == 0 && fCentralityMin == 0 && fCentralityMax == 4){ // 0-40% PbPb
3880 // a = 16.8227412106;
3881 // b = 5.8660502207;
3882 // c = -12.0978551215;
3883 // d = 2.1695068981;
3884 // e = 3.5349621182;
3885 // } else if (fModCentralityClass == 0 && fCentralityMin == 0 && fCentralityMax == 8){ // 0-80% PbPb
3886 // a = 9.4675681080;
3887 // b = 5.8114944205;
3888 // c = -10.4901523616;
3889 // d = 2.0607982712;
3890 // e = 2.9262259130;
3891 // } else if (fModCentralityClass == 0 && fCentralityMin == 4 && fCentralityMax == 8){ // 60-80% PbPb
3892 // a = 2.5985551785;
3893 // b = 5.4118895738;
3894 // c = -8.2510958428;
3895 // d = 2.2551249190;
3896 // e = 3.0700919491;
3899 // functionResultData = a*TMath::Power(mesonPt,-1*(b+c/(TMath::Power(mesonPt,d)+e)));
3904 Double_t weight = 1;
3905 if (PDGCode == 111 || PDGCode == 221){
3906 if (functionResultData != 0. && functionResultMC != 0. && isfinite(functionResultData) && isfinite(functionResultMC)){
3907 weight = functionResultData/functionResultMC;
3908 if ( kCaseGen == 3){
3909 if (!(fDoReweightHistoMCPi0 && hReweightMCHistPi0!= 0x0 && PDGCode == 111)){
3913 if (!isfinite(functionResultData)) weight = 1.;
3914 if (!isfinite(weight)) weight = 1.;
3916 } else if (PDGCode == 310 && functionResultMC != 0 && isfinite(functionResultMC)){
3917 weight = functionResultMC;
3920 // if (fModCentralityClass == 0 && fCentralityMin == 4 && fCentralityMax == 6 && PDGCode == 111){
3921 // cout << period.Data() << "\t" << kCaseGen << "\t" <<fModCentralityClass<< "\t" <<fCentralityMin<< "\t" <<fCentralityMax << "\t" << mesonPt << "\t" <<mesonMass<< "\t"<<functionResultData << "\t"<< functionResultMC << "\t" << weight <<endl;
3925 ///________________________________________________________________________
3926 AliConversionCuts* AliConversionCuts::GetStandardCuts2010PbPb(){
3927 //Create and return standard 2010 PbPb cuts
3928 AliConversionCuts *cuts=new AliConversionCuts("StandardCuts2010PbPb","StandardCuts2010PbPb");
3929 if(!cuts->InitializeCutsFromCutString("100000204209297002322000000")){
3930 cout<<"Warning: Initialization of Standardcuts2010PbPb failed"<<endl;}
3934 ///________________________________________________________________________
3935 AliConversionCuts* AliConversionCuts::GetStandardCuts2010pp(){
3936 //Create and return standard 2010 PbPb cuts
3937 AliConversionCuts *cuts=new AliConversionCuts("StandardCuts2010pp","StandardCuts2010pp");
3938 if(!cuts->InitializeCutsFromCutString("000001100209366300380000000")){
3939 cout<<"Warning: Initialization of Standardcuts2010pp failed"<<endl;}
3942 ///________________________________________________________________________
3943 void AliConversionCuts::GetCorrectEtaShiftFromPeriod(TString periodName){
3945 if(periodName.CompareTo("LHC12g") == 0 || //pilot run 2012
3946 periodName.CompareTo("LHC13b") == 0 || //mainly minimum bias
3947 periodName.CompareTo("LHC13c") == 0 || //mainly minimum bias
3948 periodName.CompareTo("LHC13d") == 0 || //mainly triggered
3949 periodName.CompareTo("LHC13e") == 0 || //mainly triggered
3950 periodName.CompareTo("LHC13c3") == 0 || //MC Starlight, anchor LHC13d+e
3951 periodName.CompareTo("LHC13c2") == 0 || //MC Starlight, coherent J/Psi, UPC muon anchor LHC13d+e
3952 periodName.CompareTo("LHC13b4") == 0 || //MC Pythia 6 (Jet-Jet), anchor LHC13b
3953 periodName.CompareTo("LHC13b2_fix_1") == 0 || //MC DPMJET, anchr LHC13b+c
3954 periodName.CompareTo("LHC13b3") == 0 || //MC HIJING, weighted to number of events per run, anchor LHC13b
3955 periodName.CompareTo("LHC13b2") == 0 || // MC DPMJET, wrong energy, anchor LHC13b
3956 periodName.CompareTo("LHC13b2_plus") == 0 || // MC DPMJET, weighted to number event per run, anchor LHC13b
3957 periodName.CompareTo("LHC13c1_bis") == 0 || // MC AMPT fast generation, pT hardbin, anchor ?
3958 periodName.CompareTo("LHC13c1") == 0 || // MC AMPT fast generation, anchor ?
3959 periodName.CompareTo("LHC13b1") == 0 || // MC DPMJET, fragments, with fixed label 0, anchor LHC12g
3960 periodName.CompareTo("LHC12g4b_fix") == 0 || // MC DPMJET, with fixed label 0, anchor LHC12g
3961 periodName.CompareTo("LHC12g1_fix") == 0 || // MC ?, with fixed label 0, anchor LHC12g
3962 periodName.CompareTo("LHC12g4c") == 0 || // MC DPMJET, shifted vertex runs, anchor LHC12g
3963 periodName.CompareTo("LHC12h6") == 0 || // MC muon cocktail, anchor LHC12g
3964 periodName.CompareTo("LHC12g4b") == 0 || // MC DPMJET 3rd iteration, anchor LHC12g
3965 periodName.CompareTo("LHC12g4a") == 0 || // MC DPMJET improved, anchor LHC12g
3966 periodName.CompareTo("LHC12g4") == 0 || // MC DPMJET, anchor LHC12g
3967 periodName.CompareTo("LHC12g5") == 0 || // MC PHOJET, anchor LHC12g
3968 periodName.CompareTo("LHC12g2") == 0 || // MC Starlight background, anchor LHC12g
3969 periodName.CompareTo("LHC12g1") == 0 ) // MC ?, anchor LHC12g
3971 printf(" Gamma Conversion Cuts %s :: pPb Run doing Eta Shift of %f \n\n",(GetCutNumber()).Data(),-0.465);
3972 SetEtaShift(-0.465);
3974 else if(periodName.CompareTo("LHC13f") == 0 ||
3975 periodName.CompareTo("LHC13c6b") == 0 ||// MC Jpsi -> mumu, anchor LHC13f
3976 periodName.CompareTo("LHC13c5") == 0 || //MC Starlight, gamma gamma UPC muon, anchor LHC13f
3977 periodName.CompareTo("LHC13c4") == 0 )//MC Starlight, coherent JPsi, UPC muon, anchor LHC13f
3979 printf(" Gamma Conversion Cuts %s :: Pbp Run doing Eta Shift of %f \n\n",(GetCutNumber()).Data(),0.465);
3980 SetEtaShift(+0.465);
3982 else printf(" Gamma Conversion Cuts %s :: Automatic Eta Shift requested but Period is not known -> No Shift \n\n",(GetCutNumber()).Data());