fixed bugs while streaming histos for weighting
[u/mrichter/AliRoot.git] / PWGGA / GammaConv / AliConversionCuts.cxx
1 /**************************************************************************
2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3  *                                                                        *
4  * Authors: Svein Lindal, Daniel Lohner                                   *
5  * Version 1.0                                                            *
6  *                                                                        *
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  **************************************************************************/
15
16 ////////////////////////////////////////////////
17 //---------------------------------------------
18 // Class handling all kinds of selection cuts for
19 // Gamma Conversion analysis
20 //---------------------------------------------
21 ////////////////////////////////////////////////
22
23 #include "AliConversionCuts.h"
24
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"
33 #include "TH1.h"
34 #include "TH2.h"
35 #include "TF1.h"
36 #include "AliStack.h"
37 #include "AliAODConversionMother.h"
38 #include "TObjString.h"
39 #include "AliAODEvent.h"
40 #include "AliESDEvent.h"
41 #include "AliCentrality.h"
42 #include "TList.h"
43 #include "TFile.h"
44 #include "AliLog.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"
53
54 class iostream;
55
56 using namespace std;
57
58 ClassImp(AliConversionCuts)
59
60
61 const char* AliConversionCuts::fgkCutNames[AliConversionCuts::kNCuts] = {
62    "HeavyIon",//0
63    "CentralityMin",//1
64    "CentralityMax",//2
65    "SelectV0AND",//3
66    "MultiplicityMethod",//4
67    "RemovePileUp",//5
68    "RejectExtraSignals",//6
69    "V0FinderType",//7
70    "EtaCut",//8
71    "MinRCut",//9
72    "SinglePtCut",//10
73    "ClsTPCCut", //11
74    "ededxSigmaCut",//12
75    "pidedxSigmaCut",//13
76    "piMomdedxSigmaCut",//14
77    "piMaxMomdedxSigmaCut",//15
78    "LowPRejectionSigmaCut",//16
79    "TOFelectronPID",//17
80    "QtMaxCut",//18
81    "Chi2GammaCut", //19
82    "PsiPair", //20
83    "DoPhotonAsymmetryCut",//21
84    "CosinePointingAngle", //22
85    "SharedElectronCuts", //23
86    "RejectToCloseV0s", //24
87    "DcaRPrimVtx", //25
88    "DcaZPrimVtx" //26
89 };
90
91
92 //________________________________________________________________________
93 AliConversionCuts::AliConversionCuts(const char *name,const char *title) :
94    AliAnalysisCuts(name,title),
95    fHistograms(NULL),
96    fHeaderList(NULL),
97    fPIDResponse(NULL),
98    fEventQuality(-1),
99    fMaxR(200),
100    fMinR(0),
101    fEtaCut(0.9),
102    fEtaCutMin(-0.1),
103    fPtCut(0.02),
104    fSinglePtCut(0),
105    fMaxZ(1000),
106    fMinClsTPC(0.),
107    fMinClsTPCToF(0.),
108    fLineCutZRSlope(0.),
109    fLineCutZValue(0),
110    fLineCutZRSlopeMin(0.),
111    fLineCutZValueMin(0),
112    fChi2CutConversion(1000),
113    fPIDProbabilityCutNegativeParticle(0),
114    fPIDProbabilityCutPositiveParticle(0),
115    fDodEdxSigmaCut(kTRUE),
116    fDoTOFsigmaCut(kFALSE),
117    fPIDTRDEfficiency(1),
118    fDoTRDPID(kFALSE),
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),
138    fQtMax(100),
139    fHighPtQtMax(0.),
140    fPtBorderForQt(0),
141    fXVertexCut(0.),
142    fYVertexCut(0.),
143    fZVertexCut(0.),
144    fNSigmaMass(0.),
145    fUseEtaMinCut(kFALSE),
146    fUseOnFlyV0Finder(kTRUE),
147    fDoPhotonAsymmetryCut(kTRUE),
148    fMinPPhotonAsymmetryCut(100.),
149    fMinPhotonAsymmetry(0.),
150    fIsHeavyIon(0),
151    fDetectorCentrality(0),
152    fModCentralityClass(0),
153    fMaxVertexZ(10),
154    fCentralityMin(0),
155    fCentralityMax(0),
156    fUseCorrectedTPCClsInfo(kFALSE),
157    fUseTOFpid(kFALSE),
158    fMultiplicityMethod(0),
159    fSpecialTrigger(0),
160    fRemovePileUp(kFALSE),
161    fOpeningAngle(0.005),
162    fPsiPairCut(10000),
163    fCosPAngleCut(10000),
164    fDoToCloseV0sCut(kFALSE),
165    fRejectExtraSignals(0),
166    fminV0Dist(200.),
167    fDoSharedElecCut(kFALSE),
168    fOfflineTriggerMask(0),
169    fHasV0AND(kTRUE),
170    fIsSDDFired(kTRUE),
171    fRandom(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),
179    fnHeaders(0),
180    fNotRejectedStart(NULL),
181    fNotRejectedEnd(NULL),
182    fGeneratorNames(NULL),
183    fCutString(NULL),
184    fUtils(NULL),
185    fEtaShift(0.0),
186    fDoEtaShift(kFALSE),
187    fDoReweightHistoMCPi0(kFALSE),
188    fDoReweightHistoMCEta(kFALSE),
189    fDoReweightHistoMCK0s(kFALSE),
190    fPathTrFReweighting(""),
191    fNameHistoReweightingPi0(""),
192    fNameHistoReweightingEta(""),
193    fNameHistoReweightingK0s(""),
194    fNameFitDataPi0(""),
195    fNameFitDataEta(""),
196    fNameFitDataK0s(""),
197    hdEdxCuts(NULL),
198    hTPCdEdxbefore(NULL),
199    hTPCdEdxafter(NULL),
200    hTPCdEdxSigbefore(NULL),
201    hTPCdEdxSigafter(NULL),
202    hTOFbefore(NULL),
203    hTOFSigbefore(NULL),
204    hTOFSigafter(NULL),
205    hTrackCuts(NULL),
206    hPhotonCuts(NULL),
207    hInvMassbefore(NULL),
208    hArmenterosbefore(NULL),
209    hInvMassafter(NULL),
210    hArmenterosafter(NULL),
211    hAcceptanceCuts(NULL),
212    hCutIndex(NULL),
213    hV0EventCuts(NULL),
214    hCentrality(NULL),
215    hCentralityVsNumberOfPrimaryTracks(NULL),
216    hVertexZ(NULL),
217    hTriggerClass(NULL),
218    hTriggerClassSelected(NULL),
219    hReweightMCHistPi0(NULL),
220    hReweightMCHistEta(NULL),
221    hReweightMCHistK0s(NULL),
222    fFitDataPi0(NULL),
223    fFitDataEta(NULL),
224    fFitDataK0s(NULL),
225    fPreSelCut(kFALSE),
226    fTriggerSelectedManually(kFALSE),
227    fSpecialTriggerName("")
228
229 {
230    InitPIDResponse();
231    for(Int_t jj=0;jj<kNCuts;jj++){fCuts[jj]=0;}
232    fCutString=new TObjString((GetCutNumber()).Data());
233
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);
238
239
240 }
241
242 //________________________________________________________________________
243 AliConversionCuts::AliConversionCuts(const AliConversionCuts &ref) :
244    AliAnalysisCuts(ref),
245    fHistograms(NULL),
246    fHeaderList(ref.fHeaderList),
247    fPIDResponse(NULL),
248    fEventQuality(ref.fEventQuality),
249    fMaxR(ref.fMaxR),
250    fMinR(ref.fMinR),
251    fEtaCut(ref.fEtaCut),
252    fEtaCutMin(ref.fEtaCutMin),
253    fPtCut(ref.fPtCut),
254    fSinglePtCut(ref.fSinglePtCut),
255    fMaxZ(ref.fMaxZ),
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),
288    fQtMax(ref.fQtMax),
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),
333    fCutString(NULL),
334    fUtils(NULL),
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),
347    hdEdxCuts(NULL),
348    hTPCdEdxbefore(NULL),
349    hTPCdEdxafter(NULL),
350    hTPCdEdxSigbefore(NULL),
351    hTPCdEdxSigafter(NULL),
352    hTOFbefore(NULL),
353    hTOFSigbefore(NULL),
354    hTOFSigafter(NULL),
355    hTrackCuts(NULL),
356    hPhotonCuts(NULL),
357    hInvMassbefore(NULL),
358    hArmenterosbefore(NULL),
359    hInvMassafter(NULL),
360    hArmenterosafter(NULL),
361    hAcceptanceCuts(NULL),
362    hCutIndex(NULL),
363    hV0EventCuts(NULL),
364    hCentrality(NULL),
365    hCentralityVsNumberOfPrimaryTracks(NULL),
366    hVertexZ(NULL),
367    hTriggerClass(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)
378 {
379    // Copy Constructor
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())
385
386 }
387
388
389 //________________________________________________________________________
390 AliConversionCuts::~AliConversionCuts() {
391    // Destructor
392    //Deleting fHistograms leads to seg fault it it's added to output collection of a task
393    // if(fHistograms)
394    //    delete fHistograms;
395    // fHistograms = NULL;
396    if(fCutString != NULL){
397       delete fCutString;
398       fCutString = NULL;
399    }
400    if(fElectronLabelArray){
401       delete fElectronLabelArray;
402       fElectronLabelArray = NULL;
403    }
404    if(fNotRejectedStart){
405       delete[] fNotRejectedStart;
406       fNotRejectedStart = NULL;
407    }
408    if(fNotRejectedEnd){
409       delete[] fNotRejectedEnd;
410       fNotRejectedEnd = NULL;
411    }
412    if(fGeneratorNames){
413       delete[] fGeneratorNames;
414       fGeneratorNames = NULL;
415    }
416    if(fUtils){
417      delete fUtils;
418      fUtils = NULL;
419    }
420
421 }
422
423 //________________________________________________________________________
424 void AliConversionCuts::InitCutHistograms(TString name, Bool_t preCut){
425
426    // Initialize Cut Histograms for QA (only initialized and filled if function is called)
427    TH1::AddDirectory(kFALSE);
428
429    if(fHistograms != NULL){
430       delete fHistograms;
431       fHistograms=NULL;
432    }
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()));
438    }
439
440    if (hReweightMCHistPi0){
441       hReweightMCHistPi0->SetName("MCInputForWeightingPi0");
442       fHistograms->Add(hReweightMCHistPi0);
443    }
444    if (hReweightMCHistEta){
445       hReweightMCHistEta->SetName("MCInputForWeightingEta");
446       fHistograms->Add(hReweightMCHistEta);
447    }
448    if (hReweightMCHistK0s){
449       hReweightMCHistK0s->SetName("MCInputForWeightingK0s");
450       fHistograms->Add(hReweightMCHistK0s);
451    }
452    if (fFitDataPi0){
453       fFitDataPi0->SetName("DataFitForWeightingPi0");
454       fHistograms->Add(fFitDataPi0);
455    }
456    if (fFitDataEta){
457       fFitDataEta->SetName("DataFitForWeightingEta");
458       fHistograms->Add(fFitDataEta);
459    }
460    if (fFitDataK0s){
461       fFitDataK0s->SetName("DataFitForWeightingK0s");
462       fHistograms->Add(fFitDataK0s);
463    }
464    // IsPhotonSelected
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);
475
476    // Track Cuts
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);
487
488    // Photon Cuts
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);
503
504    if(preCut){
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);
509    }
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);
514
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);
525
526    // dEdx Cuts
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);
539
540    TAxis *AxisBeforedEdx = NULL;
541    TAxis *AxisBeforedEdxSig = NULL;
542    TAxis *AxisBeforeTOF = NULL;
543    TAxis *AxisBeforeTOFSig = NULL;
544    if(preCut){
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();
551
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();
558
559    }
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);
562
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);
565
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);
568
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];
574    newBins[0] = from;
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);
582    if(preCut){
583       AxisBeforedEdx->Set(bins, newBins);
584       AxisBeforeTOF->Set(bins, newBins);
585       AxisBeforedEdxSig->Set(bins, newBins);
586       AxisBeforeTOFSig->Set(bins, newBins);
587    }
588    delete [] newBins;
589
590    // Event Cuts and Info
591    if(preCut){
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);
601
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);
608
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);
646    }
647    if(!preCut){
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);
684    }
685    TH1::AddDirectory(kTRUE);
686 }
687
688 //________________________________________________________________________
689 Bool_t AliConversionCuts::InitPIDResponse(){
690    // Set Pointer to AliPIDResponse
691
692    AliAnalysisManager *man=AliAnalysisManager::GetAnalysisManager();
693    if(man) {
694       AliInputEventHandler* inputHandler = (AliInputEventHandler*) (man->GetInputEventHandler());
695       fPIDResponse = (AliPIDResponse*)inputHandler->GetPIDResponse();
696       if(fPIDResponse)return kTRUE;
697
698    }
699
700
701    return kFALSE;
702 }
703 ///________________________________________________________________________
704 Bool_t AliConversionCuts::EventIsSelected(AliVEvent *fInputEvent, AliVEvent *fMCEvent){
705    // Process Event Selection
706
707    Int_t cutindex=0;
708    if(hV0EventCuts)hV0EventCuts->Fill(cutindex);
709    cutindex++;
710
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());
715       if (!mcHandler){
716          fEventQuality = 2;
717          return kFALSE;
718       }
719       if (!mcHandler->InitOk() ){
720          fEventQuality = 2;
721          return kFALSE;
722       }
723       if (!mcHandler->TreeK() ){
724          fEventQuality = 2;
725          return kFALSE;
726       }
727       if (!mcHandler->TreeTR() ) {
728          fEventQuality = 2;
729          return kFALSE;
730       }
731    }
732
733    // Event Trigger
734 //    cout << "before event trigger" << endl;
735    if(!IsTriggerSelected(fInputEvent)){
736       if(hV0EventCuts)hV0EventCuts->Fill(cutindex);
737       fEventQuality = 3;
738       return kFALSE;
739    }
740    cutindex++;
741
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);
746    }
747 //   cout << "event number " << ((AliESDEvent*)fInputEvent)->GetEventNumberInFile() << " entered"<< endl;
748
749
750    // Number of Contributors Cut
751    if(GetNumberOfContributorsVtx(fInputEvent)<=0) {
752       if(hV0EventCuts)hV0EventCuts->Fill(cutindex);
753       fEventQuality = 5;
754       return kFALSE;
755    }
756    cutindex++;
757
758    // Z Vertex Position Cut
759    if(!VertexZCut(fInputEvent)){
760       if(hV0EventCuts)hV0EventCuts->Fill(cutindex);
761       fEventQuality = 4;
762       return kFALSE;
763    }
764    cutindex++;
765
766    // Pile Up Rejection
767
768    if(fRemovePileUp){
769       if(fInputEvent->IsPileupFromSPD(3,0.8,3.,2.,5.)){
770          if(hV0EventCuts)hV0EventCuts->Fill(cutindex);
771          fEventQuality = 6;
772          return kFALSE;
773       }
774    }
775    cutindex++;
776
777    // Centrality Selection
778    if(!IsCentralitySelected(fInputEvent,fMCEvent)){
779       if(hV0EventCuts)hV0EventCuts->Fill(cutindex);
780       fEventQuality = 1;
781       return kFALSE;
782    }
783    cutindex++;
784
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());
793    fEventQuality = 0;
794    return kTRUE;
795 }
796
797 ///________________________________________________________________________
798 Bool_t AliConversionCuts::PhotonIsSelectedMC(TParticle *particle,AliStack *fMCStack,Bool_t checkForConvertedGamma){
799    // MonteCarlo Photon Selection
800
801    if(!fMCStack)return kFALSE;
802
803    if (particle->GetPdgCode() == 22){
804
805
806       if( particle->Eta() > (fEtaCut + fEtaShift) || particle->Eta() < (-fEtaCut + fEtaShift) )
807          return kFALSE;
808       if(fEtaCutMin>-0.1){
809          if( particle->Eta() < (fEtaCutMin + fEtaShift) && particle->Eta() > (-fEtaCutMin + fEtaShift) )
810             return kFALSE;
811       }
812
813       if(particle->GetMother(0) >-1 && fMCStack->Particle(particle->GetMother(0))->GetPdgCode() == 22){
814          return kFALSE; // no photon as mothers!
815       }
816
817       if(particle->GetMother(0) >= fMCStack->GetNprimary()){
818          return kFALSE; // the gamma has a mother, and it is not a primary particle
819       }
820
821       if(!checkForConvertedGamma) return kTRUE; // return in case of accepted gamma
822
823       // looking for conversion gammas (electron + positron from pairbuilding (= 5) )
824       TParticle* ePos = NULL;
825       TParticle* eNeg = NULL;
826
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){
832                   eNeg = tmpDaughter;
833                } else if(tmpDaughter->GetPdgCode() == -11){
834                   ePos = tmpDaughter;
835                }
836             }
837          }
838       }
839
840       if(ePos == NULL || eNeg == NULL){ // means we do not have two daughters from pair production
841          return kFALSE;
842       }
843
844       if(ePos->Pt()<fSinglePtCut || eNeg->Pt()<fSinglePtCut){
845          return kFALSE; // no reconstruction below the Pt cut
846       }
847
848       if( ePos->Eta() > (fEtaCut + fEtaShift) || ePos->Eta() < (-fEtaCut + fEtaShift) ||
849           eNeg->Eta() > (fEtaCut + fEtaShift) || eNeg->Eta() < (-fEtaCut + fEtaShift) )
850          return kFALSE;
851
852       if(fEtaCutMin > -0.1){
853          if( (ePos->Eta() < (fEtaCutMin + fEtaShift) && ePos->Eta() > (-fEtaCutMin + fEtaShift)) ||
854              (eNeg->Eta() < (fEtaCutMin + fEtaShift) && eNeg->Eta() > (-fEtaCutMin + fEtaShift)) )
855             return kFALSE;
856       }
857
858       if(ePos->R()>fMaxR){
859          return kFALSE; // cuts on distance from collision point
860       }
861
862       if(abs(ePos->Vz()) > fMaxZ){
863          return kFALSE;  // outside material
864       }
865       if(abs(eNeg->Vz()) > fMaxZ){
866          return kFALSE;  // outside material
867       }
868
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)){
872          return kFALSE;
873       }
874
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)){
878          return kFALSE;
879       }
880
881       return kTRUE;
882       //if(AcceptanceCut(particle,ePos,eNeg))return kTRUE;
883    }
884    return kFALSE;
885 }
886 ///________________________________________________________________________
887 Bool_t AliConversionCuts::PhotonIsSelectedAODMC(AliAODMCParticle *particle,TClonesArray *aodmcArray,Bool_t checkForConvertedGamma){
888    // MonteCarlo Photon Selection
889
890    if(!aodmcArray)return kFALSE;
891
892    if (particle->GetPdgCode() == 22){
893       if( particle->Eta() > (fEtaCut + fEtaShift) || particle->Eta() < (-fEtaCut + fEtaShift) )
894          return kFALSE;
895       if(fEtaCutMin>-0.1){
896          if( particle->Eta() < (fEtaCutMin + fEtaShift) && particle->Eta() > (-fEtaCutMin + fEtaShift) )
897             return kFALSE;
898       }
899
900       if(particle->GetMother() > -1){
901          if((static_cast<AliAODMCParticle*>(aodmcArray->At(particle->GetMother())))->GetPdgCode() == 22){
902             return kFALSE; // no photon as mothers!
903          }
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
906          }
907       }
908
909       if(!checkForConvertedGamma) return kTRUE; // return in case of accepted gamma
910
911       // looking for conversion gammas (electron + positron from pairbuilding (= 5) )
912       AliAODMCParticle* ePos = NULL;
913       AliAODMCParticle* eNeg = NULL;
914
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){
921                   eNeg = tmpDaughter;
922                } else if(tmpDaughter->GetPdgCode() == -11){
923                   ePos = tmpDaughter;
924                }
925             }
926          }
927       }
928
929       if(ePos == NULL || eNeg == NULL){ // means we do not have two daughters from pair production
930          return kFALSE;
931       }
932
933       if(ePos->Pt()<fSinglePtCut || eNeg->Pt()<fSinglePtCut){
934          return kFALSE; // no reconstruction below the Pt cut
935       }
936
937       if( ePos->Eta() > (fEtaCut + fEtaShift) || ePos->Eta() < (-fEtaCut + fEtaShift) ||
938           eNeg->Eta() > (fEtaCut + fEtaShift) || eNeg->Eta() < (-fEtaCut + fEtaShift) )
939          return kFALSE;
940
941       if(fEtaCutMin > -0.1){
942          if( (ePos->Eta() < (fEtaCutMin + fEtaShift) && ePos->Eta() > (-fEtaCutMin + fEtaShift)) ||
943              (eNeg->Eta() < (fEtaCutMin + fEtaShift) && eNeg->Eta() > (-fEtaCutMin + fEtaShift)) )
944             return kFALSE;
945       }
946
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()) );
949
950       if(rPos>fMaxR){
951          return kFALSE; // cuts on distance from collision point
952       }
953       if(abs(ePos->Zv()) > fMaxZ){
954          return kFALSE;  // outside material
955       }
956       if(abs(eNeg->Zv()) > fMaxZ){
957          return kFALSE;  // outside material
958       }
959
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)){
963          return kFALSE;
964       }
965
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)){
969          return kFALSE;
970       }
971
972       return kTRUE;
973       //if(AcceptanceCut(particle,ePos,eNeg))return kTRUE;
974    }
975    return kFALSE;
976 }
977 ///________________________________________________________________________
978 Bool_t AliConversionCuts::PhotonCuts(AliConversionPhotonBase *photon,AliVEvent *event)
979 {   // Specific Photon Cuts
980
981    Int_t cutIndex = 0;
982    if(hPhotonCuts)hPhotonCuts->Fill(cutIndex);
983    cutIndex++;
984
985    // Fill Histos before Cuts
986    if(hInvMassbefore)hInvMassbefore->Fill(photon->GetMass());
987
988    if(hArmenterosbefore)hArmenterosbefore->Fill(photon->GetArmenterosAlpha(),photon->GetArmenterosQt());
989
990    // Gamma selection based on QT from Armenteros
991    if(fDoQtGammaSelection == kTRUE){
992       if(!ArmenterosQtCut(photon)){
993          if(hPhotonCuts)hPhotonCuts->Fill(cutIndex); //1
994          return kFALSE;
995       }
996    }
997    cutIndex++; //2
998
999    // Chi Cut
1000    if(photon->GetChi2perNDF() > fChi2CutConversion || photon->GetChi2perNDF() <=0){
1001       {
1002          if(hPhotonCuts)hPhotonCuts->Fill(cutIndex); //2
1003          return kFALSE;
1004       }
1005    }
1006    cutIndex++;//3
1007
1008    // Reconstruction Acceptance Cuts
1009    if(!AcceptanceCuts(photon)){
1010       if(hPhotonCuts)hPhotonCuts->Fill(cutIndex); //3
1011       return kFALSE;
1012    }
1013
1014    cutIndex++; //4
1015    // Asymmetry Cut
1016    if(fDoPhotonAsymmetryCut == kTRUE){
1017       if(!AsymmetryCut(photon,event)){
1018          if(hPhotonCuts)hPhotonCuts->Fill(cutIndex); //4
1019          return kFALSE;
1020       }
1021    }
1022
1023    //Check the pid probability
1024    cutIndex++; //5
1025    if(!PIDProbabilityCut(photon, event)) {
1026       if(hPhotonCuts)hPhotonCuts->Fill(cutIndex); //5
1027       return kFALSE;
1028    }
1029
1030    cutIndex++; //6
1031    if(!CorrectedTPCClusterCut(photon, event)) {
1032       if(hPhotonCuts)hPhotonCuts->Fill(cutIndex); //6
1033       return kFALSE;
1034    }
1035
1036
1037    cutIndex++; //7
1038    if(!PsiPairCut(photon)) {
1039       if(hPhotonCuts)hPhotonCuts->Fill(cutIndex); //7
1040       return kFALSE;
1041    }
1042
1043    cutIndex++; //8
1044    if(!CosinePAngleCut(photon, event)) {
1045       if(hPhotonCuts)hPhotonCuts->Fill(cutIndex); //8
1046       return kFALSE;
1047    }
1048
1049    AliAODConversionPhoton* photonAOD = dynamic_cast<AliAODConversionPhoton*>(photon);
1050    if (photonAOD){
1051       photonAOD->CalculateDistanceOfClossetApproachToPrimVtx(event->GetPrimaryVertex());
1052
1053       cutIndex++; //9
1054       if(photonAOD->GetDCArToPrimVtx() > fDCARPrimVtxCut) { //DCA R cut of photon to primary vertex
1055          if(hPhotonCuts)hPhotonCuts->Fill(cutIndex); //9
1056          return kFALSE;
1057       }
1058
1059       cutIndex++; //10
1060       if(abs(photonAOD->GetDCAzToPrimVtx()) > fDCAZPrimVtxCut) { //DCA Z cut of photon to primary vertex
1061          if(hPhotonCuts)hPhotonCuts->Fill(cutIndex); //10
1062          return kFALSE;
1063       }
1064    } else {
1065       cutIndex++; //9
1066       cutIndex++; //10
1067    }
1068    cutIndex++; //11
1069    if(hPhotonCuts)hPhotonCuts->Fill(cutIndex); //11
1070
1071    // Histos after Cuts
1072    if(hInvMassafter)hInvMassafter->Fill(photon->GetMass());
1073    if(hArmenterosafter)hArmenterosafter->Fill(photon->GetArmenterosAlpha(),photon->GetArmenterosQt());
1074
1075    return kTRUE;
1076
1077 }
1078
1079 ///________________________________________________________________________
1080 Bool_t AliConversionCuts::CorrectedTPCClusterCut(AliConversionPhotonBase *photon, AliVEvent * event)
1081 {   //Cut on corrected TPC Cluster Info
1082
1083    AliVTrack * negTrack = GetTrack(event, photon->GetTrackLabelNegative());
1084    AliVTrack * posTrack = GetTrack(event, photon->GetTrackLabelPositive());
1085
1086    if(!negTrack||!posTrack)return kFALSE;
1087
1088    Double_t negclsToF=0;
1089
1090    if (!fUseCorrectedTPCClsInfo ){
1091       if(negTrack->GetTPCNclsF()!=0){
1092          negclsToF = (Double_t)negTrack->GetNcls(1)/(Double_t)negTrack->GetTPCNclsF();}// Ncluster/Nfindablecluster
1093    }
1094    else {
1095       negclsToF = negTrack->GetTPCClusterInfo(2,0,GetFirstTPCRow(photon->GetConversionRadius()));
1096    }
1097
1098    Double_t posclsToF = 0.;
1099    if (!fUseCorrectedTPCClsInfo ){
1100       if(posTrack->GetTPCNclsF()!=0){
1101          posclsToF = (Double_t)posTrack->GetNcls(1)/(Double_t)posTrack->GetTPCNclsF();
1102       }
1103    }else{
1104       posclsToF = posTrack->GetTPCClusterInfo(2,0,GetFirstTPCRow(photon->GetConversionRadius()));
1105    }
1106
1107    if( negclsToF < fMinClsTPCToF || posclsToF < fMinClsTPCToF ){
1108       return kFALSE;
1109    }
1110
1111    return kTRUE;
1112 }
1113
1114 ///________________________________________________________________________
1115 Bool_t AliConversionCuts::PhotonIsSelected(AliConversionPhotonBase *photon, AliVEvent * event)
1116 {
1117    //Selection of Reconstructed Photons
1118
1119    FillPhotonCutIndex(kPhotonIn);
1120
1121    if(event->IsA()==AliESDEvent::Class()) {
1122      if(!SelectV0Finder( ( ((AliESDEvent*)event)->GetV0(photon->GetV0Index()))->GetOnFlyStatus() ) ){
1123          FillPhotonCutIndex(kOnFly);
1124          return kFALSE;
1125       }
1126    }
1127    // else if(event->IsA()==AliAODEvent::Class()) {
1128    //    if(!SelectV0Finder( ( ((AliAODEvent*)event)->GetV0(photon->GetV0Index())) ) ){
1129    //       FillPhotonCutIndex(kOnFly);
1130    //       return kFALSE;
1131    //    }
1132    // }
1133
1134    // Get Tracks
1135    AliVTrack * negTrack = GetTrack(event, photon->GetTrackLabelNegative());
1136    AliVTrack * posTrack = GetTrack(event, photon->GetTrackLabelPositive());
1137
1138    if(!negTrack || !posTrack) {
1139       FillPhotonCutIndex(kNoTracks);
1140       return kFALSE;
1141    }
1142    photon->DeterminePhotonQuality(negTrack,posTrack);
1143    // Track Cuts
1144    if(!TracksAreSelected(negTrack, posTrack)){
1145       FillPhotonCutIndex(kTrackCuts);
1146       return kFALSE;
1147    }
1148
1149    // dEdx Cuts
1150    if(!dEdxCuts(negTrack) || !dEdxCuts(posTrack)) {
1151       FillPhotonCutIndex(kdEdxCuts);
1152       return kFALSE;
1153    }
1154
1155    // Photon Cuts
1156    if(!PhotonCuts(photon,event)){
1157       FillPhotonCutIndex(kPhotonCuts);
1158       return kFALSE;
1159    }
1160
1161    // Photon passed cuts
1162    FillPhotonCutIndex(kPhotonOut);
1163    return kTRUE;
1164 }
1165
1166 ///________________________________________________________________________
1167 Bool_t AliConversionCuts::ArmenterosQtCut(AliConversionPhotonBase *photon)
1168 {   // Armenteros Qt Cut
1169
1170    if(fDoHighPtQtGammaSelection){
1171       if(photon->GetPhotonPt() < fPtBorderForQt){
1172          if(photon->GetArmenterosQt()>fQtMax){
1173             return kFALSE;
1174          }
1175       } else {
1176          if(photon->GetArmenterosQt()>fHighPtQtMax){
1177             return kFALSE;
1178          }
1179       }
1180    } else {
1181
1182       if(photon->GetArmenterosQt()>fQtMax){
1183          return kFALSE;
1184       }
1185    }
1186    return kTRUE;
1187 }
1188
1189
1190 ///________________________________________________________________________
1191 Bool_t AliConversionCuts::AcceptanceCuts(AliConversionPhotonBase *photon) {
1192    // Exclude certain areas for photon reconstruction
1193
1194    Int_t cutIndex=0;
1195    if(hAcceptanceCuts)hAcceptanceCuts->Fill(cutIndex);
1196    cutIndex++;
1197
1198    if(photon->GetConversionRadius()>fMaxR){ // cuts on distance from collision point
1199       if(hAcceptanceCuts)hAcceptanceCuts->Fill(cutIndex);
1200       return kFALSE;
1201    }
1202    cutIndex++;
1203
1204    if(photon->GetConversionRadius()<fMinR){ // cuts on distance from collision point
1205       if(hAcceptanceCuts)hAcceptanceCuts->Fill(cutIndex);
1206       return kFALSE;
1207    }
1208    cutIndex++;
1209
1210    if(photon->GetConversionRadius() <= ((abs(photon->GetConversionZ())*fLineCutZRSlope)-fLineCutZValue)){
1211       if(hAcceptanceCuts)hAcceptanceCuts->Fill(cutIndex);
1212       return kFALSE;
1213    }
1214    else if (fUseEtaMinCut &&  photon->GetConversionRadius() >= ((abs(photon->GetConversionZ())*fLineCutZRSlopeMin)-fLineCutZValueMin )){
1215       if(hAcceptanceCuts)hAcceptanceCuts->Fill(cutIndex);
1216       return kFALSE;
1217    }
1218    cutIndex++;
1219
1220    if(abs(photon->GetConversionZ()) > fMaxZ ){ // cuts out regions where we do not reconstruct
1221       if(hAcceptanceCuts)hAcceptanceCuts->Fill(cutIndex);
1222       return kFALSE;
1223    }
1224    cutIndex++;
1225
1226
1227    if( photon->GetPhotonEta() > (fEtaCut + fEtaShift)    || photon->GetPhotonEta() < (-fEtaCut + fEtaShift) ){
1228       if(hAcceptanceCuts)hAcceptanceCuts->Fill(cutIndex);
1229       return kFALSE;
1230    }
1231    if(fEtaCutMin>-0.1){
1232       if( photon->GetPhotonEta() < (fEtaCutMin + fEtaShift) && photon->GetPhotonEta() > (-fEtaCutMin + fEtaShift) ){
1233          if(hAcceptanceCuts)hAcceptanceCuts->Fill(cutIndex);
1234          return kFALSE;
1235       }
1236    }
1237    cutIndex++;
1238
1239    if(photon->GetPhotonPt()<fPtCut){
1240       if(hAcceptanceCuts)hAcceptanceCuts->Fill(cutIndex);
1241       return kFALSE;
1242    }
1243    cutIndex++;
1244
1245    if(hAcceptanceCuts)hAcceptanceCuts->Fill(cutIndex);
1246
1247    return kTRUE;
1248 }
1249
1250
1251 ///________________________________________________________________________
1252 Bool_t AliConversionCuts::SpecificTrackCuts(AliAODTrack * negTrack, AliAODTrack * posTrack,Int_t &cutIndex) {
1253    // Track Cuts which require AOD/ESD specific implementation
1254
1255    if( !negTrack->IsOn(AliESDtrack::kTPCrefit)  || !posTrack->IsOn(AliESDtrack::kTPCrefit)   )  {
1256       if(hTrackCuts)hTrackCuts->Fill(cutIndex);
1257       return kFALSE;
1258    }
1259    cutIndex++;
1260
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);
1265       return kFALSE;
1266    }
1267    return kTRUE;
1268
1269 }
1270
1271
1272 ///________________________________________________________________________
1273 Bool_t AliConversionCuts::SpecificTrackCuts(AliESDtrack * negTrack, AliESDtrack * posTrack,Int_t &cutIndex) {
1274    // Track Cuts which require AOD/ESD specific implementation
1275
1276    if( !negTrack->IsOn(AliESDtrack::kTPCrefit)  || !posTrack->IsOn(AliESDtrack::kTPCrefit)   )  {
1277       if(hTrackCuts)hTrackCuts->Fill(cutIndex);
1278       return kFALSE;
1279    }
1280    cutIndex++;
1281
1282    if(negTrack->GetKinkIndex(0) > 0  || posTrack->GetKinkIndex(0) > 0 ) {
1283       if(hTrackCuts)hTrackCuts->Fill(cutIndex);
1284       return kFALSE;
1285    }
1286    return kTRUE;
1287 }
1288
1289
1290
1291 ///________________________________________________________________________
1292 Bool_t AliConversionCuts::TracksAreSelected(AliVTrack * negTrack, AliVTrack * posTrack) {
1293    // Track Selection for Photon Reconstruction
1294
1295    Int_t cutIndex=0;
1296    if(hTrackCuts)hTrackCuts->Fill(cutIndex);
1297    cutIndex++;
1298
1299    // avoid like sign
1300    if(negTrack->Charge() == posTrack->Charge()) {
1301       if(hTrackCuts)hTrackCuts->Fill(cutIndex);
1302       return kFALSE;
1303    }
1304    cutIndex++;
1305
1306    // Number of TPC Clusters
1307
1308
1309    if( negTrack->GetNcls(1) < fMinClsTPC || posTrack->GetNcls(1) < fMinClsTPC ) {
1310       if(hTrackCuts)hTrackCuts->Fill(cutIndex);
1311       return kFALSE;
1312    }
1313    cutIndex++;
1314
1315    // Acceptance
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);
1319       return kFALSE;
1320    }
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);
1325          return kFALSE;
1326       }
1327    }
1328    cutIndex++;
1329
1330    // Single Pt Cut
1331    if( negTrack->Pt()< fSinglePtCut || posTrack->Pt()< fSinglePtCut){
1332       if(hTrackCuts)hTrackCuts->Fill(cutIndex);
1333       return kFALSE;
1334    }
1335    cutIndex++;
1336
1337    // AOD ESD specific cuts
1338    Bool_t passCuts = kTRUE;
1339
1340    if(negTrack->IsA()==AliAODTrack::Class()) {
1341       passCuts = passCuts * SpecificTrackCuts(static_cast<AliAODTrack*>(negTrack), static_cast<AliAODTrack*>(posTrack),cutIndex);
1342    } else {
1343       passCuts = passCuts * SpecificTrackCuts(static_cast<AliESDtrack*>(negTrack), static_cast<AliESDtrack*>(posTrack),cutIndex);
1344    }
1345
1346    if(!passCuts){
1347       if(hTrackCuts)hTrackCuts->Fill(cutIndex);
1348       return kFALSE;
1349    }
1350    cutIndex++;
1351
1352    if(hTrackCuts)hTrackCuts->Fill(cutIndex);
1353
1354    return kTRUE;
1355
1356 }
1357
1358 ///________________________________________________________________________
1359 Bool_t AliConversionCuts::dEdxCuts(AliVTrack *fCurrentTrack){
1360    // Electron Identification Cuts for Photon reconstruction
1361
1362    if(!fPIDResponse){InitPIDResponse();}// Try to reinitialize PID Response
1363    if(!fPIDResponse){AliError("No PID Response"); return kTRUE;}// if still missing fatal error
1364
1365    Int_t cutIndex=0;
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());
1369    cutIndex++;
1370
1371
1372    if(fDodEdxSigmaCut == kTRUE){
1373       // TPC Electron Line
1374       if( fPIDResponse->NumberOfSigmasTPC(fCurrentTrack,AliPID::kElectron)<fPIDnSigmaBelowElectronLine ||
1375           fPIDResponse->NumberOfSigmasTPC(fCurrentTrack,AliPID::kElectron)>fPIDnSigmaAboveElectronLine){
1376
1377          if(hdEdxCuts)hdEdxCuts->Fill(cutIndex);
1378          return kFALSE;
1379       }
1380       cutIndex++;
1381
1382       // TPC Pion Line
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){
1387
1388             if(hdEdxCuts)hdEdxCuts->Fill(cutIndex);
1389             return kFALSE;
1390          }
1391       }
1392       cutIndex++;
1393
1394       // High Pt Pion rej
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){
1399
1400             if(hdEdxCuts)hdEdxCuts->Fill(cutIndex);
1401             return kFALSE;
1402          }
1403       }
1404       cutIndex++;
1405    }
1406    else{cutIndex+=3;}
1407
1408    if(fDoKaonRejectionLowP == kTRUE){
1409       if(fCurrentTrack->P()<fPIDMinPKaonRejectionLowP ){
1410          if( abs(fPIDResponse->NumberOfSigmasTPC(fCurrentTrack,AliPID::kKaon))<fPIDnSigmaAtLowPAroundKaonLine){
1411
1412             if(hdEdxCuts)hdEdxCuts->Fill(cutIndex);
1413             return kFALSE;
1414          }
1415       }
1416    }
1417    cutIndex++;
1418
1419    if(fDoProtonRejectionLowP == kTRUE){
1420       if( fCurrentTrack->P()<fPIDMinPProtonRejectionLowP ){
1421          if( abs(fPIDResponse->NumberOfSigmasTPC(fCurrentTrack,AliPID::kProton))<fPIDnSigmaAtLowPAroundProtonLine){
1422
1423             if(hdEdxCuts)hdEdxCuts->Fill(cutIndex);
1424             return kFALSE;
1425          }
1426       }
1427    }
1428    cutIndex++;
1429
1430    if(fDoPionRejectionLowP == kTRUE){
1431       if( fCurrentTrack->P()<fPIDMinPPionRejectionLowP ){
1432          if( abs(fPIDResponse->NumberOfSigmasTPC(fCurrentTrack,AliPID::kPion))<fPIDnSigmaAtLowPAroundPionLine){
1433
1434             if(hdEdxCuts)hdEdxCuts->Fill(cutIndex);
1435             return kFALSE;
1436          }
1437       }
1438    }
1439    cutIndex++;
1440
1441
1442    // cout<<"Start"<<endl;
1443    // AliPIDResponse::EDetPidStatus status=fPIDResponse->CheckPIDStatus(AliPIDResponse::kTOF,fCurrentTrack);
1444
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){
1451
1452    //    }
1453    // }
1454
1455    if((fCurrentTrack->GetStatus() & AliESDtrack::kTOFpid) && !(fCurrentTrack->GetStatus() & AliESDtrack::kTOFmismatch)){
1456       if(hTOFbefore){
1457          Double_t t0 = fPIDResponse->GetTOFResponse().GetStartTime(fCurrentTrack->P());
1458          Double_t times[5];
1459          fCurrentTrack->GetIntegratedTimes(times);
1460          Double_t TOFsignal = fCurrentTrack->GetTOFsignal();
1461          Double_t dT = TOFsignal - t0 - times[0];
1462          hTOFbefore->Fill(fCurrentTrack->P(),dT);
1463       }
1464       if(hTOFSigbefore) hTOFSigbefore->Fill(fCurrentTrack->P(),fPIDResponse->NumberOfSigmasTOF(fCurrentTrack, AliPID::kElectron));
1465       if(fUseTOFpid){
1466          if(fPIDResponse->NumberOfSigmasTOF(fCurrentTrack, AliPID::kElectron)>fTofPIDnSigmaAboveElectronLine ||
1467             fPIDResponse->NumberOfSigmasTOF(fCurrentTrack, AliPID::kElectron)<fTofPIDnSigmaBelowElectronLine ){
1468             if(hdEdxCuts)hdEdxCuts->Fill(cutIndex);
1469             return kFALSE;
1470          }
1471       }
1472       if(hTOFSigafter)hTOFSigafter->Fill(fCurrentTrack->P(),fPIDResponse->NumberOfSigmasTOF(fCurrentTrack, AliPID::kElectron));
1473    }
1474    cutIndex++;
1475
1476    // Apply TRD PID
1477    if(fDoTRDPID){
1478       if(!fPIDResponse->IdentifiedAsElectronTRD(fCurrentTrack,fPIDTRDEfficiency)){
1479          if(hdEdxCuts)hdEdxCuts->Fill(cutIndex);
1480          return kFALSE;
1481       }
1482    }
1483    cutIndex++;
1484
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());
1488    return kTRUE;
1489 }
1490
1491 ///________________________________________________________________________
1492 Bool_t AliConversionCuts::AsymmetryCut(AliConversionPhotonBase * photon,AliVEvent *event) {
1493    // Cut on Energy Assymetry
1494
1495    for(Int_t ii=0;ii<2;ii++){
1496
1497       AliVTrack *track=GetTrack(event,photon->GetTrackLabel(ii));
1498
1499       if( track->P() > fMinPPhotonAsymmetryCut ){
1500          Double_t trackNegAsy=0;
1501          if (photon->GetPhotonP()!=0.){
1502             trackNegAsy= track->P()/photon->GetPhotonP();
1503          }
1504
1505          if( trackNegAsy<fMinPhotonAsymmetry ||trackNegAsy>(1.- fMinPhotonAsymmetry)){
1506             return kFALSE;
1507          }
1508       }
1509
1510    }
1511    return kTRUE;
1512 }
1513
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)
1519
1520    AliESDEvent * esdEvent = dynamic_cast<AliESDEvent*>(event);
1521    if(esdEvent) {
1522       if(label > event->GetNumberOfTracks() ) return NULL;
1523       AliESDtrack * track = esdEvent->GetTrack(label);
1524       return track;
1525
1526    } else {
1527       AliVTrack * track = 0x0;
1528       if(((AliV0ReaderV1*)AliAnalysisManager::GetAnalysisManager()->GetTask("V0ReaderV1"))->AreAODsRelabeled()){
1529          track = dynamic_cast<AliVTrack*>(event->GetTrack(label));
1530          return track;
1531       }
1532       else{
1533          for(Int_t ii=0; ii<event->GetNumberOfTracks(); ii++) {
1534             track = dynamic_cast<AliVTrack*>(event->GetTrack(ii));
1535             if(track){
1536                if(track->GetID() == label) {
1537                   return track;
1538                }
1539             }
1540          }
1541       }
1542    }
1543    //AliDebug(5,(Form("track not found %d %d",label,event->GetNumberOfTracks()));
1544    return NULL;
1545 }
1546
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)
1552
1553    if(event) {
1554       if(label > event->GetNumberOfTracks() ) return NULL;
1555       AliESDtrack * track = event->GetTrack(label);
1556       return track;
1557    }
1558    //AliDebug(5,(Form("track not found %d %d",label,event->GetNumberOfTracks()));
1559    return NULL;
1560 }
1561
1562
1563
1564 ///________________________________________________________________________
1565 Bool_t AliConversionCuts::PIDProbabilityCut(AliConversionPhotonBase *photon, AliVEvent * event){
1566    // Cut on Electron Probability for Photon Reconstruction
1567
1568    AliESDEvent * esdEvent = dynamic_cast<AliESDEvent*>(event);
1569
1570    if(esdEvent){
1571
1572       Bool_t iResult=kFALSE;
1573
1574       Double_t *posProbArray = new Double_t[AliPID::kSPECIES];
1575       Double_t *negProbArray = new Double_t[AliPID::kSPECIES];
1576
1577       AliESDtrack* negTrack   = esdEvent->GetTrack(photon->GetTrackLabelNegative());
1578       AliESDtrack* posTrack   = esdEvent->GetTrack(photon->GetTrackLabelPositive());
1579
1580       if(negProbArray && posProbArray){
1581
1582          negTrack->GetTPCpid(negProbArray);
1583          posTrack->GetTPCpid(posProbArray);
1584
1585          if(negProbArray[AliPID::kElectron]>=fPIDProbabilityCutNegativeParticle && posProbArray[AliPID::kElectron]>=fPIDProbabilityCutPositiveParticle){
1586             iResult=kTRUE;
1587          }
1588       }
1589
1590       delete [] posProbArray;
1591       delete [] negProbArray;
1592       return iResult;
1593
1594    } else {
1595       ///Not possible for AODs
1596       return kTRUE;
1597    }
1598
1599
1600
1601
1602 }
1603
1604
1605 ///________________________________________________________________________
1606 Bool_t AliConversionCuts::AcceptanceCut(TParticle *particle, TParticle * ePos,TParticle* eNeg){
1607    // MC Acceptance Cuts
1608    //(Certain areas were excluded for photon reconstruction)
1609
1610    if(particle->R()>fMaxR){
1611       return kFALSE;}
1612
1613    if(ePos->R()>fMaxR){
1614       return kFALSE;
1615    }
1616
1617    if(ePos->R()<fMinR){
1618       return kFALSE;
1619    }
1620
1621    if( ePos->R() <= ((abs(ePos->Vz())*fLineCutZRSlope)-fLineCutZValue)){
1622       return kFALSE;
1623    }
1624    else if (fUseEtaMinCut &&  ePos->R() >= ((abs(ePos->Vz())*fLineCutZRSlopeMin)-fLineCutZValueMin )){
1625       return kFALSE;
1626    }
1627
1628    if(abs(eNeg->Vz()) > fMaxZ ){ // cuts out regions where we do not reconstruct
1629       return kFALSE;
1630    }
1631
1632    if(eNeg->Vz()!=ePos->Vz()||eNeg->R()!=ePos->R()){
1633       return kFALSE;
1634    }
1635
1636    if(abs(ePos->Vz()) > fMaxZ ){ // cuts out regions where we do not reconstruct
1637       return kFALSE;
1638    }
1639
1640
1641    if( particle->Eta() > (fEtaCut + fEtaShift) || particle->Eta() < (-fEtaCut + fEtaShift) ){
1642       return kFALSE;
1643    }
1644    if( ePos->Eta() > (fEtaCut + fEtaShift) || ePos->Eta() < (-fEtaCut + fEtaShift) ){
1645       return kFALSE;
1646    }
1647    if( eNeg->Eta() > (fEtaCut + fEtaShift) || eNeg->Eta() < (-fEtaCut + fEtaShift) ){
1648       return kFALSE;
1649    }
1650    if(fEtaCutMin>-0.1){
1651       if( particle->Eta() < (fEtaCutMin + fEtaShift) && particle->Eta() > (-fEtaCutMin + fEtaShift) ){
1652          return kFALSE;
1653       }
1654       if( ePos->Eta() < (fEtaCutMin + fEtaShift) && ePos->Eta() > (-fEtaCutMin + fEtaShift) ){
1655          return kFALSE;
1656       }
1657       if( eNeg->Eta() < (fEtaCutMin + fEtaShift) && eNeg->Eta() > (-fEtaCutMin + fEtaShift) ){
1658          return kFALSE;
1659       }
1660    }
1661
1662    if( ePos->Pt()< fSinglePtCut ||  eNeg->Pt()< fSinglePtCut){
1663       return kFALSE;
1664    }
1665
1666    if(particle->Pt()<fPtCut){
1667       return kFALSE;
1668    }
1669
1670    return kTRUE;
1671 }
1672 ///________________________________________________________________________
1673 Bool_t AliConversionCuts::UpdateCutString() {
1674    ///Update the cut string (if it has been created yet)
1675
1676    if(fCutString && fCutString->GetString().Length() == kNCuts) {
1677       fCutString->SetString(GetCutNumber());
1678    } else {
1679       return kFALSE;
1680    }
1681    return kTRUE;
1682 }
1683 ///________________________________________________________________________
1684 void AliConversionCuts::LoadReweightingHistosMCFromFile() {
1685
1686   AliInfo("Entering loading of histograms for weighting");
1687   TFile *f = TFile::Open(fPathTrFReweighting.Data());
1688   if(!f){
1689      AliError(Form("file for weighting %s not found",fPathTrFReweighting.Data()));
1690      return;
1691   }
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()));
1697   }
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() ));
1703   }
1704
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() ));
1710   }
1711
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() ));
1717
1718   }
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() ));
1724   }
1725
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() ));
1731   }
1732
1733 }
1734
1735
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();
1742    }
1743
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));
1747       return kFALSE;
1748    }
1749    if(!analysisCutSelection.IsDigit()){
1750       AliError("Cut selection contains characters");
1751       return kFALSE;
1752    }
1753
1754    const char *cutSelection = analysisCutSelection.Data();
1755 #define ASSIGNARRAY(i)  fCuts[i] = cutSelection[i] - '0'
1756    for(Int_t ii=0;ii<kNCuts;ii++){
1757       ASSIGNARRAY(ii);
1758    }
1759
1760    // Set Individual Cuts
1761    for(Int_t ii=0;ii<kNCuts;ii++){
1762       if(!SetCut(cutIds(ii),fCuts[ii]))return kFALSE;
1763    }
1764
1765    //PrintCuts();
1766
1767    return kTRUE;
1768 }
1769 ///________________________________________________________________________
1770 Bool_t AliConversionCuts::SetCut(cutIds cutID, const Int_t value) {
1771    ///Set individual cut ID
1772
1773    switch (cutID) {
1774
1775    case kv0FinderType:
1776       if( SetV0Finder(value)) {
1777          fCuts[kv0FinderType] = value;
1778          UpdateCutString();
1779          return kTRUE;
1780       } else return kFALSE;
1781
1782    case kededxSigmaCut:
1783       if( SetTPCdEdxCutElectronLine(value)) {
1784          fCuts[kededxSigmaCut] = value;
1785          UpdateCutString();
1786          return kTRUE;
1787       } else return kFALSE;
1788
1789    case kpidedxSigmaCut:
1790       if( SetTPCdEdxCutPionLine(value)) {
1791          fCuts[kpidedxSigmaCut] = value;
1792          UpdateCutString();
1793          return kTRUE;
1794       } else return kFALSE;
1795
1796    case kpiMomdedxSigmaCut:
1797       if( SetMinMomPiondEdxCut(value)) {
1798          fCuts[kpiMomdedxSigmaCut] = value;
1799          UpdateCutString();
1800          return kTRUE;
1801       } else return kFALSE;
1802
1803    case kchi2GammaCut:
1804       if( SetChi2GammaCut(value)) {
1805          fCuts[kchi2GammaCut] = value;
1806          UpdateCutString();
1807          return kTRUE;
1808       } else return kFALSE;
1809
1810    case ksinglePtCut:
1811       if( SetSinglePtCut(value)) {
1812          fCuts[ksinglePtCut] = value;
1813          UpdateCutString();
1814          return kTRUE;
1815       } else return kFALSE;
1816
1817    case kclsTPCCut:
1818       if( SetTPCClusterCut(value)) {
1819          fCuts[kclsTPCCut] = value;
1820          UpdateCutString();
1821          return kTRUE;
1822       } else return kFALSE;
1823
1824    case ketaCut:
1825       if( SetEtaCut(value)) {
1826          fCuts[ketaCut] = value;
1827          UpdateCutString();
1828          return kTRUE;
1829       } else return kFALSE;
1830
1831    case kLowPRejectionSigmaCut:
1832       if( SetLowPRejectionCuts(value)) {
1833          fCuts[kLowPRejectionSigmaCut] = value;
1834          UpdateCutString();
1835          return kTRUE;
1836       } else return kFALSE;
1837
1838    case kQtMaxCut:
1839       if( SetQtMaxCut(value)) {
1840          fCuts[kQtMaxCut] = value;
1841          UpdateCutString();
1842          return kTRUE;
1843       } else return kFALSE;
1844
1845    case kpiMaxMomdedxSigmaCut:
1846       if( SetMaxMomPiondEdxCut(value)) {
1847          fCuts[kpiMaxMomdedxSigmaCut] = value;
1848          UpdateCutString();
1849          return kTRUE;
1850       } else return kFALSE;
1851
1852    case kRCut:
1853       if( SetRCut(value)) {
1854          fCuts[kRCut] = value;
1855          UpdateCutString();
1856          return kTRUE;
1857       } else return kFALSE;
1858
1859    case kremovePileUp:
1860       if( SetRemovePileUp(value)) {
1861          fCuts[kremovePileUp] = value;
1862          UpdateCutString();
1863          return kTRUE;
1864       } else return kFALSE;
1865
1866    case kselectV0AND:
1867       if( SetSelectSpecialTrigger(value)) {
1868          fCuts[kselectV0AND] = value;
1869          UpdateCutString();
1870          return kTRUE;
1871       } else return kFALSE;
1872
1873    case kmultiplicityMethod:
1874       if( SetMultiplicityMethod(value)) {
1875          fCuts[kmultiplicityMethod] = value;
1876          UpdateCutString();
1877          return kTRUE;
1878       } else return kFALSE;
1879
1880    case kisHeavyIon:
1881       if( SetIsHeavyIon(value)) {
1882          fCuts[kisHeavyIon] = value;
1883          UpdateCutString();
1884          return kTRUE;
1885       } else return kFALSE;
1886
1887    case kCentralityMin:
1888       if( SetCentralityMin(value)) {
1889          fCuts[kCentralityMin] = value;
1890          UpdateCutString();
1891          return kTRUE;
1892       } else return kFALSE;
1893
1894    case kCentralityMax:
1895       if( SetCentralityMax(value)) {
1896          fCuts[kCentralityMax] = value;
1897          UpdateCutString();
1898          return kTRUE;
1899       } else return kFALSE;
1900
1901    case kTOFelectronPID:
1902       if( SetTOFElectronPIDCut(value)) {
1903          fCuts[kTOFelectronPID] = value;
1904          UpdateCutString();
1905          return kTRUE;
1906       } else return kFALSE;
1907
1908    case kdoPhotonAsymmetryCut:
1909       if( SetPhotonAsymmetryCut(value)) {
1910          fCuts[kdoPhotonAsymmetryCut] = value;
1911          UpdateCutString();
1912          return kTRUE;
1913       } else return kFALSE;
1914
1915    case kPsiPair:
1916       if( SetPsiPairCut(value)) {
1917          fCuts[kPsiPair] = value;
1918          UpdateCutString();
1919          return kTRUE;
1920       } else return kFALSE;
1921
1922    case kCosPAngle:
1923       if( SetCosPAngleCut(value)) {
1924          fCuts[kCosPAngle] = value;
1925          UpdateCutString();
1926          return kTRUE;
1927       } else return kFALSE;
1928
1929
1930    case kElecShare:
1931       if( SetSharedElectronCut(value)) {
1932          fCuts[kElecShare] = value;
1933          UpdateCutString();
1934          return kTRUE;
1935       } else return kFALSE;
1936
1937    case kToCloseV0s:
1938       if( SetToCloseV0sCut(value)) {
1939          fCuts[kToCloseV0s] = value;
1940          UpdateCutString();
1941          return kTRUE;
1942       } else return kFALSE;
1943
1944    case kExtraSignals:
1945       if( SetRejectExtraSignalsCut(value)) {
1946          fCuts[kExtraSignals] = value;
1947          UpdateCutString();
1948          return kTRUE;
1949       } else return kFALSE;
1950
1951    case kDcaRPrimVtx:
1952       if( SetDCARPhotonPrimVtxCut(value)) {
1953          fCuts[kDcaRPrimVtx] = value;
1954          UpdateCutString();
1955          return kTRUE;
1956       } else return kFALSE;
1957
1958    case kDcaZPrimVtx:
1959       if( SetDCAZPhotonPrimVtxCut(value)) {
1960          fCuts[kDcaZPrimVtx] = value;
1961          UpdateCutString();
1962          return kTRUE;
1963       } else return kFALSE;
1964
1965
1966    case kNCuts:
1967       AliError("Cut id out of range");
1968       return kFALSE;
1969    }
1970
1971    AliError("Cut id %d not recognized");
1972    return kFALSE;
1973
1974
1975 }
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]);
1981    }
1982 }
1983 ///________________________________________________________________________
1984 Bool_t AliConversionCuts::SetIsHeavyIon(Int_t isHeavyIon)
1985 {   // Set Cut
1986    switch(isHeavyIon){
1987    case 0:
1988       fIsHeavyIon=0;
1989       break;
1990    case 1:
1991       fIsHeavyIon=1;
1992       fDetectorCentrality=0;
1993       break;
1994    case 2:
1995       fIsHeavyIon=1;
1996       fDetectorCentrality=1;
1997       break;
1998    case 3: //allows to select centrality 0-45% in steps of 5% for V0 Multiplicity
1999       fIsHeavyIon=1;
2000       fDetectorCentrality=0;
2001       fModCentralityClass=1;
2002       break;
2003    case 4: //allows to select centrality 45-90% in steps of 5% for V0 Multiplicity
2004       fIsHeavyIon=1;
2005       fDetectorCentrality=0;
2006       fModCentralityClass=2;
2007       break;
2008    case 5: //strict cut on v0 tracks for MC
2009       fIsHeavyIon=1;
2010       fDetectorCentrality=0;
2011       fModCentralityClass=3;
2012       break;
2013    case 6: //allows to select centrality 0-45% in steps of 5% for track mult
2014       //strict cut on v0 tracks for MC
2015       fIsHeavyIon=1;
2016       fDetectorCentrality=0;
2017       fModCentralityClass=4;
2018       break;
2019    case 7: //allows to select centrality 45-90% in steps of 5% for V0 Multiplicity
2020       //strict cut on v0 tracks for MC
2021       fIsHeavyIon=1;
2022       fDetectorCentrality=0;
2023       fModCentralityClass=5;
2024       break;
2025    case 8:
2026       fIsHeavyIon=2;
2027       fDetectorCentrality=0;
2028       break;
2029    case 9:
2030       fIsHeavyIon=2;
2031       fDetectorCentrality=1;
2032       break;
2033    default:
2034       AliError(Form("SetHeavyIon not defined %d",isHeavyIon));
2035       return kFALSE;
2036    }
2037    return kTRUE;
2038 }
2039
2040 //___________________________________________________________________
2041 Bool_t AliConversionCuts::SetCentralityMin(Int_t minCentrality)
2042 {
2043    // Set Cut
2044    if(minCentrality<0||minCentrality>9){
2045       AliError(Form("minCentrality not defined %d",minCentrality));
2046       return kFALSE;
2047    }
2048
2049    fCentralityMin=minCentrality;
2050    return kTRUE;
2051 }
2052 //___________________________________________________________________
2053 Bool_t AliConversionCuts::SetCentralityMax(Int_t maxCentrality)
2054 {
2055    // Set Cut
2056    if(maxCentrality<0||maxCentrality>9){
2057       AliError(Form("maxCentrality not defined %d",maxCentrality));
2058       return kFALSE;
2059    }
2060    fCentralityMax=maxCentrality;
2061    return kTRUE;
2062 }
2063 ///________________________________________________________________________
2064 Int_t AliConversionCuts::SetSelectSpecialTrigger(Int_t selectSpecialTrigger)
2065 {// Set Cut
2066
2067    switch(selectSpecialTrigger){
2068    case 0:
2069       fSpecialTrigger=0; // dont care
2070       break;
2071    case 1:
2072       fSpecialTrigger=1; // V0AND
2073       break;
2074    case 2:
2075       fSpecialTrigger=2; // with SDD requested
2076       break;
2077    case 3:
2078       fSpecialTrigger=3; // V0AND plus with SDD requested
2079       break;
2080    // allows to run MB & 6 other different trigger classes in parallel with the same photon cut
2081    case 4:
2082       fSpecialTrigger=4; // different trigger class as MB
2083       fTriggerSelectedManually = kTRUE;
2084       break;
2085    case 5:
2086       fSpecialTrigger=4; // different trigger class as MB
2087       fTriggerSelectedManually = kTRUE;
2088       break;
2089    case 6:
2090       fSpecialTrigger=4; // different trigger class as MB
2091       fTriggerSelectedManually = kTRUE;
2092       break;
2093    case 7:
2094       fSpecialTrigger=4; // different trigger class as MB
2095       fTriggerSelectedManually = kTRUE;
2096       break;
2097     case 8:
2098       fSpecialTrigger=4; // different trigger class as MB
2099       fTriggerSelectedManually = kTRUE;
2100       break;
2101     case 9:
2102       fSpecialTrigger=4; // different trigger class as MB
2103       fTriggerSelectedManually = kTRUE;
2104       break;
2105    default:
2106       AliError("Warning: Special Trigger Not known");
2107       return kFALSE;
2108    }
2109    return kTRUE;
2110 }
2111 ///________________________________________________________________________
2112 Bool_t AliConversionCuts::SetMultiplicityMethod(Int_t multiplicityMethod)
2113 {
2114    // Set Cut
2115    fMultiplicityMethod=multiplicityMethod;
2116
2117    // 0 Photon Multiplicity
2118    // 1 TPC Track multiplicity
2119    // 2 V0 Mult
2120    // 3 SPD Mult
2121
2122    return kTRUE;
2123 }
2124 ///________________________________________________________________________
2125 Bool_t AliConversionCuts::SetRemovePileUp(Int_t removePileUp)
2126 {// Set Cut
2127    switch(removePileUp){
2128    case 0:
2129       fRemovePileUp=kFALSE;
2130       break;
2131    case 1:
2132       fRemovePileUp=kTRUE;
2133       break;
2134    default:
2135       AliError("RemovePileUpCut not defined");
2136       return kFALSE;
2137    }
2138    return kTRUE;
2139 }
2140 ///________________________________________________________________________
2141 Bool_t AliConversionCuts::SetRejectExtraSignalsCut(Int_t extraSignal) {
2142
2143    switch(extraSignal){
2144    case 0:
2145       fRejectExtraSignals = 0;
2146       break; // No Rejection
2147    case 1:
2148       fRejectExtraSignals = 1;
2149       break; // MinBias Header
2150    case 2:
2151       fRejectExtraSignals = 2;
2152       break; // User String Array
2153    case 3:
2154       fRejectExtraSignals = 3;
2155       break; // Rejection for Gamma Correction only
2156    default:
2157       AliError(Form("Extra Signal Rejection not defined %d",extraSignal));
2158       return kFALSE;
2159    }
2160    return kTRUE;
2161 }
2162 ///________________________________________________________________________
2163 Bool_t AliConversionCuts::SetV0Finder(Int_t v0FinderType)
2164 {   // Set Cut
2165    switch (v0FinderType){
2166    case 0:  // on fly V0 finder
2167       cout << "have chosen onfly V0" << endl;
2168       fUseOnFlyV0Finder=kTRUE;
2169       break;
2170    case 1:  // offline V0 finder
2171       cout << "have chosen offline V0" << endl;
2172       fUseOnFlyV0Finder=kFALSE;
2173       break;
2174    default:
2175       AliError(Form(" v0FinderType not defined %d",v0FinderType));
2176       return kFALSE;
2177    }
2178    return kTRUE;
2179 }
2180 ///________________________________________________________________________
2181 Bool_t AliConversionCuts::SetEtaCut(Int_t etaCut)
2182 {   // Set Cut
2183
2184    //Set Standard LineCutZValues
2185    fLineCutZValueMin = -2;
2186    fLineCutZValue = 7.;
2187
2188    switch(etaCut){
2189    case 0: // 0.9
2190       fEtaCut     = 0.9;
2191       fLineCutZRSlope = tan(2*atan(exp(-fEtaCut)));
2192       fEtaCutMin     = -0.1;
2193       fLineCutZRSlopeMin = 0.;
2194       break;
2195    case 1:  // 1.2  // changed from 1.2 to 0.6 on 2013.06.10
2196       fEtaCut     = 0.6;
2197       fLineCutZRSlope = tan(2*atan(exp(-fEtaCut)));
2198       fEtaCutMin     = -0.1;
2199       fLineCutZRSlopeMin = 0.;
2200       break;
2201    case 2:  // 1.4
2202       fEtaCut     = 1.4;
2203       fLineCutZRSlope = tan(2*atan(exp(-fEtaCut)));
2204       fEtaCutMin     = -0.1;
2205       fLineCutZRSlopeMin = 0.;
2206       break;
2207    case 3: // 0.8
2208       fEtaCut     = 0.8;
2209       fLineCutZRSlope = tan(2*atan(exp(-fEtaCut)));
2210       fEtaCutMin     = -0.1;
2211       fLineCutZRSlopeMin = 0.;
2212       break;
2213    case 4: // 0.75
2214       fEtaCut     = 0.75;
2215       fLineCutZRSlope = tan(2*atan(exp(-fEtaCut)));
2216       fEtaCutMin     = -0.1;
2217       fLineCutZRSlopeMin = 0.;
2218       break;
2219    case 5: // 0.5
2220       fEtaCut     = 0.5;
2221       fLineCutZRSlope = tan(2*atan(exp(-fEtaCut)));
2222       fEtaCutMin     = -0.1;
2223       fLineCutZRSlopeMin = 0.;
2224       break;
2225    case 6: // 5.
2226       fEtaCut     = 5.;
2227       fLineCutZRSlope = tan(2*atan(exp(-fEtaCut)));
2228       fEtaCutMin     = -0.1;
2229       fLineCutZRSlopeMin = 0.;
2230       break;
2231    case 7:
2232       fEtaCut     = 0.3;
2233       fLineCutZRSlope = tan(2*atan(exp(-fEtaCut)));
2234       fEtaCutMin     = -0.1;
2235       fLineCutZRSlopeMin = 0.;
2236       break;
2237    // case 8: // 0.1 - 0.8
2238    //    fEtaCut     = 0.9;
2239    //    fLineCutZRSlope = tan(2*atan(exp(-fEtaCut)));
2240    //    fEtaCutMin     = 0.1;
2241    //    fLineCutZRSlopeMin = tan(2*atan(exp(-fEtaCutMin)));
2242    //    break;
2243    case 8: // 0.4
2244       fEtaCut     = 0.4;
2245       fLineCutZRSlope = tan(2*atan(exp(-fEtaCut)));
2246       fEtaCutMin     = -0.1;
2247       fLineCutZRSlopeMin = 0.;
2248       break;
2249    case 9: // 10
2250       fEtaCut     = 10;
2251       fLineCutZRSlope = tan(2*atan(exp(-fEtaCut)));
2252       fEtaCutMin     = -0.1;
2253       fLineCutZRSlopeMin = 0.;
2254       break;
2255    default:
2256       AliError(Form(" EtaCut not defined %d",etaCut));
2257       return kFALSE;
2258    }
2259    return kTRUE;
2260 }
2261 ///________________________________________________________________________
2262 Bool_t AliConversionCuts::SetRCut(Int_t RCut){
2263    // Set Cut
2264    switch(RCut){
2265    case 0:
2266       fMinR=0;
2267       fMaxR = 180.;
2268       break;
2269    case 1:
2270       fMinR=2.8;
2271       fMaxR = 180.;
2272       break;
2273    case 2:
2274       fMinR=5.;
2275       fMaxR = 180.;
2276       break;
2277    case 3:
2278       fMaxR = 70.;
2279       fMinR = 10.;
2280       break;
2281    case 4:
2282       fMaxR = 70.;
2283       fMinR = 5.;
2284       break;
2285       // High purity cuts for PbPb (remove first layers of material)
2286    case 5:
2287       fMaxR = 180.;
2288       fMinR = 10.;
2289       break;
2290    case 6:
2291       fMaxR = 180.;
2292       fMinR = 20.;
2293       break;
2294    case 7:
2295       fMaxR = 180.;
2296       fMinR = 26.;
2297       break;
2298    case 8:
2299       fMaxR = 180.;
2300       fMinR = 35.;
2301       break;
2302    case 9:
2303       fMaxR = 35.;
2304       fMinR = 5.;
2305       break;
2306
2307    default:
2308       AliError("RCut not defined");
2309       return kFALSE;
2310    }
2311    return kTRUE;
2312 }
2313 ///________________________________________________________________________
2314 Bool_t AliConversionCuts::SetSinglePtCut(Int_t singlePtCut)
2315 {   // Set Cut
2316    switch(singlePtCut){
2317    case 0: // 0.050 GeV
2318       fSinglePtCut = 0.050;
2319       break;
2320    case 1:  // 0.100 GeV
2321       fSinglePtCut = 0.100;
2322       break;
2323    case 2:  // 0.150 GeV
2324       fSinglePtCut = 0.150;
2325       break;
2326    case 3:  // 0.200 GeV
2327       fSinglePtCut = 0.200;
2328       break;
2329    case 4:  // 0.075 GeV
2330       fSinglePtCut = 0.075;
2331       break;
2332    case 5:  // 0.125 GeV
2333       fSinglePtCut = 0.125;
2334       break;
2335    case 6:  // 0.04 GeV
2336       fSinglePtCut = 0.040;
2337       break;
2338    case 7:  // 0.0 GeV
2339       fSinglePtCut = 0.0;
2340       break;
2341    default:
2342       AliError(Form("singlePtCut not defined %d",singlePtCut));
2343       return kFALSE;
2344    }
2345    return kTRUE;
2346 }
2347 ///________________________________________________________________________
2348 Bool_t AliConversionCuts::SetTPCClusterCut(Int_t clsTPCCut)
2349 {   // Set Cut
2350    switch(clsTPCCut){
2351    case 0: // 0
2352       fMinClsTPC= 0.;
2353       break;
2354    case 1:  // 60
2355       fMinClsTPC= 60.;
2356       break;
2357    case 2:  // 80
2358       fMinClsTPC= 80.;
2359       break;
2360    case 3:  // 100
2361       fMinClsTPC= 100.;
2362       break;
2363    case 4:  // 95% of findable clusters
2364       fMinClsTPCToF= 0.95;
2365       fUseCorrectedTPCClsInfo=1;
2366       break;
2367    case 5:  // 0% of findable clusters
2368       fMinClsTPCToF= 0.0;
2369       fUseCorrectedTPCClsInfo=1;
2370       break;
2371    case 6:  // 70% of findable clusters
2372       fMinClsTPCToF= 0.7;
2373       fUseCorrectedTPCClsInfo=1;
2374       break;
2375    case 7:  // 0% of findable clusters
2376       fMinClsTPCToF= 0.35;
2377       fUseCorrectedTPCClsInfo=0;
2378       break;
2379    case 8:
2380       fMinClsTPCToF= 0.35;
2381       fUseCorrectedTPCClsInfo=1;
2382       break;
2383    case 9:
2384       fMinClsTPCToF= 0.6;
2385       fUseCorrectedTPCClsInfo=1;
2386       break;
2387    default:
2388       AliError(Form("Warning: clsTPCCut not defined %d",clsTPCCut));
2389       return kFALSE;
2390    }
2391    return kTRUE;
2392 }
2393 ///________________________________________________________________________
2394 Bool_t AliConversionCuts::SetTPCdEdxCutElectronLine(Int_t ededxSigmaCut)
2395 {   // Set Cut
2396    switch(ededxSigmaCut){
2397    case 0: // -10,10
2398       fPIDnSigmaBelowElectronLine=-10;
2399       fPIDnSigmaAboveElectronLine=10;
2400       break;
2401    case 1: // -5,5
2402       fPIDnSigmaBelowElectronLine=-5;
2403       fPIDnSigmaAboveElectronLine=5;
2404       break;
2405    case 2: // -3,5
2406       fPIDnSigmaBelowElectronLine=-3;
2407       fPIDnSigmaAboveElectronLine=5;
2408       break;
2409    case 3: // -4,5
2410       fPIDnSigmaBelowElectronLine=-4;
2411       fPIDnSigmaAboveElectronLine=5;
2412       break;
2413    case 4: // -6,7
2414       fPIDnSigmaBelowElectronLine=-6;
2415       fPIDnSigmaAboveElectronLine=7;
2416       break;
2417    case 5: // -4,4
2418       fPIDnSigmaBelowElectronLine=-4;
2419       fPIDnSigmaAboveElectronLine=4;
2420       break;
2421    case 6: // -2.5,4
2422       fPIDnSigmaBelowElectronLine=-2.5;
2423       fPIDnSigmaAboveElectronLine=4;
2424       break;
2425    case 7: // -2,3.5
2426       fPIDnSigmaBelowElectronLine=-2;
2427       fPIDnSigmaAboveElectronLine=3.5;
2428       break;
2429    default:
2430       AliError("TPCdEdxCutElectronLine not defined");
2431       return kFALSE;
2432
2433    }
2434    return kTRUE;
2435 }
2436 ///________________________________________________________________________
2437 Bool_t AliConversionCuts::SetTPCdEdxCutPionLine(Int_t pidedxSigmaCut)
2438 {   // Set Cut
2439
2440    switch(pidedxSigmaCut){
2441    case 0:  // -10
2442       fPIDnSigmaAbovePionLine=-10;
2443       fPIDnSigmaAbovePionLineHighPt=-10;
2444       break;
2445    case 1:   // 0
2446       fPIDnSigmaAbovePionLine=0;
2447       fPIDnSigmaAbovePionLineHighPt=-10;
2448       break;
2449    case 2:  // 1
2450       fPIDnSigmaAbovePionLine=1;
2451       fPIDnSigmaAbovePionLineHighPt=-10;
2452       break;
2453    case 3:  // 1
2454       fPIDnSigmaAbovePionLine=2.5;
2455       fPIDnSigmaAbovePionLineHighPt=-10;
2456       break;
2457    case 4:  // 1
2458       fPIDnSigmaAbovePionLine=0.5;
2459       fPIDnSigmaAbovePionLineHighPt=-10;
2460       break;
2461    case 5:  // 1
2462       fPIDnSigmaAbovePionLine=2.;
2463       fPIDnSigmaAbovePionLineHighPt=-10;
2464       break;
2465    case 6:  // 1
2466       fPIDnSigmaAbovePionLine=2.;
2467       fPIDnSigmaAbovePionLineHighPt=0.5;
2468       break;
2469    case 7:  // 1
2470       fPIDnSigmaAbovePionLine=3.5;
2471       fPIDnSigmaAbovePionLineHighPt=-10;
2472       break;
2473    case 8:  // 1
2474       fPIDnSigmaAbovePionLine=2.;
2475       fPIDnSigmaAbovePionLineHighPt=1.;
2476       break;
2477    case 9:
2478       fPIDnSigmaAbovePionLine=3.0; // We need a bit less tight cut on dE/dx
2479       fPIDnSigmaAbovePionLineHighPt=-10;
2480       break;
2481    default:
2482       AliError(Form("Warning: pidedxSigmaCut not defined %d",pidedxSigmaCut));
2483       return kFALSE;
2484    }
2485    return kTRUE;
2486 }
2487 ///________________________________________________________________________
2488 Bool_t AliConversionCuts::SetMinMomPiondEdxCut(Int_t piMomdedxSigmaCut)
2489 {   // Set Cut
2490    switch(piMomdedxSigmaCut){
2491    case 0:  // 0.5 GeV
2492       fPIDMinPnSigmaAbovePionLine=0.5;
2493       break;
2494    case 1:  // 1. GeV
2495       fPIDMinPnSigmaAbovePionLine=1.;
2496       break;
2497    case 2:  // 1.5 GeV
2498       fPIDMinPnSigmaAbovePionLine=1.5;
2499       break;
2500    case 3:  // 20.0 GeV
2501       fPIDMinPnSigmaAbovePionLine=20.;
2502       break;
2503    case 4:  // 50.0 GeV
2504       fPIDMinPnSigmaAbovePionLine=50.;
2505       break;
2506    case 5:  // 0.3 GeV
2507       fPIDMinPnSigmaAbovePionLine=0.3;
2508       break;
2509    case 6:  // 0.25 GeV
2510       fPIDMinPnSigmaAbovePionLine=0.25;
2511       break;
2512    case 7:  // 0.4 GeV
2513       fPIDMinPnSigmaAbovePionLine=0.4;
2514       break;
2515    case 8:  // 0.2 GeV
2516       fPIDMinPnSigmaAbovePionLine=0.2;
2517       break;
2518    default:
2519       AliError(Form("piMomdedxSigmaCut not defined %d",piMomdedxSigmaCut));
2520       return kFALSE;
2521    }
2522    return kTRUE;
2523 }
2524 ///________________________________________________________________________
2525 Bool_t AliConversionCuts::SetMaxMomPiondEdxCut(Int_t piMaxMomdedxSigmaCut)
2526 {   // Set Cut
2527    switch(piMaxMomdedxSigmaCut){
2528    case 0:  // 100. GeV
2529       fPIDMaxPnSigmaAbovePionLine=100.;
2530       break;
2531    case 1:  // 5. GeV
2532       fPIDMaxPnSigmaAbovePionLine=5.;
2533       break;
2534    case 2:  // 4. GeV
2535       fPIDMaxPnSigmaAbovePionLine=4.;
2536       break;
2537    case 3:  // 3.5 GeV
2538       fPIDMaxPnSigmaAbovePionLine=3.5;
2539       break;
2540    case 4:  // 3. GeV
2541       fPIDMaxPnSigmaAbovePionLine=3.;
2542       break;
2543    case 5:  // 7. GeV
2544       fPIDMaxPnSigmaAbovePionLine=7.;
2545       break;
2546    default:
2547       AliError(Form("piMaxMomdedxSigmaCut not defined %d",piMaxMomdedxSigmaCut));
2548       return kFALSE;
2549    }
2550    return kTRUE;
2551 }
2552 ///________________________________________________________________________
2553 Bool_t AliConversionCuts::SetLowPRejectionCuts(Int_t LowPRejectionSigmaCut)
2554 {   // Set Cut
2555    switch(LowPRejectionSigmaCut){
2556    case 0:  //
2557       fPIDnSigmaAtLowPAroundKaonLine=0;
2558       fPIDnSigmaAtLowPAroundProtonLine=0;
2559       fPIDnSigmaAtLowPAroundPionLine=0;
2560       fDoKaonRejectionLowP = kFALSE;
2561       fDoProtonRejectionLowP = kFALSE;
2562       fDoPionRejectionLowP = kFALSE;
2563       fPIDMinPPionRejectionLowP = fPIDMinPnSigmaAbovePionLine;
2564       break;
2565    case 1:  //
2566       fPIDnSigmaAtLowPAroundKaonLine=0.5;
2567       fPIDnSigmaAtLowPAroundProtonLine=0.5;
2568       fPIDnSigmaAtLowPAroundPionLine=0.5;
2569       fDoKaonRejectionLowP = kTRUE;
2570       fDoProtonRejectionLowP = kTRUE;
2571       fDoPionRejectionLowP = kTRUE;
2572       fPIDMinPPionRejectionLowP = fPIDMinPnSigmaAbovePionLine;
2573       break;
2574    case 2:  //
2575       fPIDnSigmaAtLowPAroundKaonLine=1;
2576       fPIDnSigmaAtLowPAroundProtonLine=1;
2577       fPIDnSigmaAtLowPAroundPionLine=1;
2578       fDoKaonRejectionLowP = kTRUE;
2579       fDoProtonRejectionLowP = kTRUE;
2580       fDoPionRejectionLowP = kTRUE;
2581       fPIDMinPPionRejectionLowP = fPIDMinPnSigmaAbovePionLine;
2582       break;
2583    case 3:  //
2584       fPIDnSigmaAtLowPAroundKaonLine=2.;
2585       fPIDnSigmaAtLowPAroundProtonLine=2.;
2586       fPIDnSigmaAtLowPAroundPionLine=2.;
2587       fDoKaonRejectionLowP = kTRUE;
2588       fDoProtonRejectionLowP = kTRUE;
2589       fDoPionRejectionLowP = kTRUE;
2590       fPIDMinPPionRejectionLowP = fPIDMinPnSigmaAbovePionLine;
2591       break;
2592    case 4:  //
2593       fPIDnSigmaAtLowPAroundKaonLine=0.;
2594       fPIDnSigmaAtLowPAroundProtonLine=0.;
2595       fPIDnSigmaAtLowPAroundPionLine=1;
2596       fDoKaonRejectionLowP = kFALSE;
2597       fDoProtonRejectionLowP = kFALSE;
2598       fDoPionRejectionLowP = kTRUE;
2599       fPIDMinPPionRejectionLowP = fPIDMinPnSigmaAbovePionLine;
2600       break;
2601    case 5:  //
2602       fPIDnSigmaAtLowPAroundKaonLine=0.;
2603       fPIDnSigmaAtLowPAroundProtonLine=0.;
2604       fPIDnSigmaAtLowPAroundPionLine=1.5;
2605       fDoKaonRejectionLowP = kFALSE;
2606       fDoProtonRejectionLowP = kFALSE;
2607       fDoPionRejectionLowP = kTRUE;
2608       fPIDMinPPionRejectionLowP = fPIDMinPnSigmaAbovePionLine;
2609       break;
2610    case 6:  //
2611       fPIDnSigmaAtLowPAroundKaonLine=0.;
2612       fPIDnSigmaAtLowPAroundProtonLine=0.;
2613       fPIDnSigmaAtLowPAroundPionLine=2.;
2614       fDoKaonRejectionLowP = kFALSE;
2615       fDoProtonRejectionLowP = kFALSE;
2616       fDoPionRejectionLowP = kTRUE;
2617       fPIDMinPPionRejectionLowP = fPIDMinPnSigmaAbovePionLine;
2618       break;
2619    case 7:  //
2620       fPIDnSigmaAtLowPAroundKaonLine=0.;
2621       fPIDnSigmaAtLowPAroundProtonLine=0.;
2622       fPIDnSigmaAtLowPAroundPionLine=0.5;
2623       fDoKaonRejectionLowP = kFALSE;
2624       fDoProtonRejectionLowP = kFALSE;
2625       fDoPionRejectionLowP = kTRUE;
2626       fPIDMinPPionRejectionLowP = fPIDMinPnSigmaAbovePionLine;
2627       break;
2628    default:
2629       AliError(Form("LowPRejectionSigmaCut not defined %d",LowPRejectionSigmaCut));
2630       return kFALSE;
2631    }
2632    return kTRUE;
2633 }
2634 ///________________________________________________________________________
2635 Bool_t AliConversionCuts::SetTOFElectronPIDCut(Int_t TOFelectronPID){
2636    // Set Cut
2637    switch(TOFelectronPID){
2638    case 0: // no cut
2639       fUseTOFpid = kFALSE;
2640       fTofPIDnSigmaBelowElectronLine=-100;
2641       fTofPIDnSigmaAboveElectronLine=100;
2642       break;
2643    case 1: // -7,7
2644       fUseTOFpid = kTRUE;
2645       fTofPIDnSigmaBelowElectronLine=-7;
2646       fTofPIDnSigmaAboveElectronLine=7;
2647       break;
2648    case 2: // -5,5
2649       fUseTOFpid = kTRUE;
2650       fTofPIDnSigmaBelowElectronLine=-5;
2651       fTofPIDnSigmaAboveElectronLine=5;
2652       break;
2653    case 3: // -3,5
2654       fUseTOFpid = kTRUE;
2655       fTofPIDnSigmaBelowElectronLine=-3;
2656       fTofPIDnSigmaAboveElectronLine=5;
2657       break;
2658    case 4: // -2,3
2659       fUseTOFpid = kTRUE;
2660       fTofPIDnSigmaBelowElectronLine=-2;
2661       fTofPIDnSigmaAboveElectronLine=3;
2662       break;
2663    case 5: // -3,3
2664       fUseTOFpid = kTRUE;
2665       fTofPIDnSigmaBelowElectronLine=-3;
2666       fTofPIDnSigmaAboveElectronLine=3;
2667       break;
2668    default:
2669       AliError(Form("TOFElectronCut not defined %d",TOFelectronPID));
2670       return kFALSE;
2671    }
2672    return kTRUE;
2673 }
2674 ///________________________________________________________________________
2675 Bool_t AliConversionCuts::SetQtMaxCut(Int_t QtMaxCut)
2676 {   // Set Cut
2677    switch(QtMaxCut){
2678    case 0: //
2679       fQtMax=1.;
2680       fDoQtGammaSelection=kFALSE;
2681       fDoHighPtQtGammaSelection=kFALSE;
2682       fHighPtQtMax=100.;
2683       fPtBorderForQt=100.;
2684       break;
2685    case 1:
2686       fQtMax=0.1;
2687       fDoHighPtQtGammaSelection=kFALSE;
2688       fHighPtQtMax=100.;
2689       fPtBorderForQt=100.;
2690       break;
2691    case 2:
2692       fQtMax=0.07;
2693       fDoHighPtQtGammaSelection=kFALSE;
2694       fHighPtQtMax=100.;
2695       fPtBorderForQt=100.;
2696       break;
2697    case 3:
2698       fQtMax=0.05;
2699       fDoHighPtQtGammaSelection=kFALSE;
2700       fHighPtQtMax=100.;
2701       fPtBorderForQt=100.;
2702       break;
2703    case 4:
2704       fQtMax=0.03;
2705       fDoHighPtQtGammaSelection=kFALSE;
2706       fHighPtQtMax=100.;
2707       fPtBorderForQt=100.;
2708       break;
2709    case 5:
2710       fQtMax=0.02;
2711       fDoHighPtQtGammaSelection=kFALSE;
2712       fHighPtQtMax=100.;
2713       fPtBorderForQt=100.;
2714       break;
2715    case 6:
2716       fQtMax=0.02;
2717       fDoHighPtQtGammaSelection=kTRUE;
2718       fHighPtQtMax=0.06;
2719       fPtBorderForQt=2.5;
2720       break;
2721    case 7:
2722       fQtMax=0.15;
2723       fDoHighPtQtGammaSelection=kFALSE;
2724       fHighPtQtMax=100.;
2725       fPtBorderForQt=100.;
2726       break;
2727    default:
2728       AliError(Form("Warning: QtMaxCut not defined %d",QtMaxCut));
2729       return kFALSE;
2730    }
2731    return kTRUE;
2732 }
2733 ///________________________________________________________________________
2734 Bool_t AliConversionCuts::SetChi2GammaCut(Int_t chi2GammaCut)
2735 {   // Set Cut
2736
2737    switch(chi2GammaCut){
2738    case 0: // 100
2739       fChi2CutConversion = 100.;
2740       break;
2741    case 1:  // 50
2742       fChi2CutConversion = 50.;
2743       break;
2744    case 2:  // 30
2745       fChi2CutConversion = 30.;
2746       break;
2747    case 3:
2748       fChi2CutConversion = 200.;
2749       break;
2750    case 4:
2751       fChi2CutConversion = 500.;
2752       break;
2753    case 5:
2754       fChi2CutConversion = 100000.;
2755       break;
2756    case 6:
2757       fChi2CutConversion = 5.;
2758       break;
2759    case 7:
2760       fChi2CutConversion = 10.;
2761       break;
2762    case 8:
2763       fChi2CutConversion = 20.;
2764       break;
2765    case 9:
2766       fChi2CutConversion = 15.;
2767       break;
2768    default:
2769       AliError(Form("Warning: Chi2GammaCut not defined %d",chi2GammaCut));
2770       return kFALSE;
2771    }
2772    return kTRUE;
2773 }
2774 ///________________________________________________________________________
2775 Bool_t AliConversionCuts::SetPsiPairCut(Int_t psiCut) {
2776
2777    switch(psiCut) {
2778    case 0:
2779       fPsiPairCut = 10000; //
2780       break;
2781    case 1:
2782       fPsiPairCut = 0.1; //
2783       break;
2784    case 2:
2785       fPsiPairCut = 0.05; // Standard
2786       break;
2787    case 3:
2788       fPsiPairCut = 0.035; //
2789       break;
2790    case 4:
2791       fPsiPairCut = 0.15; //
2792       break;
2793    case 5:
2794       fPsiPairCut = 0.2; //
2795       break;
2796    case 6:
2797       fPsiPairCut = 0.03; //
2798       break;
2799    case 7:
2800       fPsiPairCut = 0.025; //
2801       break;
2802    case 8:
2803       fPsiPairCut = 0.01; //
2804       break;
2805    case 9:
2806       fPsiPairCut = 0.5; //
2807       break;
2808    default:
2809       AliError(Form("PsiPairCut not defined %d",psiCut));
2810       return kFALSE;
2811    }
2812
2813    return kTRUE;
2814 }
2815 ///________________________________________________________________________
2816 Bool_t AliConversionCuts::SetPhotonAsymmetryCut(Int_t doPhotonAsymmetryCut){
2817    // Set Cut
2818    switch(doPhotonAsymmetryCut){
2819    case 0:
2820       fDoPhotonAsymmetryCut=0;
2821       fMinPPhotonAsymmetryCut=100.;
2822       fMinPhotonAsymmetry=0.;
2823       break;
2824    case 1:
2825       fDoPhotonAsymmetryCut=1;
2826       fMinPPhotonAsymmetryCut=3.5;
2827       fMinPhotonAsymmetry=0.04;
2828       break;
2829    case 2:
2830       fDoPhotonAsymmetryCut=1;
2831       fMinPPhotonAsymmetryCut=3.5;
2832       fMinPhotonAsymmetry=0.06;
2833       break;
2834    default:
2835       AliError(Form("PhotonAsymmetryCut not defined %d",doPhotonAsymmetryCut));
2836       return kFALSE;
2837    }
2838    fCuts[kdoPhotonAsymmetryCut]=doPhotonAsymmetryCut;
2839    return kTRUE;
2840 }
2841 ///________________________________________________________________________
2842 Bool_t AliConversionCuts::SetCosPAngleCut(Int_t cosCut) {
2843
2844    switch(cosCut){
2845    case 0:
2846       fCosPAngleCut = TMath::Pi(); // -1
2847       break;
2848    case 1:
2849       fCosPAngleCut = 0.1; // 0.99500
2850       break;
2851    case 2:
2852       fCosPAngleCut = 0.05; // 0.99875
2853       break;
2854    case 3:
2855       fCosPAngleCut = 0.025; // 0.99969
2856       break;
2857    case 4:
2858       fCosPAngleCut = 0.01; // 0.99995
2859       break;
2860    case 5:
2861       fCosPAngleCut = 0.2; // 0.98007
2862       break;
2863    case 6:
2864       fCosPAngleCut = 0.5; // 0.87758
2865       break;
2866    case 7:
2867       fCosPAngleCut = 0.075; // 0.73169
2868       break;
2869    default:
2870       AliError(Form("Cosine Pointing Angle cut not defined %d",cosCut));
2871       return kFALSE;
2872    }
2873
2874    return kTRUE;
2875 }
2876 ///________________________________________________________________________
2877 Bool_t AliConversionCuts::SetSharedElectronCut(Int_t sharedElec) {
2878
2879    switch(sharedElec){
2880    case 0:
2881       fDoSharedElecCut = kFALSE;
2882       break;
2883    case 1:
2884       fDoSharedElecCut = kTRUE;
2885       break;
2886    default:
2887       AliError(Form("Shared Electron Cut not defined %d",sharedElec));
2888       return kFALSE;
2889    }
2890
2891    return kTRUE;
2892 }
2893 ///________________________________________________________________________
2894 Bool_t AliConversionCuts::SetToCloseV0sCut(Int_t toClose) {
2895
2896    switch(toClose){
2897    case 0:
2898       fDoToCloseV0sCut = kFALSE;
2899       fminV0Dist = 250;
2900       break;
2901    case 1:
2902       fDoToCloseV0sCut = kTRUE;
2903       fminV0Dist = 1;
2904       break;
2905    case 2:
2906       fDoToCloseV0sCut = kTRUE;
2907       fminV0Dist = 2;
2908       break;
2909    case 3:
2910       fDoToCloseV0sCut = kTRUE;
2911       fminV0Dist = 3;
2912       break;
2913    default:
2914       AliError(Form("Shared Electron Cut not defined %d",toClose));
2915       return kFALSE;
2916    }
2917    return kTRUE;
2918 }
2919 ///________________________________________________________________________
2920 Bool_t AliConversionCuts::SetTRDElectronCut(Int_t TRDElectronCut)
2921 {   // Set Cut
2922    switch(TRDElectronCut){
2923    case 0:
2924       fDoTRDPID=kFALSE;
2925       break;
2926    case 1:
2927       fDoTRDPID=kTRUE;
2928       fPIDTRDEfficiency=0.1;
2929       break;
2930    case 8:
2931       fDoTRDPID=kTRUE;
2932       fPIDTRDEfficiency=0.8;
2933       break;
2934    case 9:
2935       fDoTRDPID=kTRUE;
2936       fPIDTRDEfficiency=0.9;
2937       break;
2938    default:
2939       AliError(Form("TRDElectronCut not defined %d",TRDElectronCut));
2940       return kFALSE;
2941    }
2942
2943    return kTRUE;
2944 }
2945
2946 ///________________________________________________________________________
2947 Bool_t AliConversionCuts::SetDCAZPhotonPrimVtxCut(Int_t DCAZPhotonPrimVtx){
2948    // Set Cut
2949    switch(DCAZPhotonPrimVtx){
2950    case 0:  //
2951       fDCAZPrimVtxCut   = 1000;
2952       break;
2953    case 1:  //
2954       fDCAZPrimVtxCut   = 10;
2955       break;
2956    case 2:  //
2957       fDCAZPrimVtxCut   = 5;
2958       break;
2959    case 3:  //
2960       fDCAZPrimVtxCut   = 4;
2961       break;
2962    case 4:  //
2963       fDCAZPrimVtxCut   = 3;
2964       break;
2965    case 5:  //
2966       fDCAZPrimVtxCut   = 2.5;
2967       break;
2968    case 6:  //
2969       fDCAZPrimVtxCut   = 2;
2970       break;
2971    case 7:  //
2972       fDCAZPrimVtxCut   = 1.5;
2973       break;
2974    case 8:  //
2975       fDCAZPrimVtxCut   = 1;
2976       break;
2977    case 9:  //
2978       fDCAZPrimVtxCut   = 0.5;
2979       break;
2980    default:
2981       cout<<"Warning: DCAZPhotonPrimVtx not defined "<<DCAZPhotonPrimVtx<<endl;
2982       return kFALSE;
2983    }
2984    return kTRUE;
2985 }
2986
2987 ///________________________________________________________________________
2988 Bool_t AliConversionCuts::SetDCARPhotonPrimVtxCut(Int_t DCARPhotonPrimVtx){
2989    // Set Cut
2990    switch(DCARPhotonPrimVtx){
2991    case 0:  //
2992       fDCARPrimVtxCut   = 1000;
2993       break;
2994    case 1:  //
2995       fDCARPrimVtxCut   = 10;
2996       break;
2997    case 2:  //
2998       fDCARPrimVtxCut   = 5;
2999       break;
3000    case 3:  //
3001       fDCARPrimVtxCut   = 4;
3002       break;
3003    case 4:  //
3004       fDCARPrimVtxCut   = 3;
3005       break;
3006    case 5:  //
3007       fDCARPrimVtxCut   = 2.5;
3008       break;
3009    case 6:  //
3010       fDCARPrimVtxCut   = 2;
3011       break;
3012    case 7:  //
3013       fDCARPrimVtxCut   = 1.5;
3014       break;
3015    case 8:  //
3016       fDCARPrimVtxCut   = 1;
3017       break;
3018    case 9:  //
3019       fDCARPrimVtxCut   = 0.5;
3020       break;
3021    default:
3022       cout<<"Warning: DCARPhotonPrimVtx not defined "<<DCARPhotonPrimVtx<<endl;
3023       return kFALSE;
3024    }
3025    return kTRUE;
3026 }
3027
3028
3029 //-------------------------------------------------------------
3030 Double_t AliConversionCuts::GetCentrality(AliVEvent *event)
3031 {   // Get Event Centrality
3032
3033    AliESDEvent *esdEvent=dynamic_cast<AliESDEvent*>(event);
3034    if(esdEvent){
3035       AliCentrality *fESDCentrality=(AliCentrality*)esdEvent->GetCentrality();
3036
3037       if(fDetectorCentrality==0){
3038          if (fIsHeavyIon==2){
3039             return fESDCentrality->GetCentralityPercentile("V0A"); // default for pPb
3040          } else{
3041             return fESDCentrality->GetCentralityPercentile("V0M"); // default
3042          }
3043       }
3044       if(fDetectorCentrality==1){
3045          return fESDCentrality->GetCentralityPercentile("CL1");
3046       }
3047    }
3048
3049    AliAODEvent *aodEvent=dynamic_cast<AliAODEvent*>(event);
3050    if(aodEvent){
3051       if(aodEvent->GetHeader()){return aodEvent->GetHeader()->GetCentrality();}
3052    }
3053
3054    return -1;
3055 }
3056 //-------------------------------------------------------------
3057 Bool_t AliConversionCuts::IsCentralitySelected(AliVEvent *event, AliVEvent *fMCEvent)
3058 {   // Centrality Selection
3059    if(!fIsHeavyIon)return kTRUE;
3060
3061    if(fCentralityMin == fCentralityMax ) return kTRUE;//0-100%
3062    else if(fCentralityMax==0) fCentralityMax=10; //CentralityRange = fCentralityMin-100%
3063
3064    Double_t centrality=GetCentrality(event);
3065    if(centrality<0)return kFALSE;
3066
3067    Int_t centralityC=0;
3068    if (fModCentralityClass == 0){
3069       centralityC= Int_t(centrality/10);
3070       if(centralityC >= fCentralityMin && centralityC < fCentralityMax)
3071          return kTRUE;
3072       else return kFALSE;
3073    }
3074    else if (fModCentralityClass ==1){
3075       centralityC= Int_t(centrality);
3076       if(centralityC >= fCentralityMin*5 && centralityC < fCentralityMax*5){
3077          return kTRUE;
3078       } else return kFALSE;
3079    }
3080    else if (fModCentralityClass ==2){
3081       centralityC= Int_t(centrality);
3082       if(centralityC >= ((fCentralityMin*5)+45) && centralityC < ((fCentralityMax*5)+45))
3083          return kTRUE;
3084       else return kFALSE;
3085    }
3086
3087    Int_t nprimaryTracks = ((AliV0ReaderV1*)AliAnalysisManager::GetAnalysisManager()->GetTask("V0ReaderV1"))->GetNumberOfPrimaryTracks();
3088    Int_t PrimaryTracks10[10][2] =
3089       {
3090          {9999,9999}, //  0
3091          {1210,2150}, // 10
3092          { 817,1435}, // 20
3093          { 536, 930}, // 30
3094          { 337, 570}, // 40
3095          { 197, 327}, // 50
3096          { 106, 173}, // 60
3097          {  51,  81}, // 70
3098          {  21,  34}, // 80
3099          {   0,   0}  // 90
3100       };
3101    Int_t PrimaryTracks5a[10][2] =
3102       {
3103          {9999,9999}, // 0
3104          {1485,2640}, // 5
3105          {1210,2150}, // 10
3106          { 995,1760}, // 15
3107          { 817,1435}, // 20
3108          { 666,1160}, // 25
3109          { 536, 930}, // 30
3110          { 428, 731}, // 35
3111          { 337, 570}, // 40
3112          { 260, 436}  // 45
3113       };
3114    Int_t PrimaryTracks5b[10][2] =
3115       {
3116          { 260, 436}, // 45
3117          { 197, 327}, // 50
3118          { 147, 239}, // 55
3119          { 106, 173}, // 60
3120          {  75, 120}, // 65
3121          {  51,  81}, // 70
3122          {  34,  53}, // 75
3123          {  21,  34}, // 80
3124          {  13,  19}, // 85
3125          {   0,   0}  // 90
3126       };
3127
3128    Int_t column = -1;
3129    if(event->IsA()==AliESDEvent::Class()) column = 0;
3130    if(event->IsA()==AliAODEvent::Class()) column = 1;
3131
3132    if (fModCentralityClass == 3){
3133       if(fMCEvent){
3134          if(nprimaryTracks > PrimaryTracks10[fCentralityMax][column] && nprimaryTracks <= PrimaryTracks10[fCentralityMin][column])
3135             return kTRUE;
3136          else return kFALSE;
3137       }
3138       else{
3139          centralityC= Int_t(centrality/10);
3140          if(centralityC >= fCentralityMin && centralityC < fCentralityMax)
3141             return kTRUE;
3142          else return kFALSE;
3143       }
3144    }
3145    else if (fModCentralityClass ==4){
3146       if(fMCEvent){
3147          if(nprimaryTracks > PrimaryTracks5a[fCentralityMax][column] && nprimaryTracks <= PrimaryTracks5a[fCentralityMin][column])
3148             return kTRUE;
3149          else return kFALSE;
3150       }
3151       else{
3152          centralityC= Int_t(centrality);
3153          if(centralityC >= fCentralityMin*5 && centralityC < fCentralityMax*5){
3154             return kTRUE;
3155          } else return kFALSE;
3156       }
3157    }
3158    else if (fModCentralityClass ==5){
3159       if(fMCEvent){
3160          if(nprimaryTracks > PrimaryTracks5b[fCentralityMax][column] && nprimaryTracks <= PrimaryTracks5b[fCentralityMin][column])
3161             return kTRUE;
3162          else return kFALSE;
3163       }
3164       else{
3165          centralityC= Int_t(centrality);
3166          if(centralityC >= ((fCentralityMin*5)+45) && centralityC < ((fCentralityMax*5)+45))
3167             return kTRUE;
3168          else return kFALSE;
3169       }
3170    }
3171
3172    return kFALSE;
3173 }
3174 ///________________________________________________________________________
3175 Bool_t AliConversionCuts::VertexZCut(AliVEvent *event){
3176    // Cut on z position of primary vertex
3177    Double_t fVertexZ=event->GetPrimaryVertex()->GetZ();
3178
3179    if(abs(fVertexZ)>fMaxVertexZ)return kFALSE;
3180
3181    if (fIsHeavyIon == 2){
3182      if(fUtils->IsFirstEventInChunk(event)) return kFALSE;
3183      if(!fUtils->IsVertexSelected2013pA(event)) return kFALSE;
3184
3185    }
3186
3187    return kTRUE;
3188 }
3189 ///________________________________________________________________________
3190
3191 Int_t AliConversionCuts::GetNumberOfContributorsVtx(AliVEvent *event){
3192    // returns number of contributors to the vertex
3193
3194    AliESDEvent *fESDEvent=dynamic_cast<AliESDEvent*>(event);
3195    if(fESDEvent){
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();
3200          }
3201       }
3202
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();
3207          }  else {
3208             AliWarning(Form("Number of contributors from bad vertex type:: %s",fESDEvent->GetPrimaryVertex()->GetName()));
3209 //            cout << "rejected " << fESDEvent->GetEventNumberInFile() << endl;
3210             return 0;
3211          }
3212       }
3213    }
3214
3215    AliAODEvent *fAODEvent=dynamic_cast<AliAODEvent*>(event);
3216    if(fAODEvent){
3217       if (fAODEvent->GetPrimaryVertex() != NULL){
3218          if(fAODEvent->GetPrimaryVertex()->GetNContributors()>0) {
3219             return fAODEvent->GetPrimaryVertex()->GetNContributors();
3220          }
3221       }
3222       if(fAODEvent->GetPrimaryVertexSPD() !=NULL){
3223          if(fAODEvent->GetPrimaryVertexSPD()->GetNContributors()>0) {
3224             return fAODEvent->GetPrimaryVertexSPD()->GetNContributors();
3225          } else {
3226             AliWarning(Form("Number of contributors from bad vertex type:: %s",fAODEvent->GetPrimaryVertex()->GetName()));
3227             return 0;
3228          }
3229       }
3230    }
3231   // cout << "rejected " << fESDEvent->GetEventNumberInFile() << endl;
3232    return 0;
3233 }
3234
3235 ///________________________________________________________________________
3236
3237 Bool_t AliConversionCuts::IsTriggerSelected(AliVEvent *fInputEvent)
3238 {
3239
3240    AliInputEventHandler *fInputHandler=(AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
3241
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;
3247          else {
3248             if (fIsHeavyIon == 1) fOfflineTriggerMask = AliVEvent::kMB | AliVEvent::kCentral | AliVEvent::kSemiCentral;
3249                else if (fIsHeavyIon == 2) fOfflineTriggerMask = AliVEvent::kINT7;
3250                else fOfflineTriggerMask = AliVEvent::kMB;
3251          }
3252       }
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;
3257
3258       if (fOfflineTriggerMask)
3259          isSelected = fOfflineTriggerMask & fInputHandler->IsEventSelected();
3260    }
3261    fIsSDDFired = !(fInputHandler->IsEventSelected() & AliVEvent::kFastOnly);
3262
3263    // Fill Histogram
3264    if(hTriggerClass){
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);
3307    }
3308
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);
3351    }
3352
3353    if(!isSelected)return kFALSE;
3354
3355    return kTRUE;
3356
3357 }
3358
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;
3369    Int_t nClsI = 63;
3370    Int_t nClsIO = 127;
3371
3372    if(radius <= radiusI){
3373       return firstTPCRow;
3374    }
3375    if(radius>radiusI && radius<=radiusO){
3376       firstTPCRow = (Int_t)((radius-radiusI)/rSizeI);
3377    }
3378    if(radius>radiusO && radius<=radiusOB){
3379       firstTPCRow = (Int_t)(nClsI+(radius-radiusO)/rSizeO);
3380    }
3381
3382    if(radius>radiusOB){
3383       firstTPCRow =(Int_t)(nClsIO+(radius-radiusOB)/rSizeOB);
3384    }
3385
3386    return firstTPCRow;
3387 }
3388
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))){
3392       return kFALSE;
3393    }
3394    return kTRUE;
3395 }
3396
3397 Double_t AliConversionCuts::GetCosineOfPointingAngle( const AliConversionPhotonBase * photon, AliVEvent * event) const{
3398    // calculates the pointing angle of the recalculated V0
3399
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]);
3407    }
3408    if(event->IsA()==AliAODEvent::Class()){
3409       momV0[0] = photon->GetPx();
3410       momV0[1] = photon->GetPy();
3411       momV0[2] = photon->GetPz();
3412    }
3413
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
3418
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];
3421
3422
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);
3426
3427    return cosinePointingAngle;
3428 }
3429
3430 ///________________________________________________________________________
3431 Bool_t AliConversionCuts::PsiPairCut(const AliConversionPhotonBase * photon) const {
3432
3433    if(photon->GetPsiPair() > fPsiPairCut){
3434       return kFALSE;}
3435    else{return kTRUE;}
3436
3437 }
3438
3439 ///________________________________________________________________________
3440 TString AliConversionCuts::GetCutNumber(){
3441    // returns TString with current cut number
3442    TString a(kNCuts);
3443    for(Int_t ii=0;ii<kNCuts;ii++){
3444       a.Append(Form("%d",fCuts[ii]));
3445    }
3446    return a;
3447 }
3448
3449 ///________________________________________________________________________
3450 void AliConversionCuts::FillElectonLabelArray(AliAODConversionPhoton* photon, Int_t nV0){
3451
3452    Int_t posLabel = photon->GetTrackLabelPositive();
3453    Int_t negLabel = photon->GetTrackLabelNegative();
3454
3455    fElectronLabelArray[nV0*2] = posLabel;
3456    fElectronLabelArray[(nV0*2)+1] = negLabel;
3457 }
3458 ///________________________________________________________________________
3459 Bool_t AliConversionCuts::RejectSharedElectronV0s(AliAODConversionPhoton* photon, Int_t nV0, Int_t nV0s){
3460
3461    Int_t posLabel = photon->GetTrackLabelPositive();
3462    Int_t negLabel = photon->GetTrackLabelNegative();
3463
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){
3468          return kFALSE;}
3469       if(fElectronLabelArray[i] == negLabel){
3470          return kFALSE;}
3471    }
3472
3473    return kTRUE;
3474 }
3475 ///________________________________________________________________________
3476 Bool_t AliConversionCuts::RejectToCloseV0s(AliAODConversionPhoton* photon, TList *photons, Int_t nV0){
3477
3478
3479    Double_t posX = photon->GetConversionX();
3480    Double_t posY = photon->GetConversionY();
3481    Double_t posZ = photon->GetConversionZ();
3482
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();
3489
3490       Double_t dist = pow((posX - posCompX),2)+pow((posY - posCompY),2)+pow((posZ - posCompZ),2);
3491
3492       if(dist < fminV0Dist*fminV0Dist){
3493          if(photon->GetChi2perNDF() < photonComp->GetChi2perNDF()) return kTRUE;
3494          else {
3495             return kFALSE;}
3496       }
3497
3498    }
3499    return kTRUE;
3500 }
3501 ///________________________________________________________________________
3502 void AliConversionCuts::GetNotRejectedParticles(Int_t rejection, TList *HeaderList, AliVEvent *MCEvent){
3503
3504
3505
3506    if(fNotRejectedStart){
3507       delete[] fNotRejectedStart;
3508       fNotRejectedStart = NULL;
3509    }
3510    if(fNotRejectedEnd){
3511       delete[] fNotRejectedEnd;
3512       fNotRejectedEnd = NULL;
3513    }
3514    if(fGeneratorNames){
3515       delete[] fGeneratorNames;
3516       fGeneratorNames = NULL;
3517    }
3518
3519    if(rejection == 0) return; // No Rejection
3520
3521    AliGenCocktailEventHeader *cHeader = 0x0;
3522    AliAODMCHeader *cHeaderAOD = 0x0;
3523    Bool_t headerFound = kFALSE;
3524
3525    if(MCEvent->IsA()==AliMCEvent::Class()){
3526       cHeader = dynamic_cast<AliGenCocktailEventHeader*>(dynamic_cast<AliMCEvent*>(MCEvent)->GenEventHeader());
3527       if(cHeader) headerFound = kTRUE;
3528    }
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;
3532    }
3533
3534    if(headerFound){
3535       TList *genHeaders = 0x0;
3536       if(cHeader) genHeaders = cHeader->GetHeaders();
3537       if(cHeaderAOD){
3538          genHeaders = cHeaderAOD->GetCocktailHeaders();
3539          if(genHeaders->GetEntries()==1){
3540             SetRejectExtraSignalsCut(0);
3541             return;
3542          }
3543       }
3544       AliGenEventHeader* gh = 0;
3545       fnHeaders = 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){
3554                   fnHeaders++;
3555                   continue;
3556                }
3557             }
3558          }
3559       }
3560
3561       fNotRejectedStart = new Int_t[fnHeaders];
3562       fNotRejectedEnd = new Int_t[fnHeaders];
3563       fGeneratorNames = new TString[fnHeaders];
3564
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();
3569          return;
3570       }
3571
3572       Int_t firstindex = 0;
3573       Int_t lastindex =  -1;
3574       Int_t nummer = 0;
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;
3586                nummer++;
3587                continue;
3588             }
3589          }
3590          firstindex = firstindex + gh->NProduced();
3591       }
3592    } else { // No Cocktail Header Found
3593       fNotRejectedStart = new Int_t[1];
3594       fNotRejectedEnd = new Int_t[1];
3595
3596       fnHeaders = 1;
3597       fNotRejectedStart[0] = 0;
3598       fNotRejectedEnd[0] = static_cast<AliMCEvent*>(MCEvent)->Stack()->GetNprimary()-1;
3599       fGeneratorNames = new TString[1];
3600       fGeneratorNames[0] = "NoCocktailGeneratorFound";
3601
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";
3608
3609       SetRejectExtraSignalsCut(0);
3610    }
3611
3612 }
3613 //_________________________________________________________________________
3614 Int_t AliConversionCuts::IsParticleFromBGEvent(Int_t index, AliStack *MCStack, AliVEvent *InputEvent){
3615
3616    // Not Accepted == kFALSE == 0
3617    //     Accepted ==  kTRUE == 1
3618    //  FirstHeader ==  kTRUE == 3
3619    if(index < 0) return 0; // No Particle
3620
3621    Int_t accepted = 0;
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);
3626       }
3627       for(Int_t i = 0;i<fnHeaders;i++){
3628          if(index >= fNotRejectedStart[i] && index <= fNotRejectedEnd[i]){
3629             accepted = 1;
3630             if(i == 0) accepted = 2; // MB Header
3631          }
3632       }
3633    }
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);
3641       }
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]){
3645             accepted = 1;
3646             if(i == 0) accepted = 2; // MB Header
3647          }
3648       }
3649    }
3650
3651    return accepted;
3652 }
3653 //_________________________________________________________________________
3654 Int_t AliConversionCuts::IsEventAcceptedByConversionCut(AliConversionCuts *ReaderCuts, AliVEvent *InputEvent, AliMCEvent *MCEvent, Bool_t isHeavyIon){
3655
3656    if ( !IsTriggerSelected(InputEvent) )
3657       return 3;
3658
3659    if(isHeavyIon && !(IsCentralitySelected(InputEvent,MCEvent)))
3660       return 1; // Check Centrality --> Not Accepted => eventQuality = 1
3661
3662    if(!isHeavyIon && GetIsFromPileup()){
3663       if(InputEvent->IsPileupFromSPD(3,0.8,3.,2.,5.)){
3664
3665          return 6; // Check Pileup --> Not Accepted => eventQuality = 6
3666       }
3667    }
3668
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
3673
3674    if( (IsSpecialTrigger() == 1 || IsSpecialTrigger() == 3) && !hasV0And)
3675       return 8; // V0AND requested but no fired
3676
3677
3678    return 0;
3679 }
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.;
3683
3684    Int_t kCaseGen = 0;
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){
3688             kCaseGen = 1;
3689          } else if (fGeneratorNames[i].CompareTo("DPMJET") == 0){
3690             kCaseGen = 2;
3691          } else if (fGeneratorNames[i].CompareTo("HIJING") == 0 ||
3692                     fGeneratorNames[i].CompareTo("Hijing") == 0 ||
3693                     fGeneratorNames[i].Contains("hijing")){
3694             kCaseGen = 3;
3695          } else if (fGeneratorNames[i].CompareTo("BOX") == 0){
3696              kCaseGen = 4;
3697          } else if (fGeneratorNames[i].CompareTo("PARAM") == 0){
3698             kCaseGen = 5;
3699          } else if (fGeneratorNames[i].CompareTo("NoCocktailGeneratorFound") == 0){
3700             kCaseGen = 6;
3701          } else if (fGeneratorNames[i].CompareTo("NoCocktailGeneratorFound_Pythia") == 0){
3702             kCaseGen = 1;
3703          } else if (fGeneratorNames[i].CompareTo("NoCocktailGeneratorFound_Phojet") == 0){
3704             kCaseGen = 2;
3705          } else if (fGeneratorNames[i].CompareTo("NoCocktailGeneratorFound_Hijing") == 0){
3706             kCaseGen = 3;
3707          }
3708          TString periodName = ((AliV0ReaderV1*)AliAnalysisManager::GetAnalysisManager()
3709                                                 ->GetTask("V0ReaderV1"))->GetPeriodName();
3710
3711          if (periodName.Contains("LHC13d2")){
3712             kCaseGen = 3;
3713          }
3714       }
3715    }
3716    if (kCaseGen == 0) return 1;
3717
3718
3719    Double_t mesonPt = 0;
3720    Double_t mesonMass = 0;
3721    Int_t PDGCode = 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();
3726    }
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();
3733    }
3734
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){
3741          dNdyMC = 2.1462;
3742          nMC = 7.06055;
3743          tMC = 0.12533;
3744       } else if ( PDGCode ==  221){
3745          dNdyMC = 0.2357;
3746          nMC = 5.9105;
3747          tMC = 0.1525;
3748       }
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){
3755          dNdyMC = 2.35978;
3756          nMC = 6.81795;
3757          tMC = 0.11492;
3758       } else if ( PDGCode ==  221){
3759          dNdyMC = 0.3690;
3760          nMC = 5.55809;
3761          tMC = 0.13387;
3762       }
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;
3767       Float_t b = 5.6661;
3768       Float_t c = -1430.5863;
3769       Float_t d = -0.6966624;
3770       Float_t e = 252.3742;
3771       if ( PDGCode ==  111){
3772          a = 0.23437;
3773          b = 5.6661;
3774          c = -1430.5863;
3775          d = -0.6966624;
3776          e = 252.3742;
3777       } else if ( PDGCode ==  221){
3778          a = 0.10399;
3779          b = 4.35311;
3780          c = -12.17723;
3781          d = -0.01172;
3782          e =1.85140;
3783       }
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);
3789       }
3790       if ( PDGCode ==  221 && fDoReweightHistoMCEta && hReweightMCHistEta!= 0x0){
3791          functionResultMC = hReweightMCHistEta->Interpolate(mesonPt);
3792       }
3793       if ( PDGCode ==  310 && fDoReweightHistoMCK0s && hReweightMCHistK0s!= 0x0){
3794          functionResultMC = hReweightMCHistK0s->Interpolate(mesonPt);
3795       }
3796    }
3797
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){
3804          dNdyData = 2.2328;
3805          nData = 7.1473;
3806          tData = 0.1346;
3807       } else if ( PDGCode ==  221){
3808          dNdyData = 0.38992; //be careful this fit is not optimal, eta in data still has problems
3809          nData = 5.72778;
3810          tData = 0.13835;
3811       }
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;
3814    } else {
3815       if ( PDGCode ==  111 && fDoReweightHistoMCPi0 && fFitDataPi0!= 0x0){
3816          functionResultData = fFitDataPi0->Eval(mesonPt);
3817       }
3818       if ( PDGCode ==  221 && fDoReweightHistoMCEta && fFitDataEta!= 0x0){
3819          functionResultData = fFitDataEta->Eval(mesonPt);
3820       }
3821       if ( PDGCode ==  310 && fDoReweightHistoMCK0s && fFitDataK0s!= 0x0){
3822          functionResultData = fFitDataK0s->Eval(mesonPt);
3823       }
3824
3825 //       Float_t a = 0.;
3826 //       Float_t b = 0.;
3827 //       Float_t c = 0.;
3828 //       Float_t d = 0.;
3829 //       Float_t e = 0.;
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;
3897 //          }
3898 //
3899 //          functionResultData = a*TMath::Power(mesonPt,-1*(b+c/(TMath::Power(mesonPt,d)+e)));
3900 //       }
3901
3902    }
3903
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)){
3910                weight = 1.;
3911             }
3912          }
3913          if (!isfinite(functionResultData)) weight = 1.;
3914          if (!isfinite(weight)) weight = 1.;
3915       }
3916    } else if (PDGCode ==  310 && functionResultMC != 0 && isfinite(functionResultMC)){
3917         weight = functionResultMC;
3918    }
3919
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;
3922 //    }
3923    return weight;
3924 }
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;}
3931     return cuts;
3932 }
3933
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;}
3940     return cuts;
3941 }
3942 ///________________________________________________________________________
3943 void AliConversionCuts::GetCorrectEtaShiftFromPeriod(TString periodName){
3944
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