]> git.uio.no Git - u/mrichter/AliRoot.git/blob - PWGGA/GammaConv/AliConversionCuts.cxx
updates to fix histogram ranges, new 2D cuts, new Analysis macro for fast QA
[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    "EvetPlane" //27
90 };
91
92
93 //________________________________________________________________________
94 AliConversionCuts::AliConversionCuts(const char *name,const char *title) :
95    AliAnalysisCuts(name,title),
96    fHistograms(NULL),
97    fHeaderList(NULL),
98    fPIDResponse(NULL),
99    fEventQuality(-1),
100    fMaxR(200),
101    fMinR(0),
102    fEtaCut(0.9),
103    fEtaCutMin(-0.1),
104    fPtCut(0.02),
105    fSinglePtCut(0),
106    fMaxZ(1000),
107    fMinClsTPC(0.),
108    fMinClsTPCToF(0.),
109    fLineCutZRSlope(0.),
110    fLineCutZValue(0),
111    fLineCutZRSlopeMin(0.),
112    fLineCutZValueMin(0),
113    fChi2CutConversion(1000),
114    fPIDProbabilityCutNegativeParticle(0),
115    fPIDProbabilityCutPositiveParticle(0),
116    fDodEdxSigmaCut(kTRUE),
117    fDoTOFsigmaCut(kFALSE),
118    fPIDTRDEfficiency(1),
119    fDoTRDPID(kFALSE),
120    fPIDnSigmaAboveElectronLine(100),
121    fPIDnSigmaBelowElectronLine(-100),
122    fTofPIDnSigmaAboveElectronLine(100),
123    fTofPIDnSigmaBelowElectronLine(-100),
124    fPIDnSigmaAbovePionLine(0),
125    fPIDnSigmaAbovePionLineHighPt(-100),
126    fPIDMinPnSigmaAbovePionLine(0),
127    fPIDMaxPnSigmaAbovePionLine(0),
128    fDoKaonRejectionLowP(kFALSE),
129    fDoProtonRejectionLowP(kFALSE),
130    fDoPionRejectionLowP(kFALSE),
131    fPIDnSigmaAtLowPAroundKaonLine(0),
132    fPIDnSigmaAtLowPAroundProtonLine(0),
133    fPIDnSigmaAtLowPAroundPionLine(0),
134    fPIDMinPKaonRejectionLowP(1.5),
135    fPIDMinPProtonRejectionLowP(2),
136    fPIDMinPPionRejectionLowP(0),
137    fDoQtGammaSelection(kTRUE),
138    fDo2DQt(kFALSE),
139    fQtMax(100),
140    fXVertexCut(0.),
141    fYVertexCut(0.),
142    fZVertexCut(0.),
143    fNSigmaMass(0.),
144    fUseEtaMinCut(kFALSE),
145    fUseOnFlyV0Finder(kTRUE),
146    fDoPhotonAsymmetryCut(kTRUE),
147    fMinPPhotonAsymmetryCut(100.),
148    fMinPhotonAsymmetry(0.),
149    fIsHeavyIon(0),
150    fDetectorCentrality(0),
151    fModCentralityClass(0),
152    fMaxVertexZ(10),
153    fCentralityMin(0),
154    fCentralityMax(0),
155    fUseCorrectedTPCClsInfo(kFALSE),
156    fUseTOFpid(kFALSE),
157    fMultiplicityMethod(0),
158    fSpecialTrigger(0),
159    fRemovePileUp(kFALSE),
160    fOpeningAngle(0.005),
161    fPsiPairCut(10000),
162    fDo2DPsiPairChi2(kFALSE),
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    fInPlaneOutOfPlane(0),
177    fConversionPointXArray(0.0),
178    fConversionPointYArray(0.0),
179    fConversionPointZArray(0.0),
180    fnHeaders(0),
181    fNotRejectedStart(NULL),
182    fNotRejectedEnd(NULL),
183    fGeneratorNames(NULL),
184    fCutString(NULL),
185    fUtils(NULL),
186    fEtaShift(0.0),
187    fDoEtaShift(kFALSE),
188    fDoReweightHistoMCPi0(kFALSE),
189    fDoReweightHistoMCEta(kFALSE),
190    fDoReweightHistoMCK0s(kFALSE),
191    fPathTrFReweighting(""),
192    fNameHistoReweightingPi0(""),
193    fNameHistoReweightingEta(""),
194    fNameHistoReweightingK0s(""),
195    fNameFitDataPi0(""),
196    fNameFitDataEta(""),
197    fNameFitDataK0s(""),
198    hdEdxCuts(NULL),
199    hTPCdEdxbefore(NULL),
200    hTPCdEdxafter(NULL),
201    hTPCdEdxSigbefore(NULL),
202    hTPCdEdxSigafter(NULL),
203    hTOFbefore(NULL),
204    hTOFSigbefore(NULL),
205    hTOFSigafter(NULL),
206    hPsiPairDeltaPhiafter(NULL),
207    hTrackCuts(NULL),
208    hPhotonCuts(NULL),
209    hInvMassbefore(NULL),
210    hArmenterosbefore(NULL),
211    hInvMassafter(NULL),
212    hArmenterosafter(NULL),
213    hAcceptanceCuts(NULL),
214    hCutIndex(NULL),
215    hV0EventCuts(NULL),
216    hCentrality(NULL),
217    hCentralityVsNumberOfPrimaryTracks(NULL),
218    hVertexZ(NULL),
219    hEventPlanePhi(NULL),
220    hTriggerClass(NULL),
221    hTriggerClassSelected(NULL),
222    hReweightMCHistPi0(NULL),
223    hReweightMCHistEta(NULL),
224    hReweightMCHistK0s(NULL),
225    fFitDataPi0(NULL),
226    fFitDataEta(NULL),
227    fFitDataK0s(NULL),
228    fPreSelCut(kFALSE),
229    fTriggerSelectedManually(kFALSE),
230    fSpecialTriggerName("")
231
232 {
233    InitPIDResponse();
234    for(Int_t jj=0;jj<kNCuts;jj++){fCuts[jj]=0;}
235    fCutString=new TObjString((GetCutNumber()).Data());
236
237    fElectronLabelArray = new Int_t[fElectronArraySize];
238    fUtils = new AliAnalysisUtils();
239    //if you do not want to apply the cut on the distance between the SPD and TRK vertex:
240    //fUtils->SetCutOnZVertexSPD(kFALSE);
241
242
243 }
244
245 //________________________________________________________________________
246 AliConversionCuts::AliConversionCuts(const AliConversionCuts &ref) :
247    AliAnalysisCuts(ref),
248    fHistograms(NULL),
249    fHeaderList(ref.fHeaderList),
250    fPIDResponse(NULL),
251    fEventQuality(ref.fEventQuality),
252    fMaxR(ref.fMaxR),
253    fMinR(ref.fMinR),
254    fEtaCut(ref.fEtaCut),
255    fEtaCutMin(ref.fEtaCutMin),
256    fPtCut(ref.fPtCut),
257    fSinglePtCut(ref.fSinglePtCut),
258    fMaxZ(ref.fMaxZ),
259    fMinClsTPC(ref.fMinClsTPC),
260    fMinClsTPCToF(ref.fMinClsTPCToF),
261    fLineCutZRSlope(ref.fLineCutZRSlope),
262    fLineCutZValue(ref.fLineCutZValue),
263    fLineCutZRSlopeMin(ref.fLineCutZRSlopeMin),
264    fLineCutZValueMin(ref.fLineCutZValueMin),
265    fChi2CutConversion(ref.fChi2CutConversion),
266    fPIDProbabilityCutNegativeParticle(ref.fPIDProbabilityCutNegativeParticle),
267    fPIDProbabilityCutPositiveParticle(ref.fPIDProbabilityCutPositiveParticle),
268    fDodEdxSigmaCut(ref. fDodEdxSigmaCut),
269    fDoTOFsigmaCut(ref.fDoTOFsigmaCut),
270    fPIDTRDEfficiency(ref.fPIDTRDEfficiency),
271    fDoTRDPID(ref.fDoTRDPID),
272    fPIDnSigmaAboveElectronLine(ref.fPIDnSigmaAboveElectronLine),
273    fPIDnSigmaBelowElectronLine(ref.fPIDnSigmaBelowElectronLine),
274    fTofPIDnSigmaAboveElectronLine(ref.fTofPIDnSigmaAboveElectronLine),
275    fTofPIDnSigmaBelowElectronLine(ref.fTofPIDnSigmaBelowElectronLine),
276    fPIDnSigmaAbovePionLine(ref.fPIDnSigmaAbovePionLine),
277    fPIDnSigmaAbovePionLineHighPt(ref.fPIDnSigmaAbovePionLineHighPt),
278    fPIDMinPnSigmaAbovePionLine(ref.fPIDMinPnSigmaAbovePionLine),
279    fPIDMaxPnSigmaAbovePionLine(ref.fPIDMaxPnSigmaAbovePionLine),
280    fDoKaonRejectionLowP(ref.fDoKaonRejectionLowP),
281    fDoProtonRejectionLowP(ref.fDoProtonRejectionLowP),
282    fDoPionRejectionLowP(ref.fDoPionRejectionLowP),
283    fPIDnSigmaAtLowPAroundKaonLine(ref.fPIDnSigmaAtLowPAroundKaonLine),
284    fPIDnSigmaAtLowPAroundProtonLine(ref.fPIDnSigmaAtLowPAroundProtonLine),
285    fPIDnSigmaAtLowPAroundPionLine(ref.fPIDnSigmaAtLowPAroundPionLine),
286    fPIDMinPKaonRejectionLowP(ref.fPIDMinPKaonRejectionLowP),
287    fPIDMinPProtonRejectionLowP(ref.fPIDMinPProtonRejectionLowP),
288    fPIDMinPPionRejectionLowP(ref.fPIDMinPPionRejectionLowP),
289    fDoQtGammaSelection(ref.fDoQtGammaSelection),
290    fDo2DQt(ref.fDo2DQt),
291    fQtMax(ref.fQtMax),
292    fXVertexCut(ref.fXVertexCut),
293    fYVertexCut(ref.fYVertexCut),
294    fZVertexCut(ref.fZVertexCut),
295    fNSigmaMass(ref.fNSigmaMass),
296    fUseEtaMinCut(ref.fUseEtaMinCut),
297    fUseOnFlyV0Finder(ref.fUseOnFlyV0Finder),
298    fDoPhotonAsymmetryCut(ref.fDoPhotonAsymmetryCut),
299    fMinPPhotonAsymmetryCut(ref.fMinPPhotonAsymmetryCut),
300    fMinPhotonAsymmetry(ref.fMinPhotonAsymmetry),
301    fIsHeavyIon(ref.fIsHeavyIon),
302    fDetectorCentrality(ref.fDetectorCentrality),
303    fModCentralityClass(ref.fModCentralityClass),
304    fMaxVertexZ(ref.fMaxVertexZ),
305    fCentralityMin(ref.fCentralityMin),
306    fCentralityMax(ref.fCentralityMax),
307    fUseCorrectedTPCClsInfo(ref.fUseCorrectedTPCClsInfo),
308    fUseTOFpid(ref.fUseTOFpid),
309    fMultiplicityMethod(ref.fMultiplicityMethod),
310    fSpecialTrigger(ref.fSpecialTrigger),
311    fRemovePileUp(ref.fRemovePileUp),
312    fOpeningAngle(ref.fOpeningAngle),
313    fPsiPairCut(ref.fPsiPairCut),
314    fDo2DPsiPairChi2(ref.fDo2DPsiPairChi2),
315    fCosPAngleCut(ref.fCosPAngleCut),
316    fDoToCloseV0sCut(ref.fDoToCloseV0sCut),
317    fRejectExtraSignals(ref.fRejectExtraSignals),
318    fminV0Dist(ref.fminV0Dist),
319    fDoSharedElecCut(ref.fDoSharedElecCut),
320    fOfflineTriggerMask(ref.fOfflineTriggerMask),
321    fHasV0AND(ref.fHasV0AND),
322    fIsSDDFired(ref.fIsSDDFired),
323    fRandom(ref.fRandom),
324    fElectronArraySize(ref.fElectronArraySize),
325    fElectronLabelArray(NULL),
326    fDCAZPrimVtxCut(ref.fDCAZPrimVtxCut),
327    fDCARPrimVtxCut(ref.fDCAZPrimVtxCut),
328    fInPlaneOutOfPlane(ref.fInPlaneOutOfPlane),
329    fConversionPointXArray(ref.fConversionPointXArray),
330    fConversionPointYArray(ref.fConversionPointYArray),
331    fConversionPointZArray(ref.fConversionPointZArray),
332    fnHeaders(ref.fnHeaders),
333    fNotRejectedStart(NULL),
334    fNotRejectedEnd(NULL),
335    fGeneratorNames(ref.fGeneratorNames),
336    fCutString(NULL),
337    fUtils(NULL),
338    fEtaShift(ref.fEtaShift),
339    fDoEtaShift(ref.fDoEtaShift),
340    fDoReweightHistoMCPi0(ref.fDoReweightHistoMCPi0),
341    fDoReweightHistoMCEta(ref.fDoReweightHistoMCEta),
342    fDoReweightHistoMCK0s(ref.fDoReweightHistoMCK0s),
343    fPathTrFReweighting(ref.fPathTrFReweighting),
344    fNameHistoReweightingPi0(ref.fNameHistoReweightingPi0),
345    fNameHistoReweightingEta(ref.fNameHistoReweightingEta),
346    fNameHistoReweightingK0s(ref.fNameHistoReweightingK0s),
347    fNameFitDataPi0(ref.fNameFitDataPi0),
348    fNameFitDataEta(ref.fNameFitDataEta),
349    fNameFitDataK0s(ref.fNameFitDataK0s),
350    hdEdxCuts(NULL),
351    hTPCdEdxbefore(NULL),
352    hTPCdEdxafter(NULL),
353    hTPCdEdxSigbefore(NULL),
354    hTPCdEdxSigafter(NULL),
355    hTOFbefore(NULL),
356    hTOFSigbefore(NULL),
357    hTOFSigafter(NULL),
358    hPsiPairDeltaPhiafter(NULL),
359    hTrackCuts(NULL),
360    hPhotonCuts(NULL),
361    hInvMassbefore(NULL),
362    hArmenterosbefore(NULL),
363    hInvMassafter(NULL),
364    hArmenterosafter(NULL),
365    hAcceptanceCuts(NULL),
366    hCutIndex(NULL),
367    hV0EventCuts(NULL),
368    hCentrality(NULL),
369    hCentralityVsNumberOfPrimaryTracks(NULL),
370    hVertexZ(NULL),
371    hEventPlanePhi(NULL),
372    hTriggerClass(NULL),
373    hTriggerClassSelected(NULL),
374    hReweightMCHistPi0(ref.hReweightMCHistPi0),
375    hReweightMCHistEta(ref.hReweightMCHistEta),
376    hReweightMCHistK0s(ref.hReweightMCHistK0s),
377    fFitDataPi0(ref.fFitDataPi0),
378    fFitDataEta(ref.fFitDataEta),
379    fFitDataK0s(ref.fFitDataK0s),
380    fPreSelCut(ref.fPreSelCut),
381    fTriggerSelectedManually(ref.fTriggerSelectedManually),
382    fSpecialTriggerName(ref.fSpecialTriggerName)
383 {
384    // Copy Constructor
385    for(Int_t jj=0;jj<kNCuts;jj++){fCuts[jj]=ref.fCuts[jj];}
386    fCutString=new TObjString((GetCutNumber()).Data());
387    fElectronLabelArray = new Int_t[fElectronArraySize];
388    fUtils = new AliAnalysisUtils();
389    // dont copy histograms (if you like histograms, call InitCutHistograms())
390
391 }
392
393
394 //________________________________________________________________________
395 AliConversionCuts::~AliConversionCuts() {
396    // Destructor
397    //Deleting fHistograms leads to seg fault it it's added to output collection of a task
398    // if(fHistograms)
399    //    delete fHistograms;
400    // fHistograms = NULL;
401    if(fCutString != NULL){
402       delete fCutString;
403       fCutString = NULL;
404    }
405    if(fElectronLabelArray){
406       delete fElectronLabelArray;
407       fElectronLabelArray = NULL;
408    }
409    if(fNotRejectedStart){
410       delete[] fNotRejectedStart;
411       fNotRejectedStart = NULL;
412    }
413    if(fNotRejectedEnd){
414       delete[] fNotRejectedEnd;
415       fNotRejectedEnd = NULL;
416    }
417    if(fGeneratorNames){
418       delete[] fGeneratorNames;
419       fGeneratorNames = NULL;
420    }
421    if(fUtils){
422      delete fUtils;
423      fUtils = NULL;
424    }
425
426 }
427
428 //________________________________________________________________________
429 void AliConversionCuts::InitCutHistograms(TString name, Bool_t preCut){
430
431    // Initialize Cut Histograms for QA (only initialized and filled if function is called)
432    TH1::AddDirectory(kFALSE);
433
434    if(fHistograms != NULL){
435       delete fHistograms;
436       fHistograms=NULL;
437    }
438    if(fHistograms==NULL){
439       fHistograms=new TList();
440       fHistograms->SetOwner(kTRUE);
441       if(name=="")fHistograms->SetName(Form("ConvCuts_%s",GetCutNumber().Data()));
442       else fHistograms->SetName(Form("%s_%s",name.Data(),GetCutNumber().Data()));
443    }
444
445    if (hReweightMCHistPi0){
446       hReweightMCHistPi0->SetName("MCInputForWeightingPi0");
447       fHistograms->Add(hReweightMCHistPi0);
448    }
449    if (hReweightMCHistEta){
450       hReweightMCHistEta->SetName("MCInputForWeightingEta");
451       fHistograms->Add(hReweightMCHistEta);
452    }
453    if (hReweightMCHistK0s){
454       hReweightMCHistK0s->SetName("MCInputForWeightingK0s");
455       fHistograms->Add(hReweightMCHistK0s);
456    }
457 //    if (fFitDataPi0){
458 //       fFitDataPi0->SetName("DataFitForWeightingPi0");
459 //       fHistograms->Add(fFitDataPi0);
460 //    }
461 //    if (fFitDataEta){
462 //       fFitDataEta->SetName("DataFitForWeightingEta");
463 //       fHistograms->Add(fFitDataEta);
464 //    }
465 //    if (fFitDataK0s){
466 //       fFitDataK0s->SetName("DataFitForWeightingK0s");
467 //       fHistograms->Add(fFitDataK0s);
468 //    }
469    // IsPhotonSelected
470    hCutIndex=new TH1F(Form("IsPhotonSelected %s",GetCutNumber().Data()),"IsPhotonSelected",11,-0.5,10.5);
471    hCutIndex->GetXaxis()->SetBinLabel(kPhotonIn+1,"in");
472    hCutIndex->GetXaxis()->SetBinLabel(kOnFly+1,"onfly");
473    hCutIndex->GetXaxis()->SetBinLabel(kNoTracks+1,"no tracks");
474    hCutIndex->GetXaxis()->SetBinLabel(kdEdxCuts+1,"dEdx");
475    hCutIndex->GetXaxis()->SetBinLabel(kTrackCuts+1,"Track cuts");
476    hCutIndex->GetXaxis()->SetBinLabel(kConvPointFail+1,"ConvPoint fail");
477    hCutIndex->GetXaxis()->SetBinLabel(kPhotonCuts+1,"PhotonCuts");
478    hCutIndex->GetXaxis()->SetBinLabel(kEventPlane+1,"EventPlane");
479    hCutIndex->GetXaxis()->SetBinLabel(kPhotonOut+1,"out");
480    fHistograms->Add(hCutIndex);
481
482    // Track Cuts
483    hTrackCuts=new TH1F(Form("TrackCuts %s",GetCutNumber().Data()),"TrackCuts",9,-0.5,8.5);
484    hTrackCuts->GetXaxis()->SetBinLabel(1,"in");
485    hTrackCuts->GetXaxis()->SetBinLabel(2,"likesign");
486    hTrackCuts->GetXaxis()->SetBinLabel(3,"ntpccl");
487    hTrackCuts->GetXaxis()->SetBinLabel(4,"acceptance");
488    hTrackCuts->GetXaxis()->SetBinLabel(5,"singlept");
489    hTrackCuts->GetXaxis()->SetBinLabel(6,"TPCrefit");
490    hTrackCuts->GetXaxis()->SetBinLabel(7,"kink");
491    hTrackCuts->GetXaxis()->SetBinLabel(8,"out");
492    fHistograms->Add(hTrackCuts);
493
494    // Photon Cuts
495    hPhotonCuts=new TH1F(Form("PhotonCuts %s",GetCutNumber().Data()),"PhotonCuts",14,-0.5,13.5);
496    hPhotonCuts->GetXaxis()->SetBinLabel(1,"in");
497    hPhotonCuts->GetXaxis()->SetBinLabel(2,"qtcut");
498    hPhotonCuts->GetXaxis()->SetBinLabel(3,"chi2");
499    hPhotonCuts->GetXaxis()->SetBinLabel(4,"acceptance");
500    hPhotonCuts->GetXaxis()->SetBinLabel(5,"asymmetry");
501    hPhotonCuts->GetXaxis()->SetBinLabel(6,"pidprob");
502    hPhotonCuts->GetXaxis()->SetBinLabel(7,"cortpcclinfo");
503    hPhotonCuts->GetXaxis()->SetBinLabel(8,"PsiPair");
504    hPhotonCuts->GetXaxis()->SetBinLabel(9,"CosPAngle");
505    hPhotonCuts->GetXaxis()->SetBinLabel(10,"DCA R");
506    hPhotonCuts->GetXaxis()->SetBinLabel(11,"DCA Z");
507    hPhotonCuts->GetXaxis()->SetBinLabel(12,"out");
508    fHistograms->Add(hPhotonCuts);
509
510    if(preCut){
511       hInvMassbefore=new TH1F(Form("InvMass_before %s",GetCutNumber().Data()),"InvMass_before",1000,0,0.3);
512       fHistograms->Add(hInvMassbefore);
513       hArmenterosbefore=new TH2F(Form("Armenteros_before %s",GetCutNumber().Data()),"Armenteros_before",200,-1,1,1000,0,1.);
514       fHistograms->Add(hArmenterosbefore);
515    }
516    hInvMassafter=new TH1F(Form("InvMass_after %s",GetCutNumber().Data()),"InvMass_after",1000,0,0.3);
517    fHistograms->Add(hInvMassafter);
518    hArmenterosafter=new TH2F(Form("Armenteros_after %s",GetCutNumber().Data()),"Armenteros_after",200,-1,1,250,0,0.25);
519    fHistograms->Add(hArmenterosafter);
520
521    hAcceptanceCuts=new TH1F(Form("PhotonAcceptanceCuts %s",GetCutNumber().Data()),"PhotonAcceptanceCuts",10,-0.5,9.5);
522    hAcceptanceCuts->GetXaxis()->SetBinLabel(1,"in");
523    hAcceptanceCuts->GetXaxis()->SetBinLabel(2,"maxR");
524    hAcceptanceCuts->GetXaxis()->SetBinLabel(3,"minR");
525    hAcceptanceCuts->GetXaxis()->SetBinLabel(4,"line");
526    hAcceptanceCuts->GetXaxis()->SetBinLabel(5,"maxZ");
527    hAcceptanceCuts->GetXaxis()->SetBinLabel(6,"eta");
528    hAcceptanceCuts->GetXaxis()->SetBinLabel(7,"minpt");
529    hAcceptanceCuts->GetXaxis()->SetBinLabel(8,"out");
530    fHistograms->Add(hAcceptanceCuts);
531
532    // dEdx Cuts
533    hdEdxCuts=new TH1F(Form("dEdxCuts %s",GetCutNumber().Data()),"dEdxCuts",10,-0.5,9.5);
534    hdEdxCuts->GetXaxis()->SetBinLabel(1,"in");
535    hdEdxCuts->GetXaxis()->SetBinLabel(2,"TPCelectron");
536    hdEdxCuts->GetXaxis()->SetBinLabel(3,"TPCpion");
537    hdEdxCuts->GetXaxis()->SetBinLabel(4,"TPCpionhighp");
538    hdEdxCuts->GetXaxis()->SetBinLabel(5,"TPCkaonlowprej");
539    hdEdxCuts->GetXaxis()->SetBinLabel(6,"TPCprotonlowprej");
540    hdEdxCuts->GetXaxis()->SetBinLabel(7,"TPCpionlowprej");
541    hdEdxCuts->GetXaxis()->SetBinLabel(8,"TOFelectron");
542    hdEdxCuts->GetXaxis()->SetBinLabel(9,"TRDelectron");
543    hdEdxCuts->GetXaxis()->SetBinLabel(10,"out");
544    fHistograms->Add(hdEdxCuts);
545
546    TAxis *AxisBeforedEdx = NULL;
547    TAxis *AxisBeforedEdxSig = NULL;
548    TAxis *AxisBeforeTOF = NULL;
549    TAxis *AxisBeforeTOFSig = NULL;
550    if(preCut){
551       hTPCdEdxbefore=new TH2F(Form("Gamma_dEdx_before %s",GetCutNumber().Data()),"dEdx Gamma before" ,150,0.03,20,800,0,200);
552       fHistograms->Add(hTPCdEdxbefore);
553       AxisBeforedEdx = hTPCdEdxbefore->GetXaxis();
554       hTPCdEdxSigbefore=new TH2F(Form("Gamma_dEdxSig_before %s",GetCutNumber().Data()),"dEdx Sigma Gamma before" ,150,0.03,20,400,-10,10);
555       fHistograms->Add(hTPCdEdxSigbefore);
556       AxisBeforedEdxSig = hTPCdEdxSigbefore->GetXaxis();
557
558       hTOFbefore=new TH2F(Form("Gamma_TOF_before %s",GetCutNumber().Data()),"TOF Gamma before" ,150,0.03,20,11000,-1000,10000);
559       fHistograms->Add(hTOFbefore);
560       AxisBeforeTOF = hTOFbefore->GetXaxis();
561       hTOFSigbefore=new TH2F(Form("Gamma_TOFSig_before %s",GetCutNumber().Data()),"TOF Sigma Gamma before" ,150,0.03,20,400,-6,10);
562       fHistograms->Add(hTOFSigbefore);
563       AxisBeforeTOFSig = hTOFSigbefore->GetXaxis();
564
565    }
566    hTPCdEdxSigafter=new TH2F(Form("Gamma_dEdxSig_after %s",GetCutNumber().Data()),"dEdx Sigma Gamma after" ,150,0.03,20,400, -10,10);
567    fHistograms->Add(hTPCdEdxSigafter);
568
569    hTPCdEdxafter=new TH2F(Form("Gamma_dEdx_after %s",GetCutNumber().Data()),"dEdx Gamma after" ,150,0.03,20,800,0,200);
570    fHistograms->Add(hTPCdEdxafter);
571
572    hTOFSigafter=new TH2F(Form("Gamma_TOFSig_after %s",GetCutNumber().Data()),"TOF Sigma Gamma after" ,150,0.03,20,400,-6,10);
573    fHistograms->Add(hTOFSigafter);
574
575    hPsiPairDeltaPhiafter=new TH2F(Form("Gamma_PsiPairDeltaPhi_after %s",GetCutNumber().Data()),"Psi Pair vs Delta Phi Gamma after" ,200,-2,2,200,-2,2);
576    fHistograms->Add(hPsiPairDeltaPhiafter);
577
578    TAxis *AxisAfter = hTPCdEdxSigafter->GetXaxis();
579    Int_t bins = AxisAfter->GetNbins();
580    Double_t from = AxisAfter->GetXmin();
581    Double_t to = AxisAfter->GetXmax();
582    Double_t *newBins = new Double_t[bins+1];
583    newBins[0] = from;
584    Double_t factor = TMath::Power(to/from, 1./bins);
585    for(Int_t i=1; i<=bins; ++i) newBins[i] = factor * newBins[i-1];
586    AxisAfter->Set(bins, newBins);
587    AxisAfter = hTOFSigafter->GetXaxis();
588    AxisAfter->Set(bins, newBins);
589    AxisAfter = hTPCdEdxafter->GetXaxis();
590    AxisAfter->Set(bins, newBins);
591    if(preCut){
592       AxisBeforedEdx->Set(bins, newBins);
593       AxisBeforeTOF->Set(bins, newBins);
594       AxisBeforedEdxSig->Set(bins, newBins);
595       AxisBeforeTOFSig->Set(bins, newBins);
596    }
597    delete [] newBins;
598
599    hCentrality=new TH1F(Form("Centrality %s",GetCutNumber().Data()),"Centrality",1000,0,100);
600    fHistograms->Add(hCentrality);
601    hCentralityVsNumberOfPrimaryTracks=new TH2F(Form("Centrality vs Primary Tracks %s",GetCutNumber().Data()),"Centrality vs Primary Tracks ",100,0,100,4000,0,4000);
602    fHistograms->Add(hCentralityVsNumberOfPrimaryTracks);
603
604    // Event Cuts and Info
605    if(preCut){
606       hV0EventCuts=new TH1F(Form("ESD_EventCuts %s",GetCutNumber().Data()),"Event Cuts",7,-0.5,6.5);
607       hV0EventCuts->GetXaxis()->SetBinLabel(1,"in");
608       hV0EventCuts->GetXaxis()->SetBinLabel(2,"OfflineTrigger");
609       hV0EventCuts->GetXaxis()->SetBinLabel(3,"nvtxcontr");
610       hV0EventCuts->GetXaxis()->SetBinLabel(4,"VertexZ");
611       hV0EventCuts->GetXaxis()->SetBinLabel(5,"pileup");
612       hV0EventCuts->GetXaxis()->SetBinLabel(6,"centrsel");
613       hV0EventCuts->GetXaxis()->SetBinLabel(7,"out");
614       fHistograms->Add(hV0EventCuts);
615
616       hVertexZ=new TH1F(Form("VertexZ %s",GetCutNumber().Data()),"VertexZ",1000,-50,50);
617       fHistograms->Add(hVertexZ);
618
619       hTriggerClass= new TH1F(Form("OfflineTrigger %s",GetCutNumber().Data()),"OfflineTrigger",35,-0.5,34.5);
620       hTriggerClass->GetXaxis()->SetBinLabel( 1,"kMB");
621       hTriggerClass->GetXaxis()->SetBinLabel( 2,"kINT7");
622       hTriggerClass->GetXaxis()->SetBinLabel( 3,"kMUON");
623       hTriggerClass->GetXaxis()->SetBinLabel( 4,"kHighMult");
624       hTriggerClass->GetXaxis()->SetBinLabel( 5,"kKEMC1");
625       hTriggerClass->GetXaxis()->SetBinLabel( 6,"kCINT5");
626       hTriggerClass->GetXaxis()->SetBinLabel( 7,"kCMUS5/kMUSPB");
627       hTriggerClass->GetXaxis()->SetBinLabel( 8,"kMUSH7/kMUSHPB");
628       hTriggerClass->GetXaxis()->SetBinLabel( 9,"kMUL7/kMuonLikePB");
629       hTriggerClass->GetXaxis()->SetBinLabel(10,"kMUU7/kMuonUnlikePB");
630       hTriggerClass->GetXaxis()->SetBinLabel(11,"kEMC7/kEMC8");
631       hTriggerClass->GetXaxis()->SetBinLabel(12,"kMUS7");
632       hTriggerClass->GetXaxis()->SetBinLabel(13,"kPHI1");
633       hTriggerClass->GetXaxis()->SetBinLabel(14,"kPHI7/kPHI8/kPHOSPb");
634       hTriggerClass->GetXaxis()->SetBinLabel(15,"kEMCEJE");
635       hTriggerClass->GetXaxis()->SetBinLabel(16,"kEMCEGA");
636       hTriggerClass->GetXaxis()->SetBinLabel(17,"kCentral");
637       hTriggerClass->GetXaxis()->SetBinLabel(18,"kSemiCentral");
638       hTriggerClass->GetXaxis()->SetBinLabel(19,"kDG5");
639       hTriggerClass->GetXaxis()->SetBinLabel(20,"kZED");
640       hTriggerClass->GetXaxis()->SetBinLabel(21,"kSPI7/kSPI");
641       hTriggerClass->GetXaxis()->SetBinLabel(22,"kINT8");
642       hTriggerClass->GetXaxis()->SetBinLabel(23,"kMuonSingleLowPt8");
643       hTriggerClass->GetXaxis()->SetBinLabel(24,"kMuonSingleHighPt8");
644       hTriggerClass->GetXaxis()->SetBinLabel(25,"kMuonLikeLowPt8");
645       hTriggerClass->GetXaxis()->SetBinLabel(26,"kMuonUnlikeLowPt8");
646       hTriggerClass->GetXaxis()->SetBinLabel(27,"kMuonUnlikeLowPt0");
647       hTriggerClass->GetXaxis()->SetBinLabel(28,"kUserDefined");
648       hTriggerClass->GetXaxis()->SetBinLabel(29,"kTRD");
649       hTriggerClass->GetXaxis()->SetBinLabel(30,"kFastOnly");
650       hTriggerClass->GetXaxis()->SetBinLabel(31,"kAnyINT");
651       hTriggerClass->GetXaxis()->SetBinLabel(32,"kAny");
652       hTriggerClass->GetXaxis()->SetBinLabel(33,"V0AND");
653       hTriggerClass->GetXaxis()->SetBinLabel(34,"NOT kFastOnly");
654       hTriggerClass->GetXaxis()->SetBinLabel(35,"failed Physics Selection");
655       fHistograms->Add(hTriggerClass);
656    }
657    if(!preCut){
658       hTriggerClassSelected= new TH1F(Form("OfflineTriggerSelected %s",GetCutNumber().Data()),"OfflineTriggerSelected",34,-0.5,33.5);
659       hTriggerClassSelected->GetXaxis()->SetBinLabel( 1,"kMB");
660       hTriggerClassSelected->GetXaxis()->SetBinLabel( 2,"kINT7");
661       hTriggerClassSelected->GetXaxis()->SetBinLabel( 3,"kMUON");
662       hTriggerClassSelected->GetXaxis()->SetBinLabel( 4,"kHighMult");
663       hTriggerClassSelected->GetXaxis()->SetBinLabel( 5,"kKEMC1");
664       hTriggerClassSelected->GetXaxis()->SetBinLabel( 6,"kCINT5");
665       hTriggerClassSelected->GetXaxis()->SetBinLabel( 7,"kCMUS5/kMUSPB");
666       hTriggerClassSelected->GetXaxis()->SetBinLabel( 8,"kMUSH7/kMUSHPB");
667       hTriggerClassSelected->GetXaxis()->SetBinLabel( 9,"kMUL7/kMuonLikePB");
668       hTriggerClassSelected->GetXaxis()->SetBinLabel(10,"kMUU7/kMuonUnlikePB");
669       hTriggerClassSelected->GetXaxis()->SetBinLabel(11,"kEMC7/kEMC8");
670       hTriggerClassSelected->GetXaxis()->SetBinLabel(12,"kMUS7");
671       hTriggerClassSelected->GetXaxis()->SetBinLabel(13,"kPHI1");
672       hTriggerClassSelected->GetXaxis()->SetBinLabel(14,"kPHI7/kPHI8/kPHOSPb");
673       hTriggerClassSelected->GetXaxis()->SetBinLabel(15,"kEMCEJE");
674       hTriggerClassSelected->GetXaxis()->SetBinLabel(16,"kEMCEGA");
675       hTriggerClassSelected->GetXaxis()->SetBinLabel(17,"kCentral");
676       hTriggerClassSelected->GetXaxis()->SetBinLabel(18,"kSemiCentral");
677       hTriggerClassSelected->GetXaxis()->SetBinLabel(19,"kDG5");
678       hTriggerClassSelected->GetXaxis()->SetBinLabel(20,"kZED");
679       hTriggerClassSelected->GetXaxis()->SetBinLabel(21,"kSPI7/kSPI");
680       hTriggerClassSelected->GetXaxis()->SetBinLabel(22,"kINT8");
681       hTriggerClassSelected->GetXaxis()->SetBinLabel(23,"kMuonSingleLowPt8");
682       hTriggerClassSelected->GetXaxis()->SetBinLabel(24,"kMuonSingleHighPt8");
683       hTriggerClassSelected->GetXaxis()->SetBinLabel(25,"kMuonLikeLowPt8");
684       hTriggerClassSelected->GetXaxis()->SetBinLabel(26,"kMuonUnlikeLowPt8");
685       hTriggerClassSelected->GetXaxis()->SetBinLabel(27,"kMuonUnlikeLowPt0");
686       hTriggerClassSelected->GetXaxis()->SetBinLabel(28,"kUserDefined");
687       hTriggerClassSelected->GetXaxis()->SetBinLabel(29,"kTRD");
688       hTriggerClassSelected->GetXaxis()->SetBinLabel(30,"kFastOnly");
689       hTriggerClassSelected->GetXaxis()->SetBinLabel(31,"kAnyINT");
690       hTriggerClassSelected->GetXaxis()->SetBinLabel(32,"kAny");
691       hTriggerClassSelected->GetXaxis()->SetBinLabel(33,"V0AND");
692       hTriggerClassSelected->GetXaxis()->SetBinLabel(34,"NOT kFastOnly");
693       fHistograms->Add(hTriggerClassSelected);
694       
695       hEventPlanePhi=new TH1F(Form("EventPlaneMinusPhotonAngle %s",GetCutNumber().Data()),"EventPlaneMinusPhotonAngle",360,-TMath::Pi(),TMath::Pi());
696       fHistograms->Add(hEventPlanePhi);
697
698       
699    }
700    TH1::AddDirectory(kTRUE);
701 }
702
703 //________________________________________________________________________
704 Bool_t AliConversionCuts::InitPIDResponse(){
705    // Set Pointer to AliPIDResponse
706
707    AliAnalysisManager *man=AliAnalysisManager::GetAnalysisManager();
708    if(man) {
709       AliInputEventHandler* inputHandler = (AliInputEventHandler*) (man->GetInputEventHandler());
710       fPIDResponse = (AliPIDResponse*)inputHandler->GetPIDResponse();
711       if(fPIDResponse)return kTRUE;
712
713    }
714
715
716    return kFALSE;
717 }
718 ///________________________________________________________________________
719 Bool_t AliConversionCuts::EventIsSelected(AliVEvent *fInputEvent, AliVEvent *fMCEvent){
720    // Process Event Selection
721
722    Int_t cutindex=0;
723    if(hV0EventCuts)hV0EventCuts->Fill(cutindex);
724    cutindex++;
725
726    // Check for MC event
727    if(fMCEvent && fInputEvent->IsA()==AliESDEvent::Class()){
728       // Check if MC event is correctly loaded
729       AliMCEventHandler* mcHandler = dynamic_cast<AliMCEventHandler*> (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
730       if (!mcHandler){
731          fEventQuality = 2;
732          return kFALSE;
733       }
734       if (!mcHandler->InitOk() ){
735          fEventQuality = 2;
736          return kFALSE;
737       }
738       if (!mcHandler->TreeK() ){
739          fEventQuality = 2;
740          return kFALSE;
741       }
742       if (!mcHandler->TreeTR() ) {
743          fEventQuality = 2;
744          return kFALSE;
745       }
746    }
747
748    // Event Trigger
749 //    cout << "before event trigger" << endl;
750    if(!IsTriggerSelected(fInputEvent)){
751       if(hV0EventCuts)hV0EventCuts->Fill(cutindex);
752       fEventQuality = 3;
753       return kFALSE;
754    }
755    cutindex++;
756
757    if(fInputEvent->IsA()==AliESDEvent::Class()){
758       AliTriggerAnalysis fTriggerAnalysis;// = new AliTriggerAnalysis;
759       fHasV0AND = fTriggerAnalysis.IsOfflineTriggerFired((AliESDEvent*)fInputEvent, AliTriggerAnalysis::kV0AND);
760       if(fHasV0AND&&hTriggerClass)hTriggerClass->Fill(32);
761    }
762 //   cout << "event number " << ((AliESDEvent*)fInputEvent)->GetEventNumberInFile() << " entered"<< endl;
763
764
765    // Number of Contributors Cut
766    if(GetNumberOfContributorsVtx(fInputEvent)<=0) {
767       if(hV0EventCuts)hV0EventCuts->Fill(cutindex);
768       fEventQuality = 5;
769       return kFALSE;
770    }
771    cutindex++;
772
773    // Z Vertex Position Cut
774    if(!VertexZCut(fInputEvent)){
775       if(hV0EventCuts)hV0EventCuts->Fill(cutindex);
776       fEventQuality = 4;
777       return kFALSE;
778    }
779    cutindex++;
780
781    // Pile Up Rejection
782
783    if(fRemovePileUp){
784       if(fInputEvent->IsPileupFromSPD(3,0.8,3.,2.,5.)){
785          if(hV0EventCuts)hV0EventCuts->Fill(cutindex);
786          fEventQuality = 6;
787          return kFALSE;
788       }
789    }
790    cutindex++;
791
792    // Centrality Selection
793    if(!IsCentralitySelected(fInputEvent,fMCEvent)){
794       if(hV0EventCuts)hV0EventCuts->Fill(cutindex);
795       fEventQuality = 1;
796       return kFALSE;
797    }
798    cutindex++;
799
800    // Fill Event Histograms
801    if(hV0EventCuts)hV0EventCuts->Fill(cutindex);
802    if(hVertexZ)hVertexZ->Fill(fInputEvent->GetPrimaryVertex()->GetZ());
803    if(hCentrality)hCentrality->Fill(GetCentrality(fInputEvent));
804    if(hCentralityVsNumberOfPrimaryTracks)
805       hCentralityVsNumberOfPrimaryTracks->Fill(GetCentrality(fInputEvent),
806                                                ((AliV0ReaderV1*)AliAnalysisManager::GetAnalysisManager()
807                                                 ->GetTask("V0ReaderV1"))->GetNumberOfPrimaryTracks());
808    fEventQuality = 0;
809    return kTRUE;
810 }
811
812 ///________________________________________________________________________
813 Bool_t AliConversionCuts::PhotonIsSelectedMC(TParticle *particle,AliStack *fMCStack,Bool_t checkForConvertedGamma){
814    // MonteCarlo Photon Selection
815
816    if(!fMCStack)return kFALSE;
817
818    if (particle->GetPdgCode() == 22){
819
820
821       if( particle->Eta() > (fEtaCut + fEtaShift) || particle->Eta() < (-fEtaCut + fEtaShift) )
822          return kFALSE;
823       if(fEtaCutMin>-0.1){
824          if( particle->Eta() < (fEtaCutMin + fEtaShift) && particle->Eta() > (-fEtaCutMin + fEtaShift) )
825             return kFALSE;
826       }
827
828       if(particle->GetMother(0) >-1 && fMCStack->Particle(particle->GetMother(0))->GetPdgCode() == 22){
829          return kFALSE; // no photon as mothers!
830       }
831
832       if(particle->GetMother(0) >= fMCStack->GetNprimary()){
833          return kFALSE; // the gamma has a mother, and it is not a primary particle
834       }
835
836       if(!checkForConvertedGamma) return kTRUE; // return in case of accepted gamma
837
838       // looking for conversion gammas (electron + positron from pairbuilding (= 5) )
839       TParticle* ePos = NULL;
840       TParticle* eNeg = NULL;
841
842       if(particle->GetNDaughters() >= 2){
843          for(Int_t daughterIndex=particle->GetFirstDaughter();daughterIndex<=particle->GetLastDaughter();daughterIndex++){
844             TParticle *tmpDaughter = fMCStack->Particle(daughterIndex);
845             if(tmpDaughter->GetUniqueID() == 5){
846                if(tmpDaughter->GetPdgCode() == 11){
847                   eNeg = tmpDaughter;
848                } else if(tmpDaughter->GetPdgCode() == -11){
849                   ePos = tmpDaughter;
850                }
851             }
852          }
853       }
854
855       if(ePos == NULL || eNeg == NULL){ // means we do not have two daughters from pair production
856          return kFALSE;
857       }
858
859       if(ePos->Pt()<fSinglePtCut || eNeg->Pt()<fSinglePtCut){
860          return kFALSE; // no reconstruction below the Pt cut
861       }
862
863       if( ePos->Eta() > (fEtaCut + fEtaShift) || ePos->Eta() < (-fEtaCut + fEtaShift) ||
864           eNeg->Eta() > (fEtaCut + fEtaShift) || eNeg->Eta() < (-fEtaCut + fEtaShift) )
865          return kFALSE;
866
867       if(fEtaCutMin > -0.1){
868          if( (ePos->Eta() < (fEtaCutMin + fEtaShift) && ePos->Eta() > (-fEtaCutMin + fEtaShift)) ||
869              (eNeg->Eta() < (fEtaCutMin + fEtaShift) && eNeg->Eta() > (-fEtaCutMin + fEtaShift)) )
870             return kFALSE;
871       }
872
873       if(ePos->R()>fMaxR){
874          return kFALSE; // cuts on distance from collision point
875       }
876
877       if(abs(ePos->Vz()) > fMaxZ){
878          return kFALSE;  // outside material
879       }
880       if(abs(eNeg->Vz()) > fMaxZ){
881          return kFALSE;  // outside material
882       }
883
884       if( ePos->R() <= ((abs(ePos->Vz()) * fLineCutZRSlope) - fLineCutZValue)){
885          return kFALSE;  // line cut to exclude regions where we do not reconstruct
886       } else if ( fEtaCutMin != -0.1 &&   ePos->R() >= ((abs(ePos->Vz()) * fLineCutZRSlopeMin) - fLineCutZValueMin)){
887          return kFALSE;
888       }
889
890       if( eNeg->R() <= ((abs(eNeg->Vz()) * fLineCutZRSlope) - fLineCutZValue)){
891          return kFALSE; // line cut to exclude regions where we do not reconstruct
892       } else if ( fEtaCutMin != -0.1 &&   eNeg->R() >= ((abs(eNeg->Vz()) * fLineCutZRSlopeMin) - fLineCutZValueMin)){
893          return kFALSE;
894       }
895
896       return kTRUE;
897       //if(AcceptanceCut(particle,ePos,eNeg))return kTRUE;
898    }
899    return kFALSE;
900 }
901 ///________________________________________________________________________
902 Bool_t AliConversionCuts::PhotonIsSelectedAODMC(AliAODMCParticle *particle,TClonesArray *aodmcArray,Bool_t checkForConvertedGamma){
903    // MonteCarlo Photon Selection
904
905    if(!aodmcArray)return kFALSE;
906
907    if (particle->GetPdgCode() == 22){
908       if( particle->Eta() > (fEtaCut + fEtaShift) || particle->Eta() < (-fEtaCut + fEtaShift) )
909          return kFALSE;
910       if(fEtaCutMin>-0.1){
911          if( particle->Eta() < (fEtaCutMin + fEtaShift) && particle->Eta() > (-fEtaCutMin + fEtaShift) )
912             return kFALSE;
913       }
914
915       if(particle->GetMother() > -1){
916          if((static_cast<AliAODMCParticle*>(aodmcArray->At(particle->GetMother())))->GetPdgCode() == 22){
917             return kFALSE; // no photon as mothers!
918          }
919          if(!(static_cast<AliAODMCParticle*>(aodmcArray->At(particle->GetMother()))->IsPrimary())){
920             return kFALSE; // the gamma has a mother, and it is not a primary particle
921          }
922       }
923
924       if(!checkForConvertedGamma) return kTRUE; // return in case of accepted gamma
925
926       // looking for conversion gammas (electron + positron from pairbuilding (= 5) )
927       AliAODMCParticle* ePos = NULL;
928       AliAODMCParticle* eNeg = NULL;
929
930       if(particle->GetNDaughters() >= 2){
931          for(Int_t daughterIndex=particle->GetDaughter(0);daughterIndex<=particle->GetDaughter(1);daughterIndex++){
932             AliAODMCParticle *tmpDaughter = static_cast<AliAODMCParticle*>(aodmcArray->At(daughterIndex));
933             if(!tmpDaughter) continue;
934             if(((tmpDaughter->GetMCProcessCode())) == 5){    // STILL A BUG IN ALIROOT >>8 HAS TPO BE REMOVED AFTER FIX
935                if(tmpDaughter->GetPdgCode() == 11){
936                   eNeg = tmpDaughter;
937                } else if(tmpDaughter->GetPdgCode() == -11){
938                   ePos = tmpDaughter;
939                }
940             }
941          }
942       }
943
944       if(ePos == NULL || eNeg == NULL){ // means we do not have two daughters from pair production
945          return kFALSE;
946       }
947
948       if(ePos->Pt()<fSinglePtCut || eNeg->Pt()<fSinglePtCut){
949          return kFALSE; // no reconstruction below the Pt cut
950       }
951
952       if( ePos->Eta() > (fEtaCut + fEtaShift) || ePos->Eta() < (-fEtaCut + fEtaShift) ||
953           eNeg->Eta() > (fEtaCut + fEtaShift) || eNeg->Eta() < (-fEtaCut + fEtaShift) )
954          return kFALSE;
955
956       if(fEtaCutMin > -0.1){
957          if( (ePos->Eta() < (fEtaCutMin + fEtaShift) && ePos->Eta() > (-fEtaCutMin + fEtaShift)) ||
958              (eNeg->Eta() < (fEtaCutMin + fEtaShift) && eNeg->Eta() > (-fEtaCutMin + fEtaShift)) )
959             return kFALSE;
960       }
961
962       Double_t rPos = sqrt( (ePos->Xv()*ePos->Xv()) + (ePos->Yv()*ePos->Yv()) );
963       Double_t rNeg = sqrt( (eNeg->Xv()*eNeg->Xv()) + (eNeg->Yv()*eNeg->Yv()) );
964
965       if(rPos>fMaxR){
966          return kFALSE; // cuts on distance from collision point
967       }
968       if(abs(ePos->Zv()) > fMaxZ){
969          return kFALSE;  // outside material
970       }
971       if(abs(eNeg->Zv()) > fMaxZ){
972          return kFALSE;  // outside material
973       }
974
975       if( rPos <= ((abs(ePos->Zv()) * fLineCutZRSlope) - fLineCutZValue)){
976          return kFALSE;  // line cut to exclude regions where we do not reconstruct
977       } else if ( fEtaCutMin != -0.1 &&   rPos >= ((abs(ePos->Zv()) * fLineCutZRSlopeMin) - fLineCutZValueMin)){
978          return kFALSE;
979       }
980
981       if( rNeg <= ((abs(eNeg->Zv()) * fLineCutZRSlope) - fLineCutZValue)){
982          return kFALSE; // line cut to exclude regions where we do not reconstruct
983       } else if ( fEtaCutMin != -0.1 &&   rNeg >= ((abs(eNeg->Zv()) * fLineCutZRSlopeMin) - fLineCutZValueMin)){
984          return kFALSE;
985       }
986
987       return kTRUE;
988       //if(AcceptanceCut(particle,ePos,eNeg))return kTRUE;
989    }
990    return kFALSE;
991 }
992 ///________________________________________________________________________
993 Bool_t AliConversionCuts::PhotonCuts(AliConversionPhotonBase *photon,AliVEvent *event)
994 {   // Specific Photon Cuts
995
996    Int_t cutIndex = 0;
997    if(hPhotonCuts)hPhotonCuts->Fill(cutIndex);
998    cutIndex++;
999
1000    // Fill Histos before Cuts
1001    if(hInvMassbefore)hInvMassbefore->Fill(photon->GetMass());
1002
1003    if(hArmenterosbefore)hArmenterosbefore->Fill(photon->GetArmenterosAlpha(),photon->GetArmenterosQt());
1004
1005    // Gamma selection based on QT from Armenteros
1006    if(fDoQtGammaSelection == kTRUE){
1007       if(!ArmenterosQtCut(photon)){
1008          if(hPhotonCuts)hPhotonCuts->Fill(cutIndex); //1
1009          return kFALSE;
1010       }
1011    }
1012    cutIndex++; //2
1013
1014    // Chi Cut
1015    if(photon->GetChi2perNDF() > fChi2CutConversion || photon->GetChi2perNDF() <=0){
1016       {
1017          if(hPhotonCuts)hPhotonCuts->Fill(cutIndex); //2
1018          return kFALSE;
1019       }
1020    }
1021    cutIndex++;//3
1022
1023    // Reconstruction Acceptance Cuts
1024    if(!AcceptanceCuts(photon)){
1025       if(hPhotonCuts)hPhotonCuts->Fill(cutIndex); //3
1026       return kFALSE;
1027    }
1028
1029    cutIndex++; //4
1030    // Asymmetry Cut
1031    if(fDoPhotonAsymmetryCut == kTRUE){
1032       if(!AsymmetryCut(photon,event)){
1033          if(hPhotonCuts)hPhotonCuts->Fill(cutIndex); //4
1034          return kFALSE;
1035       }
1036    }
1037
1038    //Check the pid probability
1039    cutIndex++; //5
1040    if(!PIDProbabilityCut(photon, event)) {
1041       if(hPhotonCuts)hPhotonCuts->Fill(cutIndex); //5
1042       return kFALSE;
1043    }
1044
1045    cutIndex++; //6
1046    if(!CorrectedTPCClusterCut(photon, event)) {
1047       if(hPhotonCuts)hPhotonCuts->Fill(cutIndex); //6
1048       return kFALSE;
1049    }
1050
1051    Double_t magField = event->GetMagneticField();
1052    if( magField  < 0.0 ){
1053       magField =  1.0;
1054    } else {
1055       magField =  -1.0;
1056    }
1057    
1058    AliVTrack * electronCandidate = GetTrack(event,photon->GetTrackLabelNegative() );
1059    AliVTrack * positronCandidate = GetTrack(event,photon->GetTrackLabelPositive() );
1060    Double_t deltaPhi = magField * TVector2::Phi_mpi_pi( electronCandidate->Phi()-positronCandidate->Phi());
1061
1062    cutIndex++; //7
1063    if(!PsiPairCut(photon)) {
1064       if(hPhotonCuts)hPhotonCuts->Fill(cutIndex); //7
1065       return kFALSE;
1066    }
1067
1068    cutIndex++; //8
1069    if(!CosinePAngleCut(photon, event)) {
1070       if(hPhotonCuts)hPhotonCuts->Fill(cutIndex); //8
1071       return kFALSE;
1072    }
1073
1074    AliAODConversionPhoton* photonAOD = dynamic_cast<AliAODConversionPhoton*>(photon);
1075    if (photonAOD){
1076       photonAOD->CalculateDistanceOfClossetApproachToPrimVtx(event->GetPrimaryVertex());
1077
1078       cutIndex++; //9
1079       if(photonAOD->GetDCArToPrimVtx() > fDCARPrimVtxCut) { //DCA R cut of photon to primary vertex
1080          if(hPhotonCuts)hPhotonCuts->Fill(cutIndex); //9
1081          return kFALSE;
1082       }
1083
1084       cutIndex++; //10
1085       if(abs(photonAOD->GetDCAzToPrimVtx()) > fDCAZPrimVtxCut) { //DCA Z cut of photon to primary vertex
1086          if(hPhotonCuts)hPhotonCuts->Fill(cutIndex); //10
1087          return kFALSE;
1088       }
1089    } else {
1090       cutIndex++; //9
1091       cutIndex++; //10
1092    }
1093    cutIndex++; //11
1094    if(hPhotonCuts)hPhotonCuts->Fill(cutIndex); //11
1095
1096    // Histos after Cuts
1097    if(hInvMassafter)hInvMassafter->Fill(photon->GetMass());
1098    if(hArmenterosafter)hArmenterosafter->Fill(photon->GetArmenterosAlpha(),photon->GetArmenterosQt());
1099    if(hPsiPairDeltaPhiafter)hPsiPairDeltaPhiafter->Fill(deltaPhi,photon->GetPsiPair());
1100    return kTRUE;
1101
1102 }
1103
1104 ///________________________________________________________________________
1105 Bool_t AliConversionCuts::CorrectedTPCClusterCut(AliConversionPhotonBase *photon, AliVEvent * event)
1106 {   //Cut on corrected TPC Cluster Info
1107
1108    AliVTrack * negTrack = GetTrack(event, photon->GetTrackLabelNegative());
1109    AliVTrack * posTrack = GetTrack(event, photon->GetTrackLabelPositive());
1110
1111    if(!negTrack||!posTrack)return kFALSE;
1112
1113    Double_t negclsToF=0;
1114
1115    if (!fUseCorrectedTPCClsInfo ){
1116       if(negTrack->GetTPCNclsF()!=0){
1117          negclsToF = (Double_t)negTrack->GetNcls(1)/(Double_t)negTrack->GetTPCNclsF();}// Ncluster/Nfindablecluster
1118    }
1119    else {
1120       negclsToF = negTrack->GetTPCClusterInfo(2,0,GetFirstTPCRow(photon->GetConversionRadius()));
1121    }
1122
1123    Double_t posclsToF = 0.;
1124    if (!fUseCorrectedTPCClsInfo ){
1125       if(posTrack->GetTPCNclsF()!=0){
1126          posclsToF = (Double_t)posTrack->GetNcls(1)/(Double_t)posTrack->GetTPCNclsF();
1127       }
1128    }else{
1129       posclsToF = posTrack->GetTPCClusterInfo(2,0,GetFirstTPCRow(photon->GetConversionRadius()));
1130    }
1131
1132    if( negclsToF < fMinClsTPCToF || posclsToF < fMinClsTPCToF ){
1133       return kFALSE;
1134    }
1135
1136    return kTRUE;
1137 }
1138
1139 ///________________________________________________________________________
1140 Bool_t AliConversionCuts::PhotonIsSelected(AliConversionPhotonBase *photon, AliVEvent * event)
1141 {
1142    //Selection of Reconstructed Photons
1143
1144    FillPhotonCutIndex(kPhotonIn);
1145
1146    if(event->IsA()==AliESDEvent::Class()) {
1147      if(!SelectV0Finder( ( ((AliESDEvent*)event)->GetV0(photon->GetV0Index()))->GetOnFlyStatus() ) ){
1148          FillPhotonCutIndex(kOnFly);
1149          return kFALSE;
1150       }
1151    }
1152    // else if(event->IsA()==AliAODEvent::Class()) {
1153    //    if(!SelectV0Finder( ( ((AliAODEvent*)event)->GetV0(photon->GetV0Index())) ) ){
1154    //       FillPhotonCutIndex(kOnFly);
1155    //       return kFALSE;
1156    //    }
1157    // }
1158
1159    // Get Tracks
1160    AliVTrack * negTrack = GetTrack(event, photon->GetTrackLabelNegative());
1161    AliVTrack * posTrack = GetTrack(event, photon->GetTrackLabelPositive());
1162
1163    if(!negTrack || !posTrack) {
1164       FillPhotonCutIndex(kNoTracks);
1165       return kFALSE;
1166    }
1167    photon->DeterminePhotonQuality(negTrack,posTrack);
1168    // Track Cuts
1169    if(!TracksAreSelected(negTrack, posTrack)){
1170       FillPhotonCutIndex(kTrackCuts);
1171       return kFALSE;
1172    }
1173
1174    // dEdx Cuts
1175    if(!dEdxCuts(negTrack) || !dEdxCuts(posTrack)) {
1176       FillPhotonCutIndex(kdEdxCuts);
1177       return kFALSE;
1178    }
1179
1180    // Photon Cuts
1181    if(!PhotonCuts(photon,event)){
1182       FillPhotonCutIndex(kPhotonCuts);
1183       return kFALSE;
1184    }
1185
1186    // Photon passed cuts
1187    FillPhotonCutIndex(kPhotonOut);
1188    return kTRUE;
1189 }
1190
1191 ///________________________________________________________________________
1192 Bool_t AliConversionCuts::ArmenterosQtCut(AliConversionPhotonBase *photon)
1193 {   // Armenteros Qt Cut
1194
1195    if(fDo2DQt){
1196       if ( !(TMath::Power(photon->GetArmenterosAlpha()/0.95,2)+TMath::Power(photon->GetArmenterosQt()/fQtMax,2) < 1) ){
1197          return kFALSE;
1198       }
1199    } else {
1200       if(photon->GetArmenterosQt()>fQtMax){
1201          return kFALSE;
1202       }
1203    }
1204    return kTRUE;
1205 }
1206
1207
1208 ///________________________________________________________________________
1209 Bool_t AliConversionCuts::AcceptanceCuts(AliConversionPhotonBase *photon) {
1210    // Exclude certain areas for photon reconstruction
1211
1212    Int_t cutIndex=0;
1213    if(hAcceptanceCuts)hAcceptanceCuts->Fill(cutIndex);
1214    cutIndex++;
1215
1216    if(photon->GetConversionRadius()>fMaxR){ // cuts on distance from collision point
1217       if(hAcceptanceCuts)hAcceptanceCuts->Fill(cutIndex);
1218       return kFALSE;
1219    }
1220    cutIndex++;
1221
1222    if(photon->GetConversionRadius()<fMinR){ // cuts on distance from collision point
1223       if(hAcceptanceCuts)hAcceptanceCuts->Fill(cutIndex);
1224       return kFALSE;
1225    }
1226    cutIndex++;
1227
1228    if(photon->GetConversionRadius() <= ((abs(photon->GetConversionZ())*fLineCutZRSlope)-fLineCutZValue)){
1229       if(hAcceptanceCuts)hAcceptanceCuts->Fill(cutIndex);
1230       return kFALSE;
1231    }
1232    else if (fUseEtaMinCut &&  photon->GetConversionRadius() >= ((abs(photon->GetConversionZ())*fLineCutZRSlopeMin)-fLineCutZValueMin )){
1233       if(hAcceptanceCuts)hAcceptanceCuts->Fill(cutIndex);
1234       return kFALSE;
1235    }
1236    cutIndex++;
1237
1238    if(abs(photon->GetConversionZ()) > fMaxZ ){ // cuts out regions where we do not reconstruct
1239       if(hAcceptanceCuts)hAcceptanceCuts->Fill(cutIndex);
1240       return kFALSE;
1241    }
1242    cutIndex++;
1243
1244
1245    if( photon->GetPhotonEta() > (fEtaCut + fEtaShift)    || photon->GetPhotonEta() < (-fEtaCut + fEtaShift) ){
1246       if(hAcceptanceCuts)hAcceptanceCuts->Fill(cutIndex);
1247       return kFALSE;
1248    }
1249    if(fEtaCutMin>-0.1){
1250       if( photon->GetPhotonEta() < (fEtaCutMin + fEtaShift) && photon->GetPhotonEta() > (-fEtaCutMin + fEtaShift) ){
1251          if(hAcceptanceCuts)hAcceptanceCuts->Fill(cutIndex);
1252          return kFALSE;
1253       }
1254    }
1255    cutIndex++;
1256
1257    if(photon->GetPhotonPt()<fPtCut){
1258       if(hAcceptanceCuts)hAcceptanceCuts->Fill(cutIndex);
1259       return kFALSE;
1260    }
1261    cutIndex++;
1262
1263    if(hAcceptanceCuts)hAcceptanceCuts->Fill(cutIndex);
1264
1265    return kTRUE;
1266 }
1267
1268
1269 ///________________________________________________________________________
1270 Bool_t AliConversionCuts::SpecificTrackCuts(AliAODTrack * negTrack, AliAODTrack * posTrack,Int_t &cutIndex) {
1271    // Track Cuts which require AOD/ESD specific implementation
1272
1273    if( !negTrack->IsOn(AliESDtrack::kTPCrefit)  || !posTrack->IsOn(AliESDtrack::kTPCrefit)   )  {
1274       if(hTrackCuts)hTrackCuts->Fill(cutIndex);
1275       return kFALSE;
1276    }
1277    cutIndex++;
1278
1279    AliAODVertex * NegVtxType=negTrack->GetProdVertex();
1280    AliAODVertex * PosVtxType=posTrack->GetProdVertex();
1281    if( (NegVtxType->GetType())==AliAODVertex::kKink || (PosVtxType->GetType())==AliAODVertex::kKink) {
1282       if(hTrackCuts)hTrackCuts->Fill(cutIndex);
1283       return kFALSE;
1284    }
1285    return kTRUE;
1286
1287 }
1288
1289
1290 ///________________________________________________________________________
1291 Bool_t AliConversionCuts::SpecificTrackCuts(AliESDtrack * negTrack, AliESDtrack * posTrack,Int_t &cutIndex) {
1292    // Track Cuts which require AOD/ESD specific implementation
1293
1294    if( !negTrack->IsOn(AliESDtrack::kTPCrefit)  || !posTrack->IsOn(AliESDtrack::kTPCrefit)   )  {
1295       if(hTrackCuts)hTrackCuts->Fill(cutIndex);
1296       return kFALSE;
1297    }
1298    cutIndex++;
1299
1300    if(negTrack->GetKinkIndex(0) > 0  || posTrack->GetKinkIndex(0) > 0 ) {
1301       if(hTrackCuts)hTrackCuts->Fill(cutIndex);
1302       return kFALSE;
1303    }
1304    return kTRUE;
1305 }
1306
1307
1308
1309 ///________________________________________________________________________
1310 Bool_t AliConversionCuts::TracksAreSelected(AliVTrack * negTrack, AliVTrack * posTrack) {
1311    // Track Selection for Photon Reconstruction
1312
1313    Int_t cutIndex=0;
1314    if(hTrackCuts)hTrackCuts->Fill(cutIndex);
1315    cutIndex++;
1316
1317    // avoid like sign
1318    if(negTrack->Charge() == posTrack->Charge()) {
1319       if(hTrackCuts)hTrackCuts->Fill(cutIndex);
1320       return kFALSE;
1321    }
1322    cutIndex++;
1323
1324    // Number of TPC Clusters
1325
1326
1327    if( negTrack->GetNcls(1) < fMinClsTPC || posTrack->GetNcls(1) < fMinClsTPC ) {
1328       if(hTrackCuts)hTrackCuts->Fill(cutIndex);
1329       return kFALSE;
1330    }
1331    cutIndex++;
1332
1333    // Acceptance
1334    if( posTrack->Eta() > (fEtaCut + fEtaShift) || posTrack->Eta() < (-fEtaCut + fEtaShift) ||
1335        negTrack->Eta() > (fEtaCut + fEtaShift) || negTrack->Eta() < (-fEtaCut + fEtaShift) ){
1336       if(hTrackCuts)hTrackCuts->Fill(cutIndex);
1337       return kFALSE;
1338    }
1339    if(fEtaCutMin>-0.1){
1340       if( (posTrack->Eta() < (fEtaCutMin + fEtaShift) && posTrack->Eta() > (-fEtaCutMin + fEtaShift)) ||
1341           (negTrack->Eta() < (fEtaCutMin + fEtaShift) && negTrack->Eta() > (-fEtaCutMin + fEtaShift)) ){
1342          if(hTrackCuts)hTrackCuts->Fill(cutIndex);
1343          return kFALSE;
1344       }
1345    }
1346    cutIndex++;
1347
1348    // Single Pt Cut
1349    if( negTrack->Pt()< fSinglePtCut || posTrack->Pt()< fSinglePtCut){
1350       if(hTrackCuts)hTrackCuts->Fill(cutIndex);
1351       return kFALSE;
1352    }
1353    cutIndex++;
1354
1355    // AOD ESD specific cuts
1356    Bool_t passCuts = kTRUE;
1357
1358    if(negTrack->IsA()==AliAODTrack::Class()) {
1359       passCuts = passCuts * SpecificTrackCuts(static_cast<AliAODTrack*>(negTrack), static_cast<AliAODTrack*>(posTrack),cutIndex);
1360    } else {
1361       passCuts = passCuts * SpecificTrackCuts(static_cast<AliESDtrack*>(negTrack), static_cast<AliESDtrack*>(posTrack),cutIndex);
1362    }
1363
1364    if(!passCuts){
1365       if(hTrackCuts)hTrackCuts->Fill(cutIndex);
1366       return kFALSE;
1367    }
1368    cutIndex++;
1369
1370    if(hTrackCuts)hTrackCuts->Fill(cutIndex);
1371
1372    return kTRUE;
1373
1374 }
1375
1376 ///________________________________________________________________________
1377 Bool_t AliConversionCuts::dEdxCuts(AliVTrack *fCurrentTrack){
1378    // Electron Identification Cuts for Photon reconstruction
1379
1380    if(!fPIDResponse){InitPIDResponse();}// Try to reinitialize PID Response
1381    if(!fPIDResponse){AliError("No PID Response"); return kTRUE;}// if still missing fatal error
1382
1383    Int_t cutIndex=0;
1384    if(hdEdxCuts)hdEdxCuts->Fill(cutIndex);
1385    if(hTPCdEdxSigbefore)hTPCdEdxSigbefore->Fill(fCurrentTrack->P(),fPIDResponse->NumberOfSigmasTPC(fCurrentTrack, AliPID::kElectron));
1386    if(hTPCdEdxbefore)hTPCdEdxbefore->Fill(fCurrentTrack->P(),fCurrentTrack->GetTPCsignal());
1387    cutIndex++;
1388
1389
1390    if(fDodEdxSigmaCut == kTRUE){
1391       // TPC Electron Line
1392       if( fPIDResponse->NumberOfSigmasTPC(fCurrentTrack,AliPID::kElectron)<fPIDnSigmaBelowElectronLine ||
1393           fPIDResponse->NumberOfSigmasTPC(fCurrentTrack,AliPID::kElectron)>fPIDnSigmaAboveElectronLine){
1394
1395          if(hdEdxCuts)hdEdxCuts->Fill(cutIndex);
1396          return kFALSE;
1397       }
1398       cutIndex++;
1399
1400       // TPC Pion Line
1401       if( fCurrentTrack->P()>fPIDMinPnSigmaAbovePionLine && fCurrentTrack->P()<fPIDMaxPnSigmaAbovePionLine ){
1402          if(fPIDResponse->NumberOfSigmasTPC(fCurrentTrack,AliPID::kElectron)>fPIDnSigmaBelowElectronLine &&
1403             fPIDResponse->NumberOfSigmasTPC(fCurrentTrack,AliPID::kElectron)<fPIDnSigmaAboveElectronLine&&
1404             fPIDResponse->NumberOfSigmasTPC(fCurrentTrack,AliPID::kPion)<fPIDnSigmaAbovePionLine){
1405
1406             if(hdEdxCuts)hdEdxCuts->Fill(cutIndex);
1407             return kFALSE;
1408          }
1409       }
1410       cutIndex++;
1411
1412       // High Pt Pion rej
1413       if( fCurrentTrack->P()>fPIDMaxPnSigmaAbovePionLine ){
1414          if(fPIDResponse->NumberOfSigmasTPC(fCurrentTrack,AliPID::kElectron)>fPIDnSigmaBelowElectronLine &&
1415             fPIDResponse->NumberOfSigmasTPC(fCurrentTrack,AliPID::kElectron)<fPIDnSigmaAboveElectronLine &&
1416             fPIDResponse->NumberOfSigmasTPC(fCurrentTrack,AliPID::kPion)<fPIDnSigmaAbovePionLineHighPt){
1417
1418             if(hdEdxCuts)hdEdxCuts->Fill(cutIndex);
1419             return kFALSE;
1420          }
1421       }
1422       cutIndex++;
1423    }
1424    else{cutIndex+=3;}
1425
1426    if(fDoKaonRejectionLowP == kTRUE){
1427       if(fCurrentTrack->P()<fPIDMinPKaonRejectionLowP ){
1428          if( abs(fPIDResponse->NumberOfSigmasTPC(fCurrentTrack,AliPID::kKaon))<fPIDnSigmaAtLowPAroundKaonLine){
1429
1430             if(hdEdxCuts)hdEdxCuts->Fill(cutIndex);
1431             return kFALSE;
1432          }
1433       }
1434    }
1435    cutIndex++;
1436
1437    if(fDoProtonRejectionLowP == kTRUE){
1438       if( fCurrentTrack->P()<fPIDMinPProtonRejectionLowP ){
1439          if( abs(fPIDResponse->NumberOfSigmasTPC(fCurrentTrack,AliPID::kProton))<fPIDnSigmaAtLowPAroundProtonLine){
1440
1441             if(hdEdxCuts)hdEdxCuts->Fill(cutIndex);
1442             return kFALSE;
1443          }
1444       }
1445    }
1446    cutIndex++;
1447
1448    if(fDoPionRejectionLowP == kTRUE){
1449       if( fCurrentTrack->P()<fPIDMinPPionRejectionLowP ){
1450          if( abs(fPIDResponse->NumberOfSigmasTPC(fCurrentTrack,AliPID::kPion))<fPIDnSigmaAtLowPAroundPionLine){
1451
1452             if(hdEdxCuts)hdEdxCuts->Fill(cutIndex);
1453             return kFALSE;
1454          }
1455       }
1456    }
1457    cutIndex++;
1458
1459
1460    // cout<<"Start"<<endl;
1461    // AliPIDResponse::EDetPidStatus status=fPIDResponse->CheckPIDStatus(AliPIDResponse::kTOF,fCurrentTrack);
1462
1463    // if( ( (status & AliVTrack::kTOFout) == AliVTrack::kTOFout ) && ( (status & AliVTrack::kTIME) == AliVTrack::kTIME ))
1464    //    {cout<<"TOF DA"<<endl;}
1465    // if(status == AliPIDResponse::kDetPidOk){
1466    //    Float_t probMis = fPIDResponse->GetTOFMismatchProbability(fCurrentTrack);
1467    //    cout<<"--> "<<probMis<<endl;
1468    //    if(probMis > 0.01){
1469
1470    //    }
1471    // }
1472
1473    if((fCurrentTrack->GetStatus() & AliESDtrack::kTOFpid) && !(fCurrentTrack->GetStatus() & AliESDtrack::kTOFmismatch)){
1474       if(hTOFbefore){
1475          Double_t t0 = fPIDResponse->GetTOFResponse().GetStartTime(fCurrentTrack->P());
1476          Double_t times[5];
1477          fCurrentTrack->GetIntegratedTimes(times);
1478          Double_t TOFsignal = fCurrentTrack->GetTOFsignal();
1479          Double_t dT = TOFsignal - t0 - times[0];
1480          hTOFbefore->Fill(fCurrentTrack->P(),dT);
1481       }
1482       if(hTOFSigbefore) hTOFSigbefore->Fill(fCurrentTrack->P(),fPIDResponse->NumberOfSigmasTOF(fCurrentTrack, AliPID::kElectron));
1483       if(fUseTOFpid){
1484          if(fPIDResponse->NumberOfSigmasTOF(fCurrentTrack, AliPID::kElectron)>fTofPIDnSigmaAboveElectronLine ||
1485             fPIDResponse->NumberOfSigmasTOF(fCurrentTrack, AliPID::kElectron)<fTofPIDnSigmaBelowElectronLine ){
1486             if(hdEdxCuts)hdEdxCuts->Fill(cutIndex);
1487             return kFALSE;
1488          }
1489       }
1490       if(hTOFSigafter)hTOFSigafter->Fill(fCurrentTrack->P(),fPIDResponse->NumberOfSigmasTOF(fCurrentTrack, AliPID::kElectron));
1491    }
1492    cutIndex++;
1493
1494    // Apply TRD PID
1495    if(fDoTRDPID){
1496       if(!fPIDResponse->IdentifiedAsElectronTRD(fCurrentTrack,fPIDTRDEfficiency)){
1497          if(hdEdxCuts)hdEdxCuts->Fill(cutIndex);
1498          return kFALSE;
1499       }
1500    }
1501    cutIndex++;
1502
1503    if(hdEdxCuts)hdEdxCuts->Fill(cutIndex);
1504    if(hTPCdEdxSigafter)hTPCdEdxSigafter->Fill(fCurrentTrack->P(),fPIDResponse->NumberOfSigmasTPC(fCurrentTrack, AliPID::kElectron));
1505    if(hTPCdEdxafter)hTPCdEdxafter->Fill(fCurrentTrack->P(),fCurrentTrack->GetTPCsignal());
1506   
1507    return kTRUE;
1508 }
1509
1510 ///________________________________________________________________________
1511 Bool_t AliConversionCuts::AsymmetryCut(AliConversionPhotonBase * photon,AliVEvent *event) {
1512    // Cut on Energy Assymetry
1513
1514    for(Int_t ii=0;ii<2;ii++){
1515
1516       AliVTrack *track=GetTrack(event,photon->GetTrackLabel(ii));
1517
1518       if( track->P() > fMinPPhotonAsymmetryCut ){
1519          Double_t trackNegAsy=0;
1520          if (photon->GetPhotonP()!=0.){
1521             trackNegAsy= track->P()/photon->GetPhotonP();
1522          }
1523
1524          if( trackNegAsy<fMinPhotonAsymmetry ||trackNegAsy>(1.- fMinPhotonAsymmetry)){
1525             return kFALSE;
1526          }
1527       }
1528
1529    }
1530    return kTRUE;
1531 }
1532
1533 ///________________________________________________________________________
1534 AliVTrack *AliConversionCuts::GetTrack(AliVEvent * event, Int_t label){
1535    //Returns pointer to the track with given ESD label
1536    //(Important for AOD implementation, since Track array in AOD data is different
1537    //from ESD array, but ESD tracklabels are stored in AOD Tracks)
1538
1539    AliESDEvent * esdEvent = dynamic_cast<AliESDEvent*>(event);
1540    if(esdEvent) {
1541       if(label > event->GetNumberOfTracks() ) return NULL;
1542       AliESDtrack * track = esdEvent->GetTrack(label);
1543       return track;
1544
1545    } else {
1546       AliVTrack * track = 0x0;
1547       if(((AliV0ReaderV1*)AliAnalysisManager::GetAnalysisManager()->GetTask("V0ReaderV1"))->AreAODsRelabeled()){
1548          track = dynamic_cast<AliVTrack*>(event->GetTrack(label));
1549          return track;
1550       }
1551       else{
1552          for(Int_t ii=0; ii<event->GetNumberOfTracks(); ii++) {
1553             track = dynamic_cast<AliVTrack*>(event->GetTrack(ii));
1554             if(track){
1555                if(track->GetID() == label) {
1556                   return track;
1557                }
1558             }
1559          }
1560       }
1561    }
1562    //AliDebug(5,(Form("track not found %d %d",label,event->GetNumberOfTracks()));
1563    return NULL;
1564 }
1565
1566 ///________________________________________________________________________
1567 AliESDtrack *AliConversionCuts::GetESDTrack(AliESDEvent * event, Int_t label){
1568    //Returns pointer to the track with given ESD label
1569    //(Important for AOD implementation, since Track array in AOD data is different
1570    //from ESD array, but ESD tracklabels are stored in AOD Tracks)
1571
1572    if(event) {
1573       if(label > event->GetNumberOfTracks() ) return NULL;
1574       AliESDtrack * track = event->GetTrack(label);
1575       return track;
1576    }
1577    //AliDebug(5,(Form("track not found %d %d",label,event->GetNumberOfTracks()));
1578    return NULL;
1579 }
1580
1581
1582
1583 ///________________________________________________________________________
1584 Bool_t AliConversionCuts::PIDProbabilityCut(AliConversionPhotonBase *photon, AliVEvent * event){
1585    // Cut on Electron Probability for Photon Reconstruction
1586
1587    AliESDEvent * esdEvent = dynamic_cast<AliESDEvent*>(event);
1588
1589    if(esdEvent){
1590
1591       Bool_t iResult=kFALSE;
1592
1593       Double_t *posProbArray = new Double_t[AliPID::kSPECIES];
1594       Double_t *negProbArray = new Double_t[AliPID::kSPECIES];
1595
1596       AliESDtrack* negTrack   = esdEvent->GetTrack(photon->GetTrackLabelNegative());
1597       AliESDtrack* posTrack   = esdEvent->GetTrack(photon->GetTrackLabelPositive());
1598
1599       if(negProbArray && posProbArray){
1600
1601          negTrack->GetTPCpid(negProbArray);
1602          posTrack->GetTPCpid(posProbArray);
1603
1604          if(negProbArray[AliPID::kElectron]>=fPIDProbabilityCutNegativeParticle && posProbArray[AliPID::kElectron]>=fPIDProbabilityCutPositiveParticle){
1605             iResult=kTRUE;
1606          }
1607       }
1608
1609       delete [] posProbArray;
1610       delete [] negProbArray;
1611       return iResult;
1612
1613    } else {
1614       ///Not possible for AODs
1615       return kTRUE;
1616    }
1617
1618
1619
1620
1621 }
1622
1623
1624 ///________________________________________________________________________
1625 Bool_t AliConversionCuts::AcceptanceCut(TParticle *particle, TParticle * ePos,TParticle* eNeg){
1626    // MC Acceptance Cuts
1627    //(Certain areas were excluded for photon reconstruction)
1628
1629    if(particle->R()>fMaxR){
1630       return kFALSE;}
1631
1632    if(ePos->R()>fMaxR){
1633       return kFALSE;
1634    }
1635
1636    if(ePos->R()<fMinR){
1637       return kFALSE;
1638    }
1639
1640    if( ePos->R() <= ((abs(ePos->Vz())*fLineCutZRSlope)-fLineCutZValue)){
1641       return kFALSE;
1642    }
1643    else if (fUseEtaMinCut &&  ePos->R() >= ((abs(ePos->Vz())*fLineCutZRSlopeMin)-fLineCutZValueMin )){
1644       return kFALSE;
1645    }
1646
1647    if(abs(eNeg->Vz()) > fMaxZ ){ // cuts out regions where we do not reconstruct
1648       return kFALSE;
1649    }
1650
1651    if(eNeg->Vz()!=ePos->Vz()||eNeg->R()!=ePos->R()){
1652       return kFALSE;
1653    }
1654
1655    if(abs(ePos->Vz()) > fMaxZ ){ // cuts out regions where we do not reconstruct
1656       return kFALSE;
1657    }
1658
1659
1660    if( particle->Eta() > (fEtaCut + fEtaShift) || particle->Eta() < (-fEtaCut + fEtaShift) ){
1661       return kFALSE;
1662    }
1663    if( ePos->Eta() > (fEtaCut + fEtaShift) || ePos->Eta() < (-fEtaCut + fEtaShift) ){
1664       return kFALSE;
1665    }
1666    if( eNeg->Eta() > (fEtaCut + fEtaShift) || eNeg->Eta() < (-fEtaCut + fEtaShift) ){
1667       return kFALSE;
1668    }
1669    if(fEtaCutMin>-0.1){
1670       if( particle->Eta() < (fEtaCutMin + fEtaShift) && particle->Eta() > (-fEtaCutMin + fEtaShift) ){
1671          return kFALSE;
1672       }
1673       if( ePos->Eta() < (fEtaCutMin + fEtaShift) && ePos->Eta() > (-fEtaCutMin + fEtaShift) ){
1674          return kFALSE;
1675       }
1676       if( eNeg->Eta() < (fEtaCutMin + fEtaShift) && eNeg->Eta() > (-fEtaCutMin + fEtaShift) ){
1677          return kFALSE;
1678       }
1679    }
1680
1681    if( ePos->Pt()< fSinglePtCut ||  eNeg->Pt()< fSinglePtCut){
1682       return kFALSE;
1683    }
1684
1685    if(particle->Pt()<fPtCut){
1686       return kFALSE;
1687    }
1688
1689    return kTRUE;
1690 }
1691 ///________________________________________________________________________
1692 Bool_t AliConversionCuts::UpdateCutString() {
1693    ///Update the cut string (if it has been created yet)
1694
1695    if(fCutString && fCutString->GetString().Length() == kNCuts) {
1696       fCutString->SetString(GetCutNumber());
1697    } else {
1698       return kFALSE;
1699    }
1700    return kTRUE;
1701 }
1702 ///________________________________________________________________________
1703 void AliConversionCuts::LoadReweightingHistosMCFromFile() {
1704
1705   AliInfo("Entering loading of histograms for weighting");
1706   TFile *f = TFile::Open(fPathTrFReweighting.Data());
1707   if(!f){
1708      AliError(Form("file for weighting %s not found",fPathTrFReweighting.Data()));
1709      return;
1710   }
1711   if (fNameHistoReweightingPi0.CompareTo("") != 0 && fDoReweightHistoMCPi0 ){
1712      TH1D *hReweightMCHistPi0temp = (TH1D*)f->Get(fNameHistoReweightingPi0.Data());
1713      hReweightMCHistPi0 = new TH1D(*hReweightMCHistPi0temp);
1714      hReweightMCHistPi0->SetDirectory(0);
1715      if (hReweightMCHistPi0) AliInfo(Form("%s has been loaded from %s", fNameHistoReweightingPi0.Data(),fPathTrFReweighting.Data() ));
1716      else AliWarning(Form("%s not found in %s", fNameHistoReweightingPi0.Data() ,fPathTrFReweighting.Data()));
1717   }
1718   if (fNameFitDataPi0.CompareTo("") != 0 && fDoReweightHistoMCPi0 ){
1719      TF1 *fFitDataPi0temp = (TF1*)f->Get(fNameFitDataPi0.Data());
1720      fFitDataPi0 = new TF1(*fFitDataPi0temp);
1721      if (fFitDataPi0) AliInfo(Form("%s has been loaded from %s", fNameFitDataPi0.Data(),fPathTrFReweighting.Data() ));
1722      else AliWarning(Form("%s not found in %s",fPathTrFReweighting.Data(), fNameFitDataPi0.Data() ));
1723   }
1724
1725   if (fNameHistoReweightingEta.CompareTo("") != 0 && fDoReweightHistoMCEta){
1726      TH1D *hReweightMCHistEtatemp = (TH1D*)f->Get(fNameHistoReweightingEta.Data());
1727      hReweightMCHistEta = new TH1D(*hReweightMCHistEtatemp);
1728      hReweightMCHistEta->SetDirectory(0);
1729      if (hReweightMCHistEta) AliInfo(Form("%s has been loaded from %s", fNameHistoReweightingEta.Data(),fPathTrFReweighting.Data() ));
1730      else AliWarning(Form("%s not found in %s", fNameHistoReweightingEta.Data(),fPathTrFReweighting.Data() ));
1731   }
1732
1733   if (fNameFitDataEta.CompareTo("") != 0 && fDoReweightHistoMCEta){
1734      TF1 *fFitDataEtatemp = (TF1*)f->Get(fNameFitDataEta.Data());
1735      fFitDataEta = new TF1(*fFitDataEtatemp);
1736      if (fFitDataEta) AliInfo(Form("%s has been loaded from %s", fNameFitDataEta.Data(),fPathTrFReweighting.Data() ));
1737      else AliWarning(Form("%s not found in %s", fNameFitDataEta.Data(),fPathTrFReweighting.Data() ));
1738
1739   }
1740   if (fNameHistoReweightingK0s.CompareTo("") != 0 && fDoReweightHistoMCK0s){
1741      TH1D *hReweightMCHistK0stemp = (TH1D*)f->Get(fNameHistoReweightingK0s.Data());
1742      hReweightMCHistK0s = new TH1D(*hReweightMCHistK0stemp);
1743      hReweightMCHistK0s->SetDirectory(0);
1744      if (hReweightMCHistK0s) AliInfo(Form("%s has been loaded from %s", fNameHistoReweightingK0s.Data(),fPathTrFReweighting.Data() ));
1745      else AliWarning(Form("%s not found in %s", fNameHistoReweightingK0s.Data(),fPathTrFReweighting.Data() ));
1746   }
1747
1748   if (fNameFitDataK0s.CompareTo("") != 0 && fDoReweightHistoMCK0s){
1749      TF1 *fFitDataK0stemp = (TF1*)f->Get(fNameFitDataK0s.Data());
1750      fFitDataK0s = new TF1(*fFitDataK0stemp);
1751      if (fFitDataK0s) AliInfo(Form("%s has been loaded from %s", fNameFitDataK0s.Data(),fPathTrFReweighting.Data() ));
1752      else AliWarning(Form("%s not found in %s", fNameFitDataK0s.Data(),fPathTrFReweighting.Data() ));
1753   }
1754   f->Close();
1755   delete f;
1756   
1757 }
1758
1759
1760 ///________________________________________________________________________
1761 Bool_t AliConversionCuts::InitializeCutsFromCutString(const TString analysisCutSelection ) {
1762    // Initialize Cuts from a given Cut string
1763    if(fDoReweightHistoMCPi0 || fDoReweightHistoMCEta || fDoReweightHistoMCK0s) {
1764       AliInfo("Weighting was enabled");
1765       LoadReweightingHistosMCFromFile();
1766    }
1767
1768    AliInfo(Form("Set Photoncut Number: %s",analysisCutSelection.Data()));
1769    if(analysisCutSelection.Length()!=kNCuts) {
1770       AliError(Form("Cut selection has the wrong length! size is %d, number of cuts is %d", analysisCutSelection.Length(), kNCuts));
1771       return kFALSE;
1772    }
1773    if(!analysisCutSelection.IsDigit()){
1774       AliError("Cut selection contains characters");
1775       return kFALSE;
1776    }
1777
1778    const char *cutSelection = analysisCutSelection.Data();
1779 #define ASSIGNARRAY(i)  fCuts[i] = cutSelection[i] - '0'
1780    for(Int_t ii=0;ii<kNCuts;ii++){
1781       ASSIGNARRAY(ii);
1782    }
1783
1784    // Set Individual Cuts
1785    for(Int_t ii=0;ii<kNCuts;ii++){
1786       if(!SetCut(cutIds(ii),fCuts[ii]))return kFALSE;
1787    }
1788
1789    //PrintCuts();
1790
1791    return kTRUE;
1792 }
1793 ///________________________________________________________________________
1794 Bool_t AliConversionCuts::SetCut(cutIds cutID, const Int_t value) {
1795    ///Set individual cut ID
1796
1797    switch (cutID) {
1798
1799    case kv0FinderType:
1800       if( SetV0Finder(value)) {
1801          fCuts[kv0FinderType] = value;
1802          UpdateCutString();
1803          return kTRUE;
1804       } else return kFALSE;
1805
1806    case kededxSigmaCut:
1807       if( SetTPCdEdxCutElectronLine(value)) {
1808          fCuts[kededxSigmaCut] = value;
1809          UpdateCutString();
1810          return kTRUE;
1811       } else return kFALSE;
1812
1813    case kpidedxSigmaCut:
1814       if( SetTPCdEdxCutPionLine(value)) {
1815          fCuts[kpidedxSigmaCut] = value;
1816          UpdateCutString();
1817          return kTRUE;
1818       } else return kFALSE;
1819
1820    case kpiMomdedxSigmaCut:
1821       if( SetMinMomPiondEdxCut(value)) {
1822          fCuts[kpiMomdedxSigmaCut] = value;
1823          UpdateCutString();
1824          return kTRUE;
1825       } else return kFALSE;
1826
1827    case kchi2GammaCut:
1828       if( SetChi2GammaCut(value)) {
1829          fCuts[kchi2GammaCut] = value;
1830          UpdateCutString();
1831          return kTRUE;
1832       } else return kFALSE;
1833
1834    case ksinglePtCut:
1835       if( SetSinglePtCut(value)) {
1836          fCuts[ksinglePtCut] = value;
1837          UpdateCutString();
1838          return kTRUE;
1839       } else return kFALSE;
1840
1841    case kclsTPCCut:
1842       if( SetTPCClusterCut(value)) {
1843          fCuts[kclsTPCCut] = value;
1844          UpdateCutString();
1845          return kTRUE;
1846       } else return kFALSE;
1847
1848    case ketaCut:
1849       if( SetEtaCut(value)) {
1850          fCuts[ketaCut] = value;
1851          UpdateCutString();
1852          return kTRUE;
1853       } else return kFALSE;
1854
1855    case kLowPRejectionSigmaCut:
1856       if( SetLowPRejectionCuts(value)) {
1857          fCuts[kLowPRejectionSigmaCut] = value;
1858          UpdateCutString();
1859          return kTRUE;
1860       } else return kFALSE;
1861
1862    case kQtMaxCut:
1863       if( SetQtMaxCut(value)) {
1864          fCuts[kQtMaxCut] = value;
1865          UpdateCutString();
1866          return kTRUE;
1867       } else return kFALSE;
1868
1869    case kpiMaxMomdedxSigmaCut:
1870       if( SetMaxMomPiondEdxCut(value)) {
1871          fCuts[kpiMaxMomdedxSigmaCut] = value;
1872          UpdateCutString();
1873          return kTRUE;
1874       } else return kFALSE;
1875
1876    case kRCut:
1877       if( SetRCut(value)) {
1878          fCuts[kRCut] = value;
1879          UpdateCutString();
1880          return kTRUE;
1881       } else return kFALSE;
1882
1883    case kremovePileUp:
1884       if( SetRemovePileUp(value)) {
1885          fCuts[kremovePileUp] = value;
1886          UpdateCutString();
1887          return kTRUE;
1888       } else return kFALSE;
1889
1890    case kselectV0AND:
1891       if( SetSelectSpecialTrigger(value)) {
1892          fCuts[kselectV0AND] = value;
1893          UpdateCutString();
1894          return kTRUE;
1895       } else return kFALSE;
1896
1897    case kmultiplicityMethod:
1898       if( SetMultiplicityMethod(value)) {
1899          fCuts[kmultiplicityMethod] = value;
1900          UpdateCutString();
1901          return kTRUE;
1902       } else return kFALSE;
1903
1904    case kisHeavyIon:
1905       if( SetIsHeavyIon(value)) {
1906          fCuts[kisHeavyIon] = value;
1907          UpdateCutString();
1908          return kTRUE;
1909       } else return kFALSE;
1910
1911    case kCentralityMin:
1912       if( SetCentralityMin(value)) {
1913          fCuts[kCentralityMin] = value;
1914          UpdateCutString();
1915          return kTRUE;
1916       } else return kFALSE;
1917
1918    case kCentralityMax:
1919       if( SetCentralityMax(value)) {
1920          fCuts[kCentralityMax] = value;
1921          UpdateCutString();
1922          return kTRUE;
1923       } else return kFALSE;
1924
1925    case kTOFelectronPID:
1926       if( SetTOFElectronPIDCut(value)) {
1927          fCuts[kTOFelectronPID] = value;
1928          UpdateCutString();
1929          return kTRUE;
1930       } else return kFALSE;
1931
1932    case kdoPhotonAsymmetryCut:
1933       if( SetPhotonAsymmetryCut(value)) {
1934          fCuts[kdoPhotonAsymmetryCut] = value;
1935          UpdateCutString();
1936          return kTRUE;
1937       } else return kFALSE;
1938
1939    case kPsiPair:
1940       if( SetPsiPairCut(value)) {
1941          fCuts[kPsiPair] = value;
1942          UpdateCutString();
1943          return kTRUE;
1944       } else return kFALSE;
1945
1946    case kCosPAngle:
1947       if( SetCosPAngleCut(value)) {
1948          fCuts[kCosPAngle] = value;
1949          UpdateCutString();
1950          return kTRUE;
1951       } else return kFALSE;
1952
1953
1954    case kElecShare:
1955       if( SetSharedElectronCut(value)) {
1956          fCuts[kElecShare] = value;
1957          UpdateCutString();
1958          return kTRUE;
1959       } else return kFALSE;
1960
1961    case kToCloseV0s:
1962       if( SetToCloseV0sCut(value)) {
1963          fCuts[kToCloseV0s] = value;
1964          UpdateCutString();
1965          return kTRUE;
1966       } else return kFALSE;
1967
1968    case kExtraSignals:
1969       if( SetRejectExtraSignalsCut(value)) {
1970          fCuts[kExtraSignals] = value;
1971          UpdateCutString();
1972          return kTRUE;
1973       } else return kFALSE;
1974
1975    case kDcaRPrimVtx:
1976       if( SetDCARPhotonPrimVtxCut(value)) {
1977          fCuts[kDcaRPrimVtx] = value;
1978          UpdateCutString();
1979          return kTRUE;
1980       } else return kFALSE;
1981
1982    case kDcaZPrimVtx:
1983       if( SetDCAZPhotonPrimVtxCut(value)) {
1984          fCuts[kDcaZPrimVtx] = value;
1985          UpdateCutString();
1986          return kTRUE;
1987       } else return kFALSE;
1988
1989    case kInPlaneOutOfPlane:
1990    if( SetInPlaneOutOfPlane(value)) {
1991       fCuts[kInPlaneOutOfPlane] = value;
1992       UpdateCutString();
1993       return kTRUE;
1994    } else return kFALSE;
1995
1996
1997       
1998
1999    case kNCuts:
2000       AliError("Cut id out of range");
2001       return kFALSE;
2002    }
2003
2004    AliError("Cut id %d not recognized");
2005    return kFALSE;
2006
2007
2008 }
2009 ///________________________________________________________________________
2010 void AliConversionCuts::PrintCuts() {
2011    // Print out current Cut Selection
2012    for(Int_t ic = 0; ic < kNCuts; ic++) {
2013       printf("%-30s : %d \n", fgkCutNames[ic], fCuts[ic]);
2014    }
2015 }
2016 ///________________________________________________________________________
2017 Bool_t AliConversionCuts::SetIsHeavyIon(Int_t isHeavyIon)
2018 {   // Set Cut
2019    switch(isHeavyIon){
2020    case 0:
2021       fIsHeavyIon=0;
2022       break;
2023    case 1:
2024       fIsHeavyIon=1;
2025       fDetectorCentrality=0;
2026       break;
2027    case 2:
2028       fIsHeavyIon=1;
2029       fDetectorCentrality=1;
2030       break;
2031    case 3: //allows to select centrality 0-45% in steps of 5% for V0 Multiplicity
2032       fIsHeavyIon=1;
2033       fDetectorCentrality=0;
2034       fModCentralityClass=1;
2035       break;
2036    case 4: //allows to select centrality 45-90% in steps of 5% for V0 Multiplicity
2037       fIsHeavyIon=1;
2038       fDetectorCentrality=0;
2039       fModCentralityClass=2;
2040       break;
2041    case 5: //strict cut on v0 tracks for MC
2042       fIsHeavyIon=1;
2043       fDetectorCentrality=0;
2044       fModCentralityClass=3;
2045       break;
2046    case 6: //allows to select centrality 0-45% in steps of 5% for track mult
2047       //strict cut on v0 tracks for MC
2048       fIsHeavyIon=1;
2049       fDetectorCentrality=0;
2050       fModCentralityClass=4;
2051       break;
2052    case 7: //allows to select centrality 45-90% in steps of 5% for V0 Multiplicity
2053       //strict cut on v0 tracks for MC
2054       fIsHeavyIon=1;
2055       fDetectorCentrality=0;
2056       fModCentralityClass=5;
2057       break;
2058    case 8:
2059       fIsHeavyIon=2;
2060       fDetectorCentrality=0;
2061       break;
2062    case 9:
2063       fIsHeavyIon=2;
2064       fDetectorCentrality=1;
2065       break;
2066    default:
2067       AliError(Form("SetHeavyIon not defined %d",isHeavyIon));
2068       return kFALSE;
2069    }
2070    return kTRUE;
2071 }
2072
2073 //___________________________________________________________________
2074 Bool_t AliConversionCuts::SetCentralityMin(Int_t minCentrality)
2075 {
2076    // Set Cut
2077    if(minCentrality<0||minCentrality>9){
2078       AliError(Form("minCentrality not defined %d",minCentrality));
2079       return kFALSE;
2080    }
2081
2082    fCentralityMin=minCentrality;
2083    return kTRUE;
2084 }
2085 //___________________________________________________________________
2086 Bool_t AliConversionCuts::SetCentralityMax(Int_t maxCentrality)
2087 {
2088    // Set Cut
2089    if(maxCentrality<0||maxCentrality>9){
2090       AliError(Form("maxCentrality not defined %d",maxCentrality));
2091       return kFALSE;
2092    }
2093    fCentralityMax=maxCentrality;
2094    return kTRUE;
2095 }
2096 ///________________________________________________________________________
2097 Int_t AliConversionCuts::SetSelectSpecialTrigger(Int_t selectSpecialTrigger)
2098 {// Set Cut
2099
2100    switch(selectSpecialTrigger){
2101    case 0:
2102       fSpecialTrigger=0; // dont care
2103       break;
2104    case 1:
2105       fSpecialTrigger=1; // V0AND
2106       break;
2107    case 2:
2108       fSpecialTrigger=2; // with SDD requested
2109       break;
2110    case 3:
2111       fSpecialTrigger=3; // V0AND plus with SDD requested
2112       break;
2113    // allows to run MB & 6 other different trigger classes in parallel with the same photon cut
2114    case 4:
2115       fSpecialTrigger=4; // different trigger class as MB
2116       fTriggerSelectedManually = kTRUE;
2117       break;
2118    case 5:
2119       fSpecialTrigger=4; // different trigger class as MB
2120       fTriggerSelectedManually = kTRUE;
2121       break;
2122    case 6:
2123       fSpecialTrigger=4; // different trigger class as MB
2124       fTriggerSelectedManually = kTRUE;
2125       break;
2126    case 7:
2127       fSpecialTrigger=4; // different trigger class as MB
2128       fTriggerSelectedManually = kTRUE;
2129       break;
2130     case 8:
2131       fSpecialTrigger=4; // different trigger class as MB
2132       fTriggerSelectedManually = kTRUE;
2133       break;
2134     case 9:
2135       fSpecialTrigger=4; // different trigger class as MB
2136       fTriggerSelectedManually = kTRUE;
2137       break;
2138    default:
2139       AliError("Warning: Special Trigger Not known");
2140       return kFALSE;
2141    }
2142    return kTRUE;
2143 }
2144 ///________________________________________________________________________
2145 Bool_t AliConversionCuts::SetMultiplicityMethod(Int_t multiplicityMethod)
2146 {
2147    // Set Cut
2148    fMultiplicityMethod=multiplicityMethod;
2149
2150    // 0 Photon Multiplicity
2151    // 1 TPC Track multiplicity
2152    // 2 V0 Mult
2153    // 3 SPD Mult
2154
2155    return kTRUE;
2156 }
2157 ///________________________________________________________________________
2158 Bool_t AliConversionCuts::SetRemovePileUp(Int_t removePileUp)
2159 {// Set Cut
2160    switch(removePileUp){
2161    case 0:
2162       fRemovePileUp=kFALSE;
2163       break;
2164    case 1:
2165       fRemovePileUp=kTRUE;
2166       break;
2167    default:
2168       AliError("RemovePileUpCut not defined");
2169       return kFALSE;
2170    }
2171    return kTRUE;
2172 }
2173 ///________________________________________________________________________
2174 Bool_t AliConversionCuts::SetRejectExtraSignalsCut(Int_t extraSignal) {
2175
2176    switch(extraSignal){
2177    case 0:
2178       fRejectExtraSignals = 0;
2179       break; // No Rejection
2180    case 1:
2181       fRejectExtraSignals = 1;
2182       break; // MinBias Header
2183    case 2:
2184       fRejectExtraSignals = 2;
2185       break; // User String Array
2186    case 3:
2187       fRejectExtraSignals = 3;
2188       break; // Rejection for Gamma Correction only
2189    default:
2190       AliError(Form("Extra Signal Rejection not defined %d",extraSignal));
2191       return kFALSE;
2192    }
2193    return kTRUE;
2194 }
2195 ///________________________________________________________________________
2196 Bool_t AliConversionCuts::SetV0Finder(Int_t v0FinderType)
2197 {   // Set Cut
2198    switch (v0FinderType){
2199    case 0:  // on fly V0 finder
2200       cout << "have chosen onfly V0" << endl;
2201       fUseOnFlyV0Finder=kTRUE;
2202       break;
2203    case 1:  // offline V0 finder
2204       cout << "have chosen offline V0" << endl;
2205       fUseOnFlyV0Finder=kFALSE;
2206       break;
2207    default:
2208       AliError(Form(" v0FinderType not defined %d",v0FinderType));
2209       return kFALSE;
2210    }
2211    return kTRUE;
2212 }
2213 ///________________________________________________________________________
2214 Bool_t AliConversionCuts::SetEtaCut(Int_t etaCut)
2215 {   // Set Cut
2216
2217    //Set Standard LineCutZValues
2218    fLineCutZValueMin = -2;
2219    fLineCutZValue = 7.;
2220
2221    switch(etaCut){
2222    case 0: // 0.9
2223       fEtaCut     = 0.9;
2224       fLineCutZRSlope = tan(2*atan(exp(-fEtaCut)));
2225       fEtaCutMin     = -0.1;
2226       fLineCutZRSlopeMin = 0.;
2227       break;
2228    case 1:  // 0.6  // changed from 1.2 to 0.6 on 2013.06.10
2229       fEtaCut     = 0.6;
2230       fLineCutZRSlope = tan(2*atan(exp(-fEtaCut)));
2231       fEtaCutMin     = -0.1;
2232       fLineCutZRSlopeMin = 0.;
2233       break;
2234    case 2:  // 1.4
2235       fEtaCut     = 1.4;
2236       fLineCutZRSlope = tan(2*atan(exp(-fEtaCut)));
2237       fEtaCutMin     = -0.1;
2238       fLineCutZRSlopeMin = 0.;
2239       break;
2240    case 3: // 0.65
2241       fEtaCut     = 0.65;
2242       fLineCutZRSlope = tan(2*atan(exp(-fEtaCut)));
2243       fEtaCutMin     = -0.1;
2244       fLineCutZRSlopeMin = 0.;
2245       break;
2246    case 4: // 0.75
2247       fEtaCut     = 0.75;
2248       fLineCutZRSlope = tan(2*atan(exp(-fEtaCut)));
2249       fEtaCutMin     = -0.1;
2250       fLineCutZRSlopeMin = 0.;
2251       break;
2252    case 5: // 0.5
2253       fEtaCut     = 0.5;
2254       fLineCutZRSlope = tan(2*atan(exp(-fEtaCut)));
2255       fEtaCutMin     = -0.1;
2256       fLineCutZRSlopeMin = 0.;
2257       break;
2258    case 6: // 5.
2259       fEtaCut     = 5.;
2260       fLineCutZRSlope = tan(2*atan(exp(-fEtaCut)));
2261       fEtaCutMin     = -0.1;
2262       fLineCutZRSlopeMin = 0.;
2263       break;
2264    case 7:
2265       if (fIsHeavyIon==1){
2266          fEtaCut     = 0.7;
2267          fLineCutZRSlope = tan(2*atan(exp(-fEtaCut)));
2268          fEtaCutMin     = -0.1;
2269          fLineCutZRSlopeMin = 0.;
2270          break;
2271       } else {   
2272          fEtaCut     = 0.3;
2273          fLineCutZRSlope = tan(2*atan(exp(-fEtaCut)));
2274          fEtaCutMin     = -0.1;
2275          fLineCutZRSlopeMin = 0.;
2276          break;
2277       }
2278    // case 8: // 0.1 - 0.8
2279    //    fEtaCut     = 0.9;
2280    //    fLineCutZRSlope = tan(2*atan(exp(-fEtaCut)));
2281    //    fEtaCutMin     = 0.1;
2282    //    fLineCutZRSlopeMin = tan(2*atan(exp(-fEtaCutMin)));
2283    //    break;
2284    case 8: // 0.4
2285       fEtaCut     = 0.4;
2286       fLineCutZRSlope = tan(2*atan(exp(-fEtaCut)));
2287       fEtaCutMin     = -0.1;
2288       fLineCutZRSlopeMin = 0.;
2289       break;
2290    case 9: // 10
2291       fEtaCut     = 10;
2292       fLineCutZRSlope = tan(2*atan(exp(-fEtaCut)));
2293       fEtaCutMin     = -0.1;
2294       fLineCutZRSlopeMin = 0.;
2295       break;
2296    default:
2297       AliError(Form(" EtaCut not defined %d",etaCut));
2298       return kFALSE;
2299    }
2300    return kTRUE;
2301 }
2302 ///________________________________________________________________________
2303 Bool_t AliConversionCuts::SetRCut(Int_t RCut){
2304    // Set Cut
2305    switch(RCut){
2306    case 0:
2307       fMinR=0;
2308       fMaxR = 180.;
2309       break;
2310    case 1:
2311       fMinR=2.8;
2312       fMaxR = 180.;
2313       break;
2314    case 2:
2315       fMinR=5.;
2316       fMaxR = 180.;
2317       break;
2318    case 3:
2319       fMaxR = 70.;
2320       fMinR = 10.;
2321       break;
2322    case 4:
2323       fMaxR = 70.;
2324       fMinR = 5.;
2325       break;
2326       // High purity cuts for PbPb (remove first layers of material)
2327    case 5:
2328       fMaxR = 180.;
2329       fMinR = 10.;
2330       break;
2331    case 6:
2332       fMaxR = 180.;
2333       fMinR = 20.;
2334       break;
2335    case 7:
2336       fMaxR = 180.;
2337       fMinR = 26.;
2338       break;
2339    case 8:
2340       fMaxR = 180.;
2341       fMinR = 12.5;
2342       break;
2343    case 9:
2344       fMaxR = 180.;
2345       fMinR = 7.5;
2346       break;
2347
2348    default:
2349       AliError("RCut not defined");
2350       return kFALSE;
2351    }
2352    return kTRUE;
2353 }
2354 ///________________________________________________________________________
2355 Bool_t AliConversionCuts::SetSinglePtCut(Int_t singlePtCut)
2356 {   // Set Cut
2357    switch(singlePtCut){
2358    case 0: // 0.050 GeV
2359       fSinglePtCut = 0.050;
2360       break;
2361    case 1:  // 0.100 GeV
2362       fSinglePtCut = 0.100;
2363       break;
2364    case 2:  // 0.150 GeV
2365       fSinglePtCut = 0.150;
2366       break;
2367    case 3:  // 0.200 GeV
2368       fSinglePtCut = 0.200;
2369       break;
2370    case 4:  // 0.075 GeV
2371       fSinglePtCut = 0.075;
2372       break;
2373    case 5:  // 0.125 GeV
2374       fSinglePtCut = 0.125;
2375       break;
2376    case 6:  // 0.04 GeV
2377       fSinglePtCut = 0.040;
2378       break;
2379    case 7:  // 0.0 GeV
2380       fSinglePtCut = 0.0;
2381       break;
2382    default:
2383       AliError(Form("singlePtCut not defined %d",singlePtCut));
2384       return kFALSE;
2385    }
2386    return kTRUE;
2387 }
2388 ///________________________________________________________________________
2389 Bool_t AliConversionCuts::SetTPCClusterCut(Int_t clsTPCCut)
2390 {   // Set Cut
2391    switch(clsTPCCut){
2392    case 0: // 0
2393       fMinClsTPC= 0.;
2394       break;
2395    case 1:  // 60
2396       fMinClsTPC= 60.;
2397       break;
2398    case 2:  // 80
2399       fMinClsTPC= 80.;
2400       break;
2401    case 3:  // 100
2402       fMinClsTPC= 100.;
2403       break;
2404    case 4:  // 95% of findable clusters
2405       fMinClsTPCToF= 0.95;
2406       fUseCorrectedTPCClsInfo=1;
2407       break;
2408    case 5:  // 0% of findable clusters
2409       fMinClsTPCToF= 0.0;
2410       fUseCorrectedTPCClsInfo=1;
2411       break;
2412    case 6:  // 70% of findable clusters
2413       fMinClsTPCToF= 0.7;
2414       fUseCorrectedTPCClsInfo=1;
2415       break;
2416    case 7:  // 0% of findable clusters
2417       fMinClsTPCToF= 0.35;
2418       fUseCorrectedTPCClsInfo=0;
2419       break;
2420    case 8:
2421       fMinClsTPCToF= 0.35;
2422       fUseCorrectedTPCClsInfo=1;
2423       break;
2424    case 9:
2425       fMinClsTPCToF= 0.6;
2426       fUseCorrectedTPCClsInfo=1;
2427       break;
2428    default:
2429       AliError(Form("Warning: clsTPCCut not defined %d",clsTPCCut));
2430       return kFALSE;
2431    }
2432    return kTRUE;
2433 }
2434 ///________________________________________________________________________
2435 Bool_t AliConversionCuts::SetTPCdEdxCutElectronLine(Int_t ededxSigmaCut)
2436 {   // Set Cut
2437    switch(ededxSigmaCut){
2438    case 0: // -10,10
2439       fPIDnSigmaBelowElectronLine=-10;
2440       fPIDnSigmaAboveElectronLine=10;
2441       break;
2442    case 1: // -5,5
2443       fPIDnSigmaBelowElectronLine=-5;
2444       fPIDnSigmaAboveElectronLine=5;
2445       break;
2446    case 2: // -3,5
2447       fPIDnSigmaBelowElectronLine=-3;
2448       fPIDnSigmaAboveElectronLine=5;
2449       break;
2450    case 3: // -4,5
2451       fPIDnSigmaBelowElectronLine=-4;
2452       fPIDnSigmaAboveElectronLine=5;
2453       break;
2454    case 4: // -6,7
2455       fPIDnSigmaBelowElectronLine=-6;
2456       fPIDnSigmaAboveElectronLine=7;
2457       break;
2458    case 5: // -4,4
2459       fPIDnSigmaBelowElectronLine=-4;
2460       fPIDnSigmaAboveElectronLine=4;
2461       break;
2462    case 6: // -2.5,4
2463       fPIDnSigmaBelowElectronLine=-2.5;
2464       fPIDnSigmaAboveElectronLine=4;
2465       break;
2466    case 7: // -2,3.5
2467       fPIDnSigmaBelowElectronLine=-2;
2468       fPIDnSigmaAboveElectronLine=3.5;
2469       break;
2470    default:
2471       AliError("TPCdEdxCutElectronLine not defined");
2472       return kFALSE;
2473
2474    }
2475    return kTRUE;
2476 }
2477 ///________________________________________________________________________
2478 Bool_t AliConversionCuts::SetTPCdEdxCutPionLine(Int_t pidedxSigmaCut)
2479 {   // Set Cut
2480
2481    switch(pidedxSigmaCut){
2482    case 0:  // -10
2483       fPIDnSigmaAbovePionLine=-10;
2484       fPIDnSigmaAbovePionLineHighPt=-10;
2485       break;
2486    case 1:   // 0
2487       fPIDnSigmaAbovePionLine=0;
2488       fPIDnSigmaAbovePionLineHighPt=-10;
2489       break;
2490    case 2:  // 1
2491       fPIDnSigmaAbovePionLine=1;
2492       fPIDnSigmaAbovePionLineHighPt=-10;
2493       break;
2494    case 3:  // 1
2495       fPIDnSigmaAbovePionLine=2.5;
2496       fPIDnSigmaAbovePionLineHighPt=-10;
2497       break;
2498    case 4:  // 1
2499       fPIDnSigmaAbovePionLine=0.5;
2500       fPIDnSigmaAbovePionLineHighPt=-10;
2501       break;
2502    case 5:  // 1
2503       fPIDnSigmaAbovePionLine=2.;
2504       fPIDnSigmaAbovePionLineHighPt=-10;
2505       break;
2506    case 6:  // 1
2507       fPIDnSigmaAbovePionLine=2.;
2508       fPIDnSigmaAbovePionLineHighPt=0.5;
2509       break;
2510    case 7:  // 1
2511       fPIDnSigmaAbovePionLine=3.5;
2512       fPIDnSigmaAbovePionLineHighPt=-10;
2513       break;
2514    case 8:  // 1
2515       fPIDnSigmaAbovePionLine=2.;
2516       fPIDnSigmaAbovePionLineHighPt=1.;
2517       break;
2518    case 9:
2519       fPIDnSigmaAbovePionLine=3.0; // We need a bit less tight cut on dE/dx
2520       fPIDnSigmaAbovePionLineHighPt=-10;
2521       break;
2522    default:
2523       AliError(Form("Warning: pidedxSigmaCut not defined %d",pidedxSigmaCut));
2524       return kFALSE;
2525    }
2526    return kTRUE;
2527 }
2528 ///________________________________________________________________________
2529 Bool_t AliConversionCuts::SetMinMomPiondEdxCut(Int_t piMomdedxSigmaCut)
2530 {   // Set Cut
2531    switch(piMomdedxSigmaCut){
2532    case 0:  // 0.5 GeV
2533       fPIDMinPnSigmaAbovePionLine=0.5;
2534       break;
2535    case 1:  // 1. GeV
2536       fPIDMinPnSigmaAbovePionLine=1.;
2537       break;
2538    case 2:  // 1.5 GeV
2539       fPIDMinPnSigmaAbovePionLine=1.5;
2540       break;
2541    case 3:  // 20.0 GeV
2542       fPIDMinPnSigmaAbovePionLine=20.;
2543       break;
2544    case 4:  // 50.0 GeV
2545       fPIDMinPnSigmaAbovePionLine=50.;
2546       break;
2547    case 5:  // 0.3 GeV
2548       fPIDMinPnSigmaAbovePionLine=0.3;
2549       break;
2550    case 6:  // 0.25 GeV
2551       fPIDMinPnSigmaAbovePionLine=0.25;
2552       break;
2553    case 7:  // 0.4 GeV
2554       fPIDMinPnSigmaAbovePionLine=0.4;
2555       break;
2556    case 8:  // 0.2 GeV
2557       fPIDMinPnSigmaAbovePionLine=0.2;
2558       break;
2559    default:
2560       AliError(Form("piMomdedxSigmaCut not defined %d",piMomdedxSigmaCut));
2561       return kFALSE;
2562    }
2563    return kTRUE;
2564 }
2565 ///________________________________________________________________________
2566 Bool_t AliConversionCuts::SetMaxMomPiondEdxCut(Int_t piMaxMomdedxSigmaCut)
2567 {   // Set Cut
2568    switch(piMaxMomdedxSigmaCut){
2569    case 0:  // 100. GeV
2570       fPIDMaxPnSigmaAbovePionLine=100.;
2571       break;
2572    case 1:  // 5. GeV
2573       fPIDMaxPnSigmaAbovePionLine=5.;
2574       break;
2575    case 2:  // 4. GeV
2576       fPIDMaxPnSigmaAbovePionLine=4.;
2577       break;
2578    case 3:  // 3.5 GeV
2579       fPIDMaxPnSigmaAbovePionLine=3.5;
2580       break;
2581    case 4:  // 3. GeV
2582       fPIDMaxPnSigmaAbovePionLine=3.;
2583       break;
2584    case 5:  // 7. GeV
2585       fPIDMaxPnSigmaAbovePionLine=7.;
2586       break;
2587    default:
2588       AliError(Form("piMaxMomdedxSigmaCut not defined %d",piMaxMomdedxSigmaCut));
2589       return kFALSE;
2590    }
2591    return kTRUE;
2592 }
2593 ///________________________________________________________________________
2594 Bool_t AliConversionCuts::SetLowPRejectionCuts(Int_t LowPRejectionSigmaCut)
2595 {   // Set Cut
2596    switch(LowPRejectionSigmaCut){
2597    case 0:  //
2598       fPIDnSigmaAtLowPAroundKaonLine=0;
2599       fPIDnSigmaAtLowPAroundProtonLine=0;
2600       fPIDnSigmaAtLowPAroundPionLine=0;
2601       fDoKaonRejectionLowP = kFALSE;
2602       fDoProtonRejectionLowP = kFALSE;
2603       fDoPionRejectionLowP = kFALSE;
2604       fPIDMinPPionRejectionLowP = fPIDMinPnSigmaAbovePionLine;
2605       break;
2606    case 1:  //
2607       fPIDnSigmaAtLowPAroundKaonLine=0.5;
2608       fPIDnSigmaAtLowPAroundProtonLine=0.5;
2609       fPIDnSigmaAtLowPAroundPionLine=0.5;
2610       fDoKaonRejectionLowP = kTRUE;
2611       fDoProtonRejectionLowP = kTRUE;
2612       fDoPionRejectionLowP = kTRUE;
2613       fPIDMinPPionRejectionLowP = fPIDMinPnSigmaAbovePionLine;
2614       break;
2615    case 2:  //
2616       fPIDnSigmaAtLowPAroundKaonLine=1;
2617       fPIDnSigmaAtLowPAroundProtonLine=1;
2618       fPIDnSigmaAtLowPAroundPionLine=1;
2619       fDoKaonRejectionLowP = kTRUE;
2620       fDoProtonRejectionLowP = kTRUE;
2621       fDoPionRejectionLowP = kTRUE;
2622       fPIDMinPPionRejectionLowP = fPIDMinPnSigmaAbovePionLine;
2623       break;
2624    case 3:  //
2625       fPIDnSigmaAtLowPAroundKaonLine=2.;
2626       fPIDnSigmaAtLowPAroundProtonLine=2.;
2627       fPIDnSigmaAtLowPAroundPionLine=2.;
2628       fDoKaonRejectionLowP = kTRUE;
2629       fDoProtonRejectionLowP = kTRUE;
2630       fDoPionRejectionLowP = kTRUE;
2631       fPIDMinPPionRejectionLowP = fPIDMinPnSigmaAbovePionLine;
2632       break;
2633    case 4:  //
2634       fPIDnSigmaAtLowPAroundKaonLine=0.;
2635       fPIDnSigmaAtLowPAroundProtonLine=0.;
2636       fPIDnSigmaAtLowPAroundPionLine=1;
2637       fDoKaonRejectionLowP = kFALSE;
2638       fDoProtonRejectionLowP = kFALSE;
2639       fDoPionRejectionLowP = kTRUE;
2640       fPIDMinPPionRejectionLowP = fPIDMinPnSigmaAbovePionLine;
2641       break;
2642    case 5:  //
2643       fPIDnSigmaAtLowPAroundKaonLine=0.;
2644       fPIDnSigmaAtLowPAroundProtonLine=0.;
2645       fPIDnSigmaAtLowPAroundPionLine=1.5;
2646       fDoKaonRejectionLowP = kFALSE;
2647       fDoProtonRejectionLowP = kFALSE;
2648       fDoPionRejectionLowP = kTRUE;
2649       fPIDMinPPionRejectionLowP = fPIDMinPnSigmaAbovePionLine;
2650       break;
2651    case 6:  //
2652       fPIDnSigmaAtLowPAroundKaonLine=0.;
2653       fPIDnSigmaAtLowPAroundProtonLine=0.;
2654       fPIDnSigmaAtLowPAroundPionLine=2.;
2655       fDoKaonRejectionLowP = kFALSE;
2656       fDoProtonRejectionLowP = kFALSE;
2657       fDoPionRejectionLowP = kTRUE;
2658       fPIDMinPPionRejectionLowP = fPIDMinPnSigmaAbovePionLine;
2659       break;
2660    case 7:  //
2661       fPIDnSigmaAtLowPAroundKaonLine=0.;
2662       fPIDnSigmaAtLowPAroundProtonLine=0.;
2663       fPIDnSigmaAtLowPAroundPionLine=0.5;
2664       fDoKaonRejectionLowP = kFALSE;
2665       fDoProtonRejectionLowP = kFALSE;
2666       fDoPionRejectionLowP = kTRUE;
2667       fPIDMinPPionRejectionLowP = fPIDMinPnSigmaAbovePionLine;
2668       break;
2669    default:
2670       AliError(Form("LowPRejectionSigmaCut not defined %d",LowPRejectionSigmaCut));
2671       return kFALSE;
2672    }
2673    return kTRUE;
2674 }
2675 ///________________________________________________________________________
2676 Bool_t AliConversionCuts::SetTOFElectronPIDCut(Int_t TOFelectronPID){
2677    // Set Cut
2678    switch(TOFelectronPID){
2679    case 0: // no cut
2680       fUseTOFpid = kFALSE;
2681       fTofPIDnSigmaBelowElectronLine=-100;
2682       fTofPIDnSigmaAboveElectronLine=100;
2683       break;
2684    case 1: // -7,7
2685       fUseTOFpid = kTRUE;
2686       fTofPIDnSigmaBelowElectronLine=-7;
2687       fTofPIDnSigmaAboveElectronLine=7;
2688       break;
2689    case 2: // -5,5
2690       fUseTOFpid = kTRUE;
2691       fTofPIDnSigmaBelowElectronLine=-5;
2692       fTofPIDnSigmaAboveElectronLine=5;
2693       break;
2694    case 3: // -3,5
2695       fUseTOFpid = kTRUE;
2696       fTofPIDnSigmaBelowElectronLine=-3;
2697       fTofPIDnSigmaAboveElectronLine=5;
2698       break;
2699    case 4: // -2,3
2700       fUseTOFpid = kTRUE;
2701       fTofPIDnSigmaBelowElectronLine=-2;
2702       fTofPIDnSigmaAboveElectronLine=3;
2703       break;
2704    case 5: // -3,3
2705       fUseTOFpid = kTRUE;
2706       fTofPIDnSigmaBelowElectronLine=-3;
2707       fTofPIDnSigmaAboveElectronLine=3;
2708       break;
2709    default:
2710       AliError(Form("TOFElectronCut not defined %d",TOFelectronPID));
2711       return kFALSE;
2712    }
2713    return kTRUE;
2714 }
2715 ///________________________________________________________________________
2716 Bool_t AliConversionCuts::SetQtMaxCut(Int_t QtMaxCut)
2717 {   // Set Cut
2718    switch(QtMaxCut){
2719    case 0: //
2720       fQtMax=1.;
2721       fDoQtGammaSelection=kFALSE;
2722       fDo2DQt=kFALSE;
2723       break;
2724    case 1:
2725       fQtMax=0.1;
2726       fDo2DQt=kFALSE;
2727       break;
2728    case 2:
2729       fQtMax=0.07;
2730       fDo2DQt=kFALSE;
2731       break;
2732    case 3:
2733       fQtMax=0.05;
2734       fDo2DQt=kFALSE;
2735       break;
2736    case 4:
2737       fQtMax=0.03;
2738       fDo2DQt=kFALSE;
2739       break;
2740    case 5:
2741       fQtMax=0.02;
2742       fDo2DQt=kFALSE;
2743       break;
2744    case 6:
2745       fQtMax=0.02;
2746       fDo2DQt=kTRUE;
2747       break;
2748    case 7:
2749       fQtMax=0.15;
2750       fDo2DQt=kFALSE;
2751       break;
2752    case 8:
2753       fQtMax=0.05;
2754       fDo2DQt=kTRUE;
2755       break;   
2756    case 9:
2757       fQtMax=0.03;
2758       fDo2DQt=kTRUE;
2759       break;      
2760    default:
2761       AliError(Form("Warning: QtMaxCut not defined %d",QtMaxCut));
2762       return kFALSE;
2763    }
2764    return kTRUE;
2765 }
2766 ///________________________________________________________________________
2767 Bool_t AliConversionCuts::SetChi2GammaCut(Int_t chi2GammaCut)
2768 {   // Set Cut
2769
2770    switch(chi2GammaCut){
2771    case 0: // 100
2772       fChi2CutConversion = 100.;
2773       break;
2774    case 1:  // 50
2775       fChi2CutConversion = 50.;
2776       break;
2777    case 2:  // 30
2778       fChi2CutConversion = 30.;
2779       break;
2780    case 3:
2781       fChi2CutConversion = 200.;
2782       break;
2783    case 4:
2784       fChi2CutConversion = 500.;
2785       break;
2786    case 5:
2787       fChi2CutConversion = 100000.;
2788       break;
2789    case 6:
2790       fChi2CutConversion = 5.;
2791       break;
2792    case 7:
2793       fChi2CutConversion = 10.;
2794       break;
2795    case 8:
2796       fChi2CutConversion = 20.;
2797       break;
2798    case 9:
2799       fChi2CutConversion = 15.;
2800       break;
2801    default:
2802       AliError(Form("Warning: Chi2GammaCut not defined %d",chi2GammaCut));
2803       return kFALSE;
2804    }
2805    return kTRUE;
2806 }
2807 ///________________________________________________________________________
2808 Bool_t AliConversionCuts::SetPsiPairCut(Int_t psiCut) {
2809
2810    switch(psiCut) {
2811    case 0:
2812       fPsiPairCut = 10000; //
2813       break;
2814    case 1:
2815       fPsiPairCut = 0.1; //
2816       break;
2817    case 2:
2818       fPsiPairCut = 0.05; // Standard
2819       break;
2820    case 3:
2821       fPsiPairCut = 0.035; //
2822       break;
2823    case 4:
2824       fPsiPairCut = 0.2; //
2825       break;   
2826    case 5:
2827       fPsiPairCut = 0.1; //
2828       fDo2DPsiPairChi2 = kTRUE;
2829       break;
2830    case 6:
2831       fPsiPairCut = 0.05; //
2832       fDo2DPsiPairChi2 = kTRUE;
2833       break;
2834    case 7:
2835       fPsiPairCut = 0.035; //
2836       fDo2DPsiPairChi2 = kTRUE;
2837       break;
2838    case 8:
2839       fPsiPairCut = 0.2; //
2840       fDo2DPsiPairChi2 = kTRUE; //
2841       break;
2842    case 9:
2843       fPsiPairCut = 0.5; //
2844       break;
2845    default:
2846       AliError(Form("PsiPairCut not defined %d",psiCut));
2847       return kFALSE;
2848    }
2849
2850    return kTRUE;
2851 }
2852 ///________________________________________________________________________
2853 Bool_t AliConversionCuts::SetPhotonAsymmetryCut(Int_t doPhotonAsymmetryCut){
2854    // Set Cut
2855    switch(doPhotonAsymmetryCut){
2856    case 0:
2857       fDoPhotonAsymmetryCut=0;
2858       fMinPPhotonAsymmetryCut=100.;
2859       fMinPhotonAsymmetry=0.;
2860       break;
2861    case 1:
2862       fDoPhotonAsymmetryCut=1;
2863       fMinPPhotonAsymmetryCut=3.5;
2864       fMinPhotonAsymmetry=0.04;
2865       break;
2866    case 2:
2867       fDoPhotonAsymmetryCut=1;
2868       fMinPPhotonAsymmetryCut=3.5;
2869       fMinPhotonAsymmetry=0.06;
2870       break;
2871    case 3:
2872       fDoPhotonAsymmetryCut=1;
2873       fMinPPhotonAsymmetryCut=0.0;
2874       fMinPhotonAsymmetry=0.05;
2875       break; 
2876    default:
2877       AliError(Form("PhotonAsymmetryCut not defined %d",doPhotonAsymmetryCut));
2878       return kFALSE;
2879    }
2880    fCuts[kdoPhotonAsymmetryCut]=doPhotonAsymmetryCut;
2881    return kTRUE;
2882 }
2883 ///________________________________________________________________________
2884 Bool_t AliConversionCuts::SetCosPAngleCut(Int_t cosCut) {
2885
2886    switch(cosCut){
2887    case 0:
2888       fCosPAngleCut = -1; 
2889       break;
2890    case 1:
2891       fCosPAngleCut = 0; 
2892       break;
2893    case 2:
2894       fCosPAngleCut = 0.5; 
2895       break;
2896    case 3:
2897       fCosPAngleCut = 0.75; 
2898       break;
2899    case 4:
2900       fCosPAngleCut = 0.85; 
2901       break;
2902    case 5:
2903       fCosPAngleCut = 0.88; 
2904       break;
2905    case 6:
2906       fCosPAngleCut = 0.9;
2907       break;
2908    case 7:
2909       fCosPAngleCut = 0.95;
2910       break;
2911    default:
2912       AliError(Form("Cosine Pointing Angle cut not defined %d",cosCut));
2913       return kFALSE;
2914    }
2915
2916    return kTRUE;
2917 }
2918 ///________________________________________________________________________
2919 Bool_t AliConversionCuts::SetSharedElectronCut(Int_t sharedElec) {
2920
2921    switch(sharedElec){
2922    case 0:
2923       fDoSharedElecCut = kFALSE;
2924       break;
2925    case 1:
2926       fDoSharedElecCut = kTRUE;
2927       break;
2928    default:
2929       AliError(Form("Shared Electron Cut not defined %d",sharedElec));
2930       return kFALSE;
2931    }
2932
2933    return kTRUE;
2934 }
2935 ///________________________________________________________________________
2936 Bool_t AliConversionCuts::SetToCloseV0sCut(Int_t toClose) {
2937
2938    switch(toClose){
2939    case 0:
2940       fDoToCloseV0sCut = kFALSE;
2941       fminV0Dist = 250;
2942       break;
2943    case 1:
2944       fDoToCloseV0sCut = kTRUE;
2945       fminV0Dist = 1;
2946       break;
2947    case 2:
2948       fDoToCloseV0sCut = kTRUE;
2949       fminV0Dist = 2;
2950       break;
2951    case 3:
2952       fDoToCloseV0sCut = kTRUE;
2953       fminV0Dist = 3;
2954       break;
2955    default:
2956       AliError(Form("Shared Electron Cut not defined %d",toClose));
2957       return kFALSE;
2958    }
2959    return kTRUE;
2960 }
2961 ///________________________________________________________________________
2962 Bool_t AliConversionCuts::SetTRDElectronCut(Int_t TRDElectronCut)
2963 {   // Set Cut
2964    switch(TRDElectronCut){
2965    case 0:
2966       fDoTRDPID=kFALSE;
2967       break;
2968    case 1:
2969       fDoTRDPID=kTRUE;
2970       fPIDTRDEfficiency=0.1;
2971       break;
2972    case 8:
2973       fDoTRDPID=kTRUE;
2974       fPIDTRDEfficiency=0.8;
2975       break;
2976    case 9:
2977       fDoTRDPID=kTRUE;
2978       fPIDTRDEfficiency=0.9;
2979       break;
2980    default:
2981       AliError(Form("TRDElectronCut not defined %d",TRDElectronCut));
2982       return kFALSE;
2983    }
2984
2985    return kTRUE;
2986 }
2987
2988 ///________________________________________________________________________
2989 Bool_t AliConversionCuts::SetDCAZPhotonPrimVtxCut(Int_t DCAZPhotonPrimVtx){
2990    // Set Cut
2991    switch(DCAZPhotonPrimVtx){
2992    case 0:  //
2993       fDCAZPrimVtxCut   = 1000;
2994       break;
2995    case 1:  //
2996       fDCAZPrimVtxCut   = 10;
2997       break;
2998    case 2:  //
2999       fDCAZPrimVtxCut   = 5;
3000       break;
3001    case 3:  //
3002       fDCAZPrimVtxCut   = 4;
3003       break;
3004    case 4:  //
3005       fDCAZPrimVtxCut   = 3;
3006       break;
3007    case 5:  //
3008       fDCAZPrimVtxCut   = 2.5;
3009       break;
3010    case 6:  //
3011       fDCAZPrimVtxCut   = 2;
3012       break;
3013    case 7:  //
3014       fDCAZPrimVtxCut   = 1.5;
3015       break;
3016    case 8:  //
3017       fDCAZPrimVtxCut   = 1;
3018       break;
3019    case 9:  //
3020       fDCAZPrimVtxCut   = 0.5;
3021       break;
3022    default:
3023       cout<<"Warning: DCAZPhotonPrimVtx not defined "<<DCAZPhotonPrimVtx<<endl;
3024       return kFALSE;
3025    }
3026    return kTRUE;
3027 }
3028
3029 ///________________________________________________________________________
3030 Bool_t AliConversionCuts::SetDCARPhotonPrimVtxCut(Int_t DCARPhotonPrimVtx){
3031    // Set Cut
3032    switch(DCARPhotonPrimVtx){
3033    case 0:  //
3034       fDCARPrimVtxCut   = 1000;
3035       break;
3036    case 1:  //
3037       fDCARPrimVtxCut   = 10;
3038       break;
3039    case 2:  //
3040       fDCARPrimVtxCut   = 5;
3041       break;
3042    case 3:  //
3043       fDCARPrimVtxCut   = 4;
3044       break;
3045    case 4:  //
3046       fDCARPrimVtxCut   = 3;
3047       break;
3048    case 5:  //
3049       fDCARPrimVtxCut   = 2.5;
3050       break;
3051    case 6:  //
3052       fDCARPrimVtxCut   = 2;
3053       break;
3054    case 7:  //
3055       fDCARPrimVtxCut   = 1.5;
3056       break;
3057    case 8:  //
3058       fDCARPrimVtxCut   = 1;
3059       break;
3060    case 9:  //
3061       fDCARPrimVtxCut   = 0.5;
3062       break;
3063    default:
3064       cout<<"Warning: DCARPhotonPrimVtx not defined "<<DCARPhotonPrimVtx<<endl;
3065       return kFALSE;
3066    }
3067    return kTRUE;
3068 }
3069
3070 ///________________________________________________________________________
3071 Bool_t AliConversionCuts::SetInPlaneOutOfPlane(Int_t inOutPlane){
3072    // Set Cut
3073    switch(inOutPlane){
3074    case 0:  //
3075       fInPlaneOutOfPlane = 0; // No Event Plane
3076       break;
3077    case 1:  //
3078       fInPlaneOutOfPlane = 1; // In-Plane
3079       break;
3080    case 2:  //
3081       fInPlaneOutOfPlane = 2; // Out-Of-Plane
3082       break;
3083    default:
3084       cout<<"Warning: In-Plane or Out-Of-Plane not defined "<<inOutPlane<<endl;
3085       return kFALSE;
3086    }
3087    return kTRUE;
3088 }
3089
3090
3091 //-------------------------------------------------------------
3092 Double_t AliConversionCuts::GetCentrality(AliVEvent *event)
3093 {   // Get Event Centrality
3094
3095    AliESDEvent *esdEvent=dynamic_cast<AliESDEvent*>(event);
3096    if(esdEvent){
3097       AliCentrality *fESDCentrality=(AliCentrality*)esdEvent->GetCentrality();
3098
3099       if(fDetectorCentrality==0){
3100          if (fIsHeavyIon==2){
3101             return fESDCentrality->GetCentralityPercentile("V0A"); // default for pPb
3102          } else{
3103             return fESDCentrality->GetCentralityPercentile("V0M"); // default
3104          }
3105       }
3106       if(fDetectorCentrality==1){
3107          return fESDCentrality->GetCentralityPercentile("CL1");
3108       }
3109    }
3110
3111    AliAODEvent *aodEvent=dynamic_cast<AliAODEvent*>(event);
3112    if(aodEvent){
3113       if(aodEvent->GetHeader()){return aodEvent->GetHeader()->GetCentrality();}
3114    }
3115
3116    return -1;
3117 }
3118 //-------------------------------------------------------------
3119 Bool_t AliConversionCuts::IsCentralitySelected(AliVEvent *event, AliVEvent *fMCEvent)
3120 {   // Centrality Selection
3121    if(!fIsHeavyIon)return kTRUE;
3122
3123    if(fCentralityMin == fCentralityMax ) return kTRUE;//0-100%
3124    else if(fCentralityMax==0) fCentralityMax=10; //CentralityRange = fCentralityMin-100%
3125
3126    Double_t centrality=GetCentrality(event);
3127    if(centrality<0)return kFALSE;
3128
3129    Int_t centralityC=0;
3130    if (fModCentralityClass == 0){
3131       centralityC= Int_t(centrality/10);
3132       if(centralityC >= fCentralityMin && centralityC < fCentralityMax)
3133          return kTRUE;
3134       else return kFALSE;
3135    }
3136    else if (fModCentralityClass ==1){
3137       centralityC= Int_t(centrality);
3138       if(centralityC >= fCentralityMin*5 && centralityC < fCentralityMax*5){
3139          return kTRUE;
3140       } else return kFALSE;
3141    }
3142    else if (fModCentralityClass ==2){
3143       centralityC= Int_t(centrality);
3144       if(centralityC >= ((fCentralityMin*5)+45) && centralityC < ((fCentralityMax*5)+45))
3145          return kTRUE;
3146       else return kFALSE;
3147    }
3148
3149    Int_t nprimaryTracks = ((AliV0ReaderV1*)AliAnalysisManager::GetAnalysisManager()->GetTask("V0ReaderV1"))->GetNumberOfPrimaryTracks();
3150    Int_t PrimaryTracks10[10][2] =
3151       {
3152          {9999,9999}, //  0
3153          {1210,2150}, // 10
3154          { 817,1435}, // 20
3155          { 536, 930}, // 30
3156          { 337, 570}, // 40
3157          { 197, 327}, // 50
3158          { 106, 173}, // 60
3159          {  51,  81}, // 70
3160          {  21,  34}, // 80
3161          {   0,   0}  // 90
3162       };
3163    Int_t PrimaryTracks5a[10][2] =
3164       {
3165          {9999,9999}, // 0
3166          {1485,2640}, // 5
3167          {1210,2150}, // 10
3168          { 995,1760}, // 15
3169          { 817,1435}, // 20
3170          { 666,1160}, // 25
3171          { 536, 930}, // 30
3172          { 428, 731}, // 35
3173          { 337, 570}, // 40
3174          { 260, 436}  // 45
3175       };
3176    Int_t PrimaryTracks5b[10][2] =
3177       {
3178          { 260, 436}, // 45
3179          { 197, 327}, // 50
3180          { 147, 239}, // 55
3181          { 106, 173}, // 60
3182          {  75, 120}, // 65
3183          {  51,  81}, // 70
3184          {  34,  53}, // 75
3185          {  21,  34}, // 80
3186          {  13,  19}, // 85
3187          {   0,   0}  // 90
3188       };
3189
3190    Int_t column = -1;
3191    if(event->IsA()==AliESDEvent::Class()) column = 0;
3192    if(event->IsA()==AliAODEvent::Class()) column = 1;
3193
3194    if (fModCentralityClass == 3){
3195       if(fMCEvent){
3196          if(nprimaryTracks > PrimaryTracks10[fCentralityMax][column] && nprimaryTracks <= PrimaryTracks10[fCentralityMin][column])
3197             return kTRUE;
3198          else return kFALSE;
3199       }
3200       else{
3201          centralityC= Int_t(centrality/10);
3202          if(centralityC >= fCentralityMin && centralityC < fCentralityMax)
3203             return kTRUE;
3204          else return kFALSE;
3205       }
3206    }
3207    else if (fModCentralityClass ==4){
3208       if(fMCEvent){
3209          if(nprimaryTracks > PrimaryTracks5a[fCentralityMax][column] && nprimaryTracks <= PrimaryTracks5a[fCentralityMin][column])
3210             return kTRUE;
3211          else return kFALSE;
3212       }
3213       else{
3214          centralityC= Int_t(centrality);
3215          if(centralityC >= fCentralityMin*5 && centralityC < fCentralityMax*5){
3216             return kTRUE;
3217          } else return kFALSE;
3218       }
3219    }
3220    else if (fModCentralityClass ==5){
3221       if(fMCEvent){
3222          if(nprimaryTracks > PrimaryTracks5b[fCentralityMax][column] && nprimaryTracks <= PrimaryTracks5b[fCentralityMin][column])
3223             return kTRUE;
3224          else return kFALSE;
3225       }
3226       else{
3227          centralityC= Int_t(centrality);
3228          if(centralityC >= ((fCentralityMin*5)+45) && centralityC < ((fCentralityMax*5)+45))
3229             return kTRUE;
3230          else return kFALSE;
3231       }
3232    }
3233
3234    return kFALSE;
3235 }
3236 ///________________________________________________________________________
3237 Bool_t AliConversionCuts::VertexZCut(AliVEvent *event){
3238    // Cut on z position of primary vertex
3239    Double_t fVertexZ=event->GetPrimaryVertex()->GetZ();
3240
3241    if(abs(fVertexZ)>fMaxVertexZ)return kFALSE;
3242
3243    if (fIsHeavyIon == 2){
3244      if(fUtils->IsFirstEventInChunk(event)) return kFALSE;
3245      if(!fUtils->IsVertexSelected2013pA(event)) return kFALSE;
3246      if(fUtils->IsPileUpEvent(event)) return kFALSE;
3247    }
3248
3249    return kTRUE;
3250 }
3251 ///________________________________________________________________________
3252
3253 Int_t AliConversionCuts::GetNumberOfContributorsVtx(AliVEvent *event){
3254    // returns number of contributors to the vertex
3255
3256    AliESDEvent *fESDEvent=dynamic_cast<AliESDEvent*>(event);
3257    if(fESDEvent){
3258       if (fESDEvent->GetPrimaryVertex() != NULL){
3259          if(fESDEvent->GetPrimaryVertex()->GetNContributors()>0) {
3260 //     cout << "accepted global" << fESDEvent->GetEventNumberInFile() << " with NCont: " << fESDEvent->GetPrimaryVertex()->GetNContributors() << endl;
3261             return fESDEvent->GetPrimaryVertex()->GetNContributors();
3262          }
3263       }
3264
3265       if(fESDEvent->GetPrimaryVertexSPD() !=NULL){
3266          if(fESDEvent->GetPrimaryVertexSPD()->GetNContributors()>0) {
3267 //     cout << "accepted SPD" << fESDEvent->GetEventNumberInFile() << " with NCont: " << fESDEvent->GetPrimaryVertexSPD()->GetNContributors() << endl;
3268             return fESDEvent->GetPrimaryVertexSPD()->GetNContributors();
3269          }  else {
3270             AliWarning(Form("Number of contributors from bad vertex type:: %s",fESDEvent->GetPrimaryVertex()->GetName()));
3271 //            cout << "rejected " << fESDEvent->GetEventNumberInFile() << endl;
3272             return 0;
3273          }
3274       }
3275    }
3276
3277    AliAODEvent *fAODEvent=dynamic_cast<AliAODEvent*>(event);
3278    if(fAODEvent){
3279       if (fAODEvent->GetPrimaryVertex() != NULL){
3280          if(fAODEvent->GetPrimaryVertex()->GetNContributors()>0) {
3281             return fAODEvent->GetPrimaryVertex()->GetNContributors();
3282          }
3283       }
3284       if(fAODEvent->GetPrimaryVertexSPD() !=NULL){
3285          if(fAODEvent->GetPrimaryVertexSPD()->GetNContributors()>0) {
3286             return fAODEvent->GetPrimaryVertexSPD()->GetNContributors();
3287          } else {
3288             AliWarning(Form("Number of contributors from bad vertex type:: %s",fAODEvent->GetPrimaryVertex()->GetName()));
3289             return 0;
3290          }
3291       }
3292    }
3293   // cout << "rejected " << fESDEvent->GetEventNumberInFile() << endl;
3294    return 0;
3295 }
3296
3297 ///________________________________________________________________________
3298
3299 Bool_t AliConversionCuts::IsTriggerSelected(AliVEvent *fInputEvent)
3300 {
3301
3302    AliInputEventHandler *fInputHandler=(AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
3303
3304    UInt_t isSelected = AliVEvent::kAny;
3305    if (fInputHandler==NULL) return kFALSE;
3306    if( fInputHandler->GetEventSelection() || fInputEvent->IsA()==AliAODEvent::Class()) {
3307       if (!fTriggerSelectedManually){
3308          if (fPreSelCut) fOfflineTriggerMask = AliVEvent::kAny;
3309          else {
3310             if (fIsHeavyIon == 1) fOfflineTriggerMask = AliVEvent::kMB | AliVEvent::kCentral | AliVEvent::kSemiCentral;
3311                else if (fIsHeavyIon == 2) fOfflineTriggerMask = AliVEvent::kINT7;
3312                else fOfflineTriggerMask = AliVEvent::kMB;
3313          }
3314       }
3315       // Get the actual offline trigger mask for the event and AND it with the
3316       // requested mask. If no mask requested select by default the event.
3317 //       if (fPreSelCut) cout << "Trigger selected from outside: "<< fTriggerSelectedManually <<"\t Offline Trigger mask for Precut: " << fOfflineTriggerMask << endl;
3318 //       else cout << "Trigger selected from outside: "<< fTriggerSelectedManually <<"\t Offline Trigger mask: " << fOfflineTriggerMask << endl;
3319
3320       if (fOfflineTriggerMask)
3321          isSelected = fOfflineTriggerMask & fInputHandler->IsEventSelected();
3322    }
3323    fIsSDDFired = !(fInputHandler->IsEventSelected() & AliVEvent::kFastOnly);
3324
3325    // Fill Histogram
3326    if(hTriggerClass){
3327       if (fIsSDDFired) hTriggerClass->Fill(33);
3328       if (fInputHandler->IsEventSelected() & AliVEvent::kMB)hTriggerClass->Fill(0);
3329       if (fInputHandler->IsEventSelected() & AliVEvent::kINT7)hTriggerClass->Fill(1);
3330       if (fInputHandler->IsEventSelected() & AliVEvent::kMUON)hTriggerClass->Fill(2);
3331       if (fInputHandler->IsEventSelected() & AliVEvent::kHighMult)hTriggerClass->Fill(3);
3332       if (fInputHandler->IsEventSelected() & AliVEvent::kEMC1)hTriggerClass->Fill(4);
3333       if (fInputHandler->IsEventSelected() & AliVEvent::kCINT5)hTriggerClass->Fill(5);
3334       if (fInputHandler->IsEventSelected() & AliVEvent::kCMUS5)hTriggerClass->Fill(6);
3335       if (fInputHandler->IsEventSelected() & AliVEvent::kMUSPB)hTriggerClass->Fill(6);
3336       if (fInputHandler->IsEventSelected() & AliVEvent::kMUSH7)hTriggerClass->Fill(7);
3337       if (fInputHandler->IsEventSelected() & AliVEvent::kMUSHPB)hTriggerClass->Fill(7);
3338       if (fInputHandler->IsEventSelected() & AliVEvent::kMUL7)hTriggerClass->Fill(8);
3339       if (fInputHandler->IsEventSelected() & AliVEvent::kMuonLikePB)hTriggerClass->Fill(8);
3340       if (fInputHandler->IsEventSelected() & AliVEvent::kMUU7)hTriggerClass->Fill(9);
3341       if (fInputHandler->IsEventSelected() & AliVEvent::kMuonUnlikePB)hTriggerClass->Fill(9);
3342       if (fInputHandler->IsEventSelected() & AliVEvent::kEMC7)hTriggerClass->Fill(10);
3343       if (fInputHandler->IsEventSelected() & AliVEvent::kEMC8)hTriggerClass->Fill(10);
3344       if (fInputHandler->IsEventSelected() & AliVEvent::kMUS7)hTriggerClass->Fill(11);
3345       if (fInputHandler->IsEventSelected() & AliVEvent::kPHI1)hTriggerClass->Fill(12);
3346       if (fInputHandler->IsEventSelected() & AliVEvent::kPHI7)hTriggerClass->Fill(13);
3347       if (fInputHandler->IsEventSelected() & AliVEvent::kPHI8)hTriggerClass->Fill(13);
3348       if (fInputHandler->IsEventSelected() & AliVEvent::kPHOSPb)hTriggerClass->Fill(13);
3349       if (fInputHandler->IsEventSelected() & AliVEvent::kEMCEJE)hTriggerClass->Fill(14);
3350       if (fInputHandler->IsEventSelected() & AliVEvent::kEMCEGA)hTriggerClass->Fill(15);
3351       if (fInputHandler->IsEventSelected() & AliVEvent::kCentral)hTriggerClass->Fill(16);
3352       if (fInputHandler->IsEventSelected() & AliVEvent::kSemiCentral)hTriggerClass->Fill(17);
3353       if (fInputHandler->IsEventSelected() & AliVEvent::kDG5)hTriggerClass->Fill(18);
3354       if (fInputHandler->IsEventSelected() & AliVEvent::kZED)hTriggerClass->Fill(19);
3355       if (fInputHandler->IsEventSelected() & AliVEvent::kSPI7)hTriggerClass->Fill(20);
3356       if (fInputHandler->IsEventSelected() & AliVEvent::kSPI)hTriggerClass->Fill(20);
3357       if (fInputHandler->IsEventSelected() & AliVEvent::kINT8)hTriggerClass->Fill(21);
3358       if (fInputHandler->IsEventSelected() & AliVEvent::kMuonSingleLowPt8)hTriggerClass->Fill(22);
3359       if (fInputHandler->IsEventSelected() & AliVEvent::kMuonSingleHighPt8)hTriggerClass->Fill(23);
3360       if (fInputHandler->IsEventSelected() & AliVEvent::kMuonLikeLowPt8)hTriggerClass->Fill(24);
3361       if (fInputHandler->IsEventSelected() & AliVEvent::kMuonUnlikeLowPt8)hTriggerClass->Fill(25);
3362       if (fInputHandler->IsEventSelected() & AliVEvent::kMuonUnlikeLowPt0)hTriggerClass->Fill(26);
3363       if (fInputHandler->IsEventSelected() & AliVEvent::kUserDefined)hTriggerClass->Fill(27);
3364       if (fInputHandler->IsEventSelected() & AliVEvent::kTRD)hTriggerClass->Fill(28);
3365       if (fInputHandler->IsEventSelected() & AliVEvent::kFastOnly)hTriggerClass->Fill(29);
3366       if (fInputHandler->IsEventSelected() & AliVEvent::kAnyINT)hTriggerClass->Fill(30);
3367       if (fInputHandler->IsEventSelected() & AliVEvent::kAny)hTriggerClass->Fill(31);
3368       if (!fInputHandler->IsEventSelected()) hTriggerClass->Fill(34);
3369    }
3370
3371    if(hTriggerClassSelected && isSelected){
3372       if (!fIsSDDFired) hTriggerClassSelected->Fill(33);
3373       if (fInputHandler->IsEventSelected() & AliVEvent::kMB)hTriggerClassSelected->Fill(0);
3374       if (fInputHandler->IsEventSelected() & AliVEvent::kINT7)hTriggerClassSelected->Fill(1);
3375       if (fInputHandler->IsEventSelected() & AliVEvent::kMUON)hTriggerClassSelected->Fill(2);
3376       if (fInputHandler->IsEventSelected() & AliVEvent::kHighMult)hTriggerClassSelected->Fill(3);
3377       if (fInputHandler->IsEventSelected() & AliVEvent::kEMC1)hTriggerClassSelected->Fill(4);
3378       if (fInputHandler->IsEventSelected() & AliVEvent::kCINT5)hTriggerClassSelected->Fill(5);
3379       if (fInputHandler->IsEventSelected() & AliVEvent::kCMUS5)hTriggerClassSelected->Fill(6);
3380       if (fInputHandler->IsEventSelected() & AliVEvent::kMUSPB)hTriggerClassSelected->Fill(6);
3381       if (fInputHandler->IsEventSelected() & AliVEvent::kMUSH7)hTriggerClassSelected->Fill(7);
3382       if (fInputHandler->IsEventSelected() & AliVEvent::kMUSHPB)hTriggerClassSelected->Fill(7);
3383       if (fInputHandler->IsEventSelected() & AliVEvent::kMUL7)hTriggerClassSelected->Fill(8);
3384       if (fInputHandler->IsEventSelected() & AliVEvent::kMuonLikePB)hTriggerClassSelected->Fill(8);
3385       if (fInputHandler->IsEventSelected() & AliVEvent::kMUU7)hTriggerClassSelected->Fill(9);
3386       if (fInputHandler->IsEventSelected() & AliVEvent::kMuonUnlikePB)hTriggerClassSelected->Fill(9);
3387       if (fInputHandler->IsEventSelected() & AliVEvent::kEMC7)hTriggerClassSelected->Fill(10);
3388       if (fInputHandler->IsEventSelected() & AliVEvent::kEMC8)hTriggerClassSelected->Fill(10);
3389       if (fInputHandler->IsEventSelected() & AliVEvent::kMUS7)hTriggerClassSelected->Fill(11);
3390       if (fInputHandler->IsEventSelected() & AliVEvent::kPHI1)hTriggerClassSelected->Fill(12);
3391       if (fInputHandler->IsEventSelected() & AliVEvent::kPHI7)hTriggerClassSelected->Fill(13);
3392       if (fInputHandler->IsEventSelected() & AliVEvent::kPHI8)hTriggerClassSelected->Fill(13);
3393       if (fInputHandler->IsEventSelected() & AliVEvent::kPHOSPb)hTriggerClassSelected->Fill(13);
3394       if (fInputHandler->IsEventSelected() & AliVEvent::kEMCEJE)hTriggerClassSelected->Fill(14);
3395       if (fInputHandler->IsEventSelected() & AliVEvent::kEMCEGA)hTriggerClassSelected->Fill(15);
3396       if (fInputHandler->IsEventSelected() & AliVEvent::kCentral)hTriggerClassSelected->Fill(16);
3397       if (fInputHandler->IsEventSelected() & AliVEvent::kSemiCentral)hTriggerClassSelected->Fill(17);
3398       if (fInputHandler->IsEventSelected() & AliVEvent::kDG5)hTriggerClassSelected->Fill(18);
3399       if (fInputHandler->IsEventSelected() & AliVEvent::kZED)hTriggerClassSelected->Fill(19);
3400       if (fInputHandler->IsEventSelected() & AliVEvent::kSPI7)hTriggerClassSelected->Fill(20);
3401       if (fInputHandler->IsEventSelected() & AliVEvent::kSPI)hTriggerClassSelected->Fill(20);
3402       if (fInputHandler->IsEventSelected() & AliVEvent::kINT8)hTriggerClassSelected->Fill(21);
3403       if (fInputHandler->IsEventSelected() & AliVEvent::kMuonSingleLowPt8)hTriggerClassSelected->Fill(22);
3404       if (fInputHandler->IsEventSelected() & AliVEvent::kMuonSingleHighPt8)hTriggerClassSelected->Fill(23);
3405       if (fInputHandler->IsEventSelected() & AliVEvent::kMuonLikeLowPt8)hTriggerClassSelected->Fill(24);
3406       if (fInputHandler->IsEventSelected() & AliVEvent::kMuonUnlikeLowPt8)hTriggerClassSelected->Fill(25);
3407       if (fInputHandler->IsEventSelected() & AliVEvent::kMuonUnlikeLowPt0)hTriggerClassSelected->Fill(26);
3408       if (fInputHandler->IsEventSelected() & AliVEvent::kUserDefined)hTriggerClassSelected->Fill(27);
3409       if (fInputHandler->IsEventSelected() & AliVEvent::kTRD)hTriggerClassSelected->Fill(28);
3410       if (fInputHandler->IsEventSelected() & AliVEvent::kFastOnly)hTriggerClassSelected->Fill(29);
3411       if (fInputHandler->IsEventSelected() & AliVEvent::kAnyINT)hTriggerClassSelected->Fill(30);
3412       if (fInputHandler->IsEventSelected() & AliVEvent::kAny)hTriggerClassSelected->Fill(31);
3413    }
3414
3415    if(!isSelected)return kFALSE;
3416
3417    return kTRUE;
3418
3419 }
3420
3421 ///________________________________________________________________________
3422 Int_t AliConversionCuts::GetFirstTPCRow(Double_t radius){
3423    // Get first TPC row
3424    Int_t firstTPCRow = 0;
3425    Double_t radiusI = 84.8;
3426    Double_t radiusO = 134.6;
3427    Double_t radiusOB = 198.;
3428    Double_t rSizeI = 0.75;
3429    Double_t rSizeO = 1.;
3430    Double_t rSizeOB = 1.5;
3431    Int_t nClsI = 63;
3432    Int_t nClsIO = 127;
3433
3434    if(radius <= radiusI){
3435       return firstTPCRow;
3436    }
3437    if(radius>radiusI && radius<=radiusO){
3438       firstTPCRow = (Int_t)((radius-radiusI)/rSizeI);
3439    }
3440    if(radius>radiusO && radius<=radiusOB){
3441       firstTPCRow = (Int_t)(nClsI+(radius-radiusO)/rSizeO);
3442    }
3443
3444    if(radius>radiusOB){
3445       firstTPCRow =(Int_t)(nClsIO+(radius-radiusOB)/rSizeOB);
3446    }
3447
3448    return firstTPCRow;
3449 }
3450
3451 Bool_t AliConversionCuts::CosinePAngleCut(const AliConversionPhotonBase * photon, AliVEvent * event) const {
3452    ///Check if passes cosine of pointing angle cut
3453    if(GetCosineOfPointingAngle(photon, event) < fCosPAngleCut){
3454       return kFALSE;
3455    }
3456    return kTRUE;
3457 }
3458
3459 Double_t AliConversionCuts::GetCosineOfPointingAngle( const AliConversionPhotonBase * photon, AliVEvent * event) const{
3460    // calculates the pointing angle of the recalculated V0
3461
3462    Double_t momV0[3] = {0,0,0};
3463    if(event->IsA()==AliESDEvent::Class()){
3464       AliESDEvent *esdEvent = dynamic_cast<AliESDEvent*>(event);
3465       if(!esdEvent) return -999;
3466       AliESDv0 *v0 = esdEvent->GetV0(photon->GetV0Index());
3467       if(!v0) return -999;
3468       v0->GetPxPyPz(momV0[0],momV0[1],momV0[2]);
3469    }
3470    if(event->IsA()==AliAODEvent::Class()){
3471       momV0[0] = photon->GetPx();
3472       momV0[1] = photon->GetPy();
3473       momV0[2] = photon->GetPz();
3474    }
3475
3476    //Double_t momV0[3] = { photon->GetPx(), photon->GetPy(), photon->GetPz() }; //momentum of the V0
3477    Double_t PosV0[3] = { photon->GetConversionX() - event->GetPrimaryVertex()->GetX(),
3478                          photon->GetConversionY() - event->GetPrimaryVertex()->GetY(),
3479                          photon->GetConversionZ() - event->GetPrimaryVertex()->GetZ() }; //Recalculated V0 Position vector
3480
3481    Double_t momV02 = momV0[0]*momV0[0] + momV0[1]*momV0[1] + momV0[2]*momV0[2];
3482    Double_t PosV02 = PosV0[0]*PosV0[0] + PosV0[1]*PosV0[1] + PosV0[2]*PosV0[2];
3483
3484
3485    Double_t cosinePointingAngle = -999;
3486    if(momV02*PosV02 > 0.0)
3487       cosinePointingAngle = (PosV0[0]*momV0[0] +  PosV0[1]*momV0[1] + PosV0[2]*momV0[2] ) / TMath::Sqrt(momV02 * PosV02);
3488
3489    return cosinePointingAngle;
3490 }
3491
3492 ///________________________________________________________________________
3493 Bool_t AliConversionCuts::PsiPairCut(const AliConversionPhotonBase * photon) const {
3494
3495    if (fDo2DPsiPairChi2){
3496       if (abs(photon->GetPsiPair()) < -fPsiPairCut/fChi2CutConversion*photon->GetChi2perNDF() + fPsiPairCut ){  
3497          return kTRUE;
3498       } else {
3499          return kFALSE;
3500       }    
3501    } else {
3502       if(abs(photon->GetPsiPair()) > fPsiPairCut){
3503          return kFALSE;}
3504       else{return kTRUE;}
3505    } 
3506 }
3507
3508 ///________________________________________________________________________
3509 TString AliConversionCuts::GetCutNumber(){
3510    // returns TString with current cut number
3511    TString a(kNCuts);
3512    for(Int_t ii=0;ii<kNCuts;ii++){
3513       a.Append(Form("%d",fCuts[ii]));
3514    }
3515    return a;
3516 }
3517
3518 ///________________________________________________________________________
3519 void AliConversionCuts::FillElectonLabelArray(AliAODConversionPhoton* photon, Int_t nV0){
3520
3521    Int_t posLabel = photon->GetTrackLabelPositive();
3522    Int_t negLabel = photon->GetTrackLabelNegative();
3523
3524    fElectronLabelArray[nV0*2] = posLabel;
3525    fElectronLabelArray[(nV0*2)+1] = negLabel;
3526 }
3527 ///________________________________________________________________________
3528 Bool_t AliConversionCuts::RejectSharedElectronV0s(AliAODConversionPhoton* photon, Int_t nV0, Int_t nV0s){
3529
3530    Int_t posLabel = photon->GetTrackLabelPositive();
3531    Int_t negLabel = photon->GetTrackLabelNegative();
3532
3533    for(Int_t i = 0; i<nV0s*2;i++){
3534       if(i==nV0*2)     continue;
3535       if(i==(nV0*2)+1) continue;
3536       if(fElectronLabelArray[i] == posLabel){
3537          return kFALSE;}
3538       if(fElectronLabelArray[i] == negLabel){
3539          return kFALSE;}
3540    }
3541
3542    return kTRUE;
3543 }
3544 ///________________________________________________________________________
3545 Bool_t AliConversionCuts::RejectToCloseV0s(AliAODConversionPhoton* photon, TList *photons, Int_t nV0){
3546
3547
3548    Double_t posX = photon->GetConversionX();
3549    Double_t posY = photon->GetConversionY();
3550    Double_t posZ = photon->GetConversionZ();
3551
3552    for(Int_t i = 0;i<photons->GetEntries();i++){
3553       if(nV0 == i) continue;
3554       AliAODConversionPhoton *photonComp = (AliAODConversionPhoton*) photons->At(i);
3555       Double_t posCompX = photonComp->GetConversionX();
3556       Double_t posCompY = photonComp->GetConversionY();
3557       Double_t posCompZ = photonComp->GetConversionZ();
3558
3559       Double_t dist = pow((posX - posCompX),2)+pow((posY - posCompY),2)+pow((posZ - posCompZ),2);
3560
3561       if(dist < fminV0Dist*fminV0Dist){
3562          if(photon->GetChi2perNDF() < photonComp->GetChi2perNDF()) return kTRUE;
3563          else {
3564             return kFALSE;}
3565       }
3566
3567    }
3568    return kTRUE;
3569 }
3570 ///________________________________________________________________________
3571 void AliConversionCuts::GetNotRejectedParticles(Int_t rejection, TList *HeaderList, AliVEvent *MCEvent){
3572
3573
3574
3575    if(fNotRejectedStart){
3576       delete[] fNotRejectedStart;
3577       fNotRejectedStart = NULL;
3578    }
3579    if(fNotRejectedEnd){
3580       delete[] fNotRejectedEnd;
3581       fNotRejectedEnd = NULL;
3582    }
3583    if(fGeneratorNames){
3584       delete[] fGeneratorNames;
3585       fGeneratorNames = NULL;
3586    }
3587
3588    if(rejection == 0) return; // No Rejection
3589
3590    AliGenCocktailEventHeader *cHeader = 0x0;
3591    AliAODMCHeader *cHeaderAOD = 0x0;
3592    Bool_t headerFound = kFALSE;
3593
3594    if(MCEvent->IsA()==AliMCEvent::Class()){
3595       cHeader = dynamic_cast<AliGenCocktailEventHeader*>(dynamic_cast<AliMCEvent*>(MCEvent)->GenEventHeader());
3596       if(cHeader) headerFound = kTRUE;
3597    }
3598    if(MCEvent->IsA()==AliAODEvent::Class()){ // MCEvent is a AODEvent in case of AOD
3599       cHeaderAOD = dynamic_cast<AliAODMCHeader*>(MCEvent->FindListObject(AliAODMCHeader::StdBranchName()));
3600       if(cHeaderAOD) headerFound = kTRUE;
3601    }
3602
3603    if(headerFound){
3604       TList *genHeaders = 0x0;
3605       if(cHeader) genHeaders = cHeader->GetHeaders();
3606       if(cHeaderAOD){
3607          genHeaders = cHeaderAOD->GetCocktailHeaders();
3608          if(genHeaders->GetEntries()==1){
3609             SetRejectExtraSignalsCut(0);
3610             return;
3611          }
3612       }
3613       AliGenEventHeader* gh = 0;
3614       fnHeaders = 0;
3615       if(rejection == 1 || rejection == 3) fnHeaders = 1; // MinBiasHeader
3616       if(rejection == 2){ // TList of Headers Names
3617          for(Int_t i = 0; i<genHeaders->GetEntries();i++){
3618             gh = (AliGenEventHeader*)genHeaders->At(i);
3619             TString GeneratorName = gh->GetName();
3620             for(Int_t j = 0; j<HeaderList->GetEntries();j++){
3621                TString GeneratorInList = ((TObjString*)HeaderList->At(j))->GetString();
3622                if(GeneratorName.CompareTo(GeneratorInList) == 0){
3623                   fnHeaders++;
3624                   continue;
3625                }
3626             }
3627          }
3628       }
3629
3630       fNotRejectedStart = new Int_t[fnHeaders];
3631       fNotRejectedEnd = new Int_t[fnHeaders];
3632       fGeneratorNames = new TString[fnHeaders];
3633
3634       if(rejection == 1 || rejection == 3){
3635          fNotRejectedStart[0] = 0;
3636          fNotRejectedEnd[0] = ((AliGenEventHeader*)genHeaders->At(0))->NProduced()-1;
3637          fGeneratorNames[0] = ((AliGenEventHeader*)genHeaders->At(0))->GetName();
3638          return;
3639       }
3640
3641       Int_t firstindex = 0;
3642       Int_t lastindex =  -1;
3643       Int_t nummer = 0;
3644       for(Int_t i = 0; i<genHeaders->GetEntries();i++){
3645          gh = (AliGenEventHeader*)genHeaders->At(i);
3646          TString GeneratorName = gh->GetName();
3647          lastindex = lastindex + gh->NProduced();
3648          for(Int_t j = 0; j<HeaderList->GetEntries();j++){
3649             TString GeneratorInList = ((TObjString*)HeaderList->At(j))->GetString();
3650             if(GeneratorName.CompareTo(GeneratorInList) == 0){
3651                fNotRejectedStart[nummer] = firstindex;
3652                fNotRejectedEnd[nummer] = lastindex;
3653                fGeneratorNames[nummer] = GeneratorName;
3654                //cout << "Number of particles produced for: " << i << "\t" << GeneratorName.Data() << "\t" << lastindex-firstindex+1 << endl;
3655                nummer++;
3656                continue;
3657             }
3658          }
3659          firstindex = firstindex + gh->NProduced();
3660       }
3661    } else { // No Cocktail Header Found
3662       fNotRejectedStart = new Int_t[1];
3663       fNotRejectedEnd = new Int_t[1];
3664
3665       fnHeaders = 1;
3666       fNotRejectedStart[0] = 0;
3667       fNotRejectedEnd[0] = static_cast<AliMCEvent*>(MCEvent)->Stack()->GetNprimary()-1;
3668       fGeneratorNames = new TString[1];
3669       fGeneratorNames[0] = "NoCocktailGeneratorFound";
3670
3671       AliGenPythiaEventHeader *mcHeaderPythia = dynamic_cast<AliGenPythiaEventHeader*>(static_cast<AliMCEvent*>(MCEvent)->GenEventHeader());
3672       if (mcHeaderPythia) fGeneratorNames[0] = "NoCocktailGeneratorFound_Pythia";
3673       AliGenDPMjetEventHeader *mcHeaderPhojet = dynamic_cast<AliGenDPMjetEventHeader*>(static_cast<AliMCEvent*>(MCEvent)->GenEventHeader());
3674       if (mcHeaderPhojet) fGeneratorNames[0] = "NoCocktailGeneratorFound_Phojet";
3675       AliGenHijingEventHeader *mcHeaderHijing = dynamic_cast<AliGenHijingEventHeader*>(static_cast<AliMCEvent*>(MCEvent)->GenEventHeader());
3676       if (mcHeaderHijing) fGeneratorNames[0] = "NoCocktailGeneratorFound_Hijing";
3677
3678       SetRejectExtraSignalsCut(0);
3679    }
3680
3681 }
3682 //_________________________________________________________________________
3683 Int_t AliConversionCuts::IsParticleFromBGEvent(Int_t index, AliStack *MCStack, AliVEvent *InputEvent){
3684
3685    // Not Accepted == kFALSE == 0
3686    //     Accepted ==  kTRUE == 1
3687    //  FirstHeader ==  kTRUE == 3
3688    if(index < 0) return 0; // No Particle
3689
3690    Int_t accepted = 0;
3691    if(!InputEvent || InputEvent->IsA()==AliESDEvent::Class()){
3692       if( index >= MCStack->GetNprimary()){ // Secondary Particle
3693          if( ((TParticle*)MCStack->Particle(index))->GetMother(0) < 0) return 1; // Secondary Particle without Mother??
3694          return IsParticleFromBGEvent(((TParticle*)MCStack->Particle(index))->GetMother(0),MCStack,InputEvent);
3695       }
3696       for(Int_t i = 0;i<fnHeaders;i++){
3697          if(index >= fNotRejectedStart[i] && index <= fNotRejectedEnd[i]){
3698             accepted = 1;
3699             if(i == 0) accepted = 2; // MB Header
3700          }
3701       }
3702    }
3703    else if(InputEvent->IsA()==AliAODEvent::Class()){
3704       TClonesArray *AODMCTrackArray = dynamic_cast<TClonesArray*>(InputEvent->FindListObject(AliAODMCParticle::StdBranchName()));
3705       AliAODMCParticle *aodMCParticle = static_cast<AliAODMCParticle*>(AODMCTrackArray->At(index));
3706       if(!aodMCParticle) return 1; // Photon Without a Mother ? --> Accepted
3707       if(!aodMCParticle->IsPrimary()){
3708          if( aodMCParticle->GetMother() < 0) return 1;// Secondary Particle without Mother??
3709          return IsParticleFromBGEvent(aodMCParticle->GetMother(),MCStack,InputEvent);
3710       }
3711       index = abs(static_cast<AliAODMCParticle*>(AODMCTrackArray->At(index))->GetLabel());
3712       for(Int_t i = 0;i<fnHeaders;i++){
3713          if(index >= fNotRejectedStart[i] && index <= fNotRejectedEnd[i]){
3714             accepted = 1;
3715             if(i == 0) accepted = 2; // MB Header
3716          }
3717       }
3718    }
3719
3720    return accepted;
3721 }
3722 //_________________________________________________________________________
3723 Int_t AliConversionCuts::IsEventAcceptedByConversionCut(AliConversionCuts *ReaderCuts, AliVEvent *InputEvent, AliMCEvent *MCEvent, Bool_t isHeavyIon){
3724
3725    if ( !IsTriggerSelected(InputEvent) )
3726       return 3;
3727
3728    if(isHeavyIon && !(IsCentralitySelected(InputEvent,MCEvent)))
3729       return 1; // Check Centrality --> Not Accepted => eventQuality = 1
3730       
3731       
3732    if(!isHeavyIon && GetIsFromPileup()){
3733       if(InputEvent->IsPileupFromSPD(3,0.8,3.,2.,5.)){
3734
3735          return 6; // Check Pileup --> Not Accepted => eventQuality = 6
3736       }
3737    }
3738
3739    Bool_t hasV0And = ReaderCuts->HasV0AND();
3740    Bool_t isSDDFired = ReaderCuts->IsSDDFired();
3741    if( (IsSpecialTrigger() == 2 || IsSpecialTrigger() == 3) && !isSDDFired && !MCEvent)
3742       return 7; // With SDD requested but no fired
3743
3744    if( (IsSpecialTrigger() == 1 || IsSpecialTrigger() == 3) && !hasV0And)
3745       return 8; // V0AND requested but no fired
3746
3747    if(hCentrality)hCentrality->Fill(GetCentrality(InputEvent));
3748    if(hCentralityVsNumberOfPrimaryTracks)
3749       hCentralityVsNumberOfPrimaryTracks->Fill(GetCentrality(InputEvent),
3750                                              ((AliV0ReaderV1*)AliAnalysisManager::GetAnalysisManager()
3751                                                 ->GetTask("V0ReaderV1"))->GetNumberOfPrimaryTracks());     
3752
3753    return 0;
3754 }
3755 //_________________________________________________________________________
3756 Float_t AliConversionCuts::GetWeightForMeson(TString period, Int_t index, AliStack *MCStack, AliVEvent *InputEvent){
3757    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.Contains("LHC13d2"))) return 1.;
3758
3759    Int_t kCaseGen = 0;
3760    for (Int_t i = 0; i < fnHeaders; i++){
3761       if (index >= fNotRejectedStart[i] && index < fNotRejectedEnd[i]+1){
3762          if (fGeneratorNames[i].CompareTo("Pythia") == 0){
3763             kCaseGen = 1;
3764          } else if (fGeneratorNames[i].CompareTo("DPMJET") == 0){
3765             kCaseGen = 2;
3766          } else if (fGeneratorNames[i].CompareTo("HIJING") == 0 ||
3767                     fGeneratorNames[i].CompareTo("Hijing") == 0 ||
3768                     fGeneratorNames[i].Contains("hijing")){
3769             kCaseGen = 3;
3770          } else if (fGeneratorNames[i].CompareTo("BOX") == 0){
3771              kCaseGen = 4;
3772          } else if (fGeneratorNames[i].CompareTo("PARAM") == 0){
3773             kCaseGen = 5;
3774          } else if (fGeneratorNames[i].CompareTo("NoCocktailGeneratorFound") == 0){
3775             kCaseGen = 6;
3776          } else if (fGeneratorNames[i].CompareTo("NoCocktailGeneratorFound_Pythia") == 0){
3777             kCaseGen = 1;
3778          } else if (fGeneratorNames[i].CompareTo("NoCocktailGeneratorFound_Phojet") == 0){
3779             kCaseGen = 2;
3780          } else if (fGeneratorNames[i].CompareTo("NoCocktailGeneratorFound_Hijing") == 0){
3781             kCaseGen = 3;
3782          }
3783          if (period.Contains("LHC13d2")){
3784             kCaseGen = 3;
3785          }
3786       }
3787    }
3788    if (kCaseGen == 0) return 1;
3789
3790
3791    Double_t mesonPt = 0;
3792    Double_t mesonMass = 0;
3793    Int_t PDGCode = 0;
3794    if(!InputEvent || InputEvent->IsA()==AliESDEvent::Class()){
3795       mesonPt = ((TParticle*)MCStack->Particle(index))->Pt();
3796       mesonMass = ((TParticle*)MCStack->Particle(index))->GetCalcMass();
3797       PDGCode = ((TParticle*)MCStack->Particle(index))->GetPdgCode();
3798    }
3799    else if(InputEvent->IsA()==AliAODEvent::Class()){
3800       TClonesArray *AODMCTrackArray = dynamic_cast<TClonesArray*>(InputEvent->FindListObject(AliAODMCParticle::StdBranchName()));
3801       AliAODMCParticle *aodMCParticle = static_cast<AliAODMCParticle*>(AODMCTrackArray->At(index));
3802       mesonPt = aodMCParticle->Pt();
3803       mesonMass = aodMCParticle->GetCalcMass();
3804       PDGCode = aodMCParticle->GetPdgCode();
3805    }
3806
3807    Float_t functionResultMC = 1.;
3808    if (kCaseGen == 1){ // Pythia 6
3809       Float_t dNdyMC = 2.1462;
3810       Float_t nMC = 7.06055;
3811       Float_t tMC = 0.12533;
3812       if ( PDGCode ==  111){
3813          dNdyMC = 2.1462;
3814          nMC = 7.06055;
3815          tMC = 0.12533;
3816       } else if ( PDGCode ==  221){
3817          dNdyMC = 0.2357;
3818          nMC = 5.9105;
3819          tMC = 0.1525;
3820       }
3821       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);
3822    } else if (kCaseGen == 2){ // Phojet
3823       Float_t dNdyMC = 2.35978;
3824       Float_t nMC = 6.81795;
3825       Float_t tMC = 0.11492;
3826       if ( PDGCode ==  111){
3827          dNdyMC = 2.35978;
3828          nMC = 6.81795;
3829          tMC = 0.11492;
3830       } else if ( PDGCode ==  221){
3831          dNdyMC = 0.3690;
3832          nMC = 5.55809;
3833          tMC = 0.13387;
3834       }
3835       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);
3836    } else if (kCaseGen == 4){ // BOX generators pp
3837 //       functionResultMC = 1./sqrt(1.-mesonMass*mesonMass/((mesonMass*mesonMass+mesonPt*mesonPt)*cosh(mesonY)*cosh(mesonY)));
3838       Float_t a = 0.23437;
3839       Float_t b = 5.6661;
3840       Float_t c = -1430.5863;
3841       Float_t d = -0.6966624;
3842       Float_t e = 252.3742;
3843       if ( PDGCode ==  111){
3844          a = 0.23437;
3845          b = 5.6661;
3846          c = -1430.5863;
3847          d = -0.6966624;
3848          e = 252.3742;
3849       } else if ( PDGCode ==  221){
3850          a = 0.10399;
3851          b = 4.35311;
3852          c = -12.17723;
3853          d = -0.01172;
3854          e =1.85140;
3855       }
3856       functionResultMC = a*TMath::Power(mesonPt,-1.*(b+c/(TMath::Power(mesonPt,d)+e)))*1./mesonPt *1./1.6 *1./(2.* TMath::Pi());
3857 //       cout << functionResultMC << endl;
3858    } else if (kCaseGen == 3 ){ // HIJING
3859       if ( PDGCode ==  111 && fDoReweightHistoMCPi0 && hReweightMCHistPi0!= 0x0){
3860          functionResultMC = hReweightMCHistPi0->Interpolate(mesonPt);
3861       }
3862       if ( PDGCode ==  221 && fDoReweightHistoMCEta && hReweightMCHistEta!= 0x0){
3863          functionResultMC = hReweightMCHistEta->Interpolate(mesonPt);
3864       }
3865       if ( PDGCode ==  310 && fDoReweightHistoMCK0s && hReweightMCHistK0s!= 0x0){
3866          functionResultMC = hReweightMCHistK0s->Interpolate(mesonPt);
3867       }
3868    }
3869
3870    Float_t functionResultData = 1;
3871    if (kCaseGen == 1 || kCaseGen == 2 || kCaseGen == 4 ){
3872       Float_t dNdyData = 2.2328;
3873       Float_t nData = 7.1473;
3874       Float_t tData = 0.1346;
3875       if ( PDGCode ==  111){
3876          dNdyData = 2.2328;
3877          nData = 7.1473;
3878          tData = 0.1346;
3879       } else if ( PDGCode ==  221){
3880          dNdyData = 0.38992; //be careful this fit is not optimal, eta in data still has problems
3881          nData = 5.72778;
3882          tData = 0.13835;
3883       }
3884       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);
3885 //       cout << functionResultData << endl;
3886    } else {
3887       if ( PDGCode ==  111 && fDoReweightHistoMCPi0 && fFitDataPi0!= 0x0){
3888          functionResultData = fFitDataPi0->Eval(mesonPt);
3889       }
3890       if ( PDGCode ==  221 && fDoReweightHistoMCEta && fFitDataEta!= 0x0){
3891          functionResultData = fFitDataEta->Eval(mesonPt);
3892       }
3893       if ( PDGCode ==  310 && fDoReweightHistoMCK0s && fFitDataK0s!= 0x0){
3894          functionResultData = fFitDataK0s->Eval(mesonPt);
3895       }
3896
3897    }
3898
3899    Double_t weight = 1;
3900    if (PDGCode ==  111 || PDGCode ==  221){
3901       if (functionResultData != 0. && functionResultMC != 0. && isfinite(functionResultData) && isfinite(functionResultMC)){
3902          weight = functionResultData/functionResultMC;
3903          if ( kCaseGen == 3){
3904             if (!(fDoReweightHistoMCPi0 && hReweightMCHistPi0!= 0x0 && PDGCode ==  111)){
3905                weight = 1.;
3906             }
3907          }
3908          if (!isfinite(functionResultData)) weight = 1.;
3909          if (!isfinite(weight)) weight = 1.;
3910       }
3911    } else if (PDGCode ==  310 && functionResultMC != 0 && isfinite(functionResultMC)){
3912         weight = functionResultMC;
3913    }
3914
3915 //    if (fModCentralityClass == 0 && fCentralityMin == 4 && fCentralityMax == 6 && PDGCode ==  111){
3916 //        cout << period.Data() << "\t" << kCaseGen << "\t" <<fModCentralityClass<< "\t" <<fCentralityMin<< "\t" <<fCentralityMax << "\t" << mesonPt << "\t" <<mesonMass<< "\t"<<functionResultData << "\t"<< functionResultMC << "\t" << weight <<endl;
3917 //    }
3918    return weight;
3919 }
3920 ///________________________________________________________________________
3921 AliConversionCuts* AliConversionCuts::GetStandardCuts2010PbPb(){
3922     //Create and return standard 2010 PbPb cuts
3923     AliConversionCuts *cuts=new AliConversionCuts("StandardCuts2010PbPb","StandardCuts2010PbPb");
3924     if(!cuts->InitializeCutsFromCutString("100000204209297002322000000")){
3925    cout<<"Warning: Initialization of Standardcuts2010PbPb failed"<<endl;}
3926     return cuts;
3927 }
3928
3929 ///________________________________________________________________________
3930 AliConversionCuts* AliConversionCuts::GetStandardCuts2010pp(){
3931     //Create and return standard 2010 PbPb cuts
3932     AliConversionCuts *cuts=new AliConversionCuts("StandardCuts2010pp","StandardCuts2010pp");
3933     if(!cuts->InitializeCutsFromCutString("000001100209366300380000000")){
3934    cout<<"Warning: Initialization of Standardcuts2010pp failed"<<endl;}
3935     return cuts;
3936 }
3937 ///________________________________________________________________________
3938 void AliConversionCuts::GetCorrectEtaShiftFromPeriod(TString periodName){
3939
3940    if(periodName.CompareTo("LHC12g") == 0 || //pilot run 2012
3941       periodName.CompareTo("LHC13b") == 0 || //mainly minimum bias
3942       periodName.CompareTo("LHC13c") == 0 || //mainly minimum bias
3943       periodName.CompareTo("LHC13d") == 0 || //mainly triggered
3944       periodName.CompareTo("LHC13e") == 0 || //mainly triggered
3945       periodName.CompareTo("LHC13c3") == 0 || //MC Starlight, anchor LHC13d+e
3946       periodName.CompareTo("LHC13c2") == 0 || //MC Starlight, coherent J/Psi, UPC muon anchor LHC13d+e
3947       periodName.CompareTo("LHC13b4") == 0 || //MC Pythia 6 (Jet-Jet), anchor LHC13b
3948       periodName.CompareTo("LHC13b2_fix_1") == 0 || //MC DPMJET, anchr LHC13b+c
3949       periodName.CompareTo("LHC13b2_efix_p1") == 0 || //MC DPMJET, anchr LHC13b+c
3950       periodName.CompareTo("LHC13b2_efix_p2") == 0 || //MC DPMJET, anchr LHC13b+c
3951       periodName.CompareTo("LHC13b2_efix_p3") == 0 || //MC DPMJET, anchr LHC13b+c
3952       periodName.CompareTo("LHC13b2_efix_p4") == 0 || //MC DPMJET, anchr LHC13b+c
3953       periodName.CompareTo("LHC13b3") == 0 || //MC HIJING, weighted to number of events per run, anchor LHC13b
3954       periodName.CompareTo("LHC13b2") == 0 ||  // MC DPMJET, wrong energy, anchor LHC13b
3955       periodName.CompareTo("LHC13b2_plus") == 0 || // MC DPMJET, weighted to number event per run, anchor LHC13b
3956       periodName.CompareTo("LHC13c1_bis") == 0 || // MC AMPT fast generation, pT hardbin, anchor ?
3957       periodName.CompareTo("LHC13c1") == 0 || // MC AMPT fast generation, anchor ?
3958       periodName.CompareTo("LHC13b1") == 0 || // MC DPMJET, fragments, with fixed label 0, anchor LHC12g
3959       periodName.CompareTo("LHC12g4b_fix") == 0 || // MC DPMJET, with fixed label 0, anchor LHC12g
3960       periodName.CompareTo("LHC12g1_fix") == 0 || // MC ?, with fixed label 0, anchor LHC12g
3961       periodName.CompareTo("LHC12g4c") == 0 || // MC DPMJET, shifted vertex runs, anchor LHC12g
3962       periodName.CompareTo("LHC12h6") == 0 || // MC muon cocktail, anchor LHC12g
3963       periodName.CompareTo("LHC12g4b") == 0 || // MC DPMJET 3rd iteration, anchor LHC12g
3964       periodName.CompareTo("LHC12g4a") == 0 || // MC DPMJET improved, anchor LHC12g
3965       periodName.CompareTo("LHC12g4") == 0 || // MC DPMJET, anchor LHC12g
3966       periodName.CompareTo("LHC12g5") == 0 || // MC PHOJET, anchor LHC12g
3967       periodName.CompareTo("LHC12g2") == 0 || // MC Starlight background, anchor LHC12g
3968       periodName.CompareTo("LHC12g1") == 0 ) // MC ?, anchor LHC12g
3969       {
3970          printf(" Gamma Conversion Cuts %s :: pPb Run doing Eta Shift of %f \n\n",(GetCutNumber()).Data(),-0.465);
3971          SetEtaShift(-0.465);
3972       }
3973    else if(periodName.CompareTo("LHC13f") == 0 ||
3974            periodName.CompareTo("LHC13c6b") == 0 ||// MC Jpsi -> mumu, anchor LHC13f
3975            periodName.CompareTo("LHC13c5") == 0 || //MC Starlight, gamma gamma UPC muon, anchor LHC13f
3976            periodName.CompareTo("LHC13c4") == 0 )//MC Starlight, coherent JPsi, UPC muon, anchor LHC13f
3977       {
3978          printf(" Gamma Conversion Cuts %s :: Pbp Run doing Eta Shift of %f \n\n",(GetCutNumber()).Data(),0.465);
3979          SetEtaShift(+0.465);
3980       }
3981    else printf(" Gamma Conversion Cuts %s :: Automatic Eta Shift requested but Period is not known -> No Shift \n\n",(GetCutNumber()).Data());
3982 }
3983
3984 ///________________________________________________________________________
3985 Bool_t AliConversionCuts::InPlaneOutOfPlaneCut(Double_t photonPhi, Double_t eventPlaneAngle, Bool_t fill){
3986    
3987    //GetPhotonPhi() 0-2 Pi  //eventPlaneAngle -1pi-1pi
3988    eventPlaneAngle=eventPlaneAngle+TMath::Pi();
3989    Double_t gammaToEPAngle = eventPlaneAngle-photonPhi;
3990    if(gammaToEPAngle < 0) gammaToEPAngle=gammaToEPAngle+2*TMath::Pi();
3991    gammaToEPAngle = gammaToEPAngle-TMath::Pi(); // angle from -pi +pi
3992
3993    if(!fInPlaneOutOfPlane){
3994       if(fill&&hEventPlanePhi)hEventPlanePhi->Fill(gammaToEPAngle);
3995       return kTRUE;
3996    }
3997    else if(fInPlaneOutOfPlane == 1){
3998       if(abs(gammaToEPAngle)<=0.25*TMath::Pi() || abs(gammaToEPAngle)>=0.75*TMath::Pi()){
3999          if(fill&&hEventPlanePhi)hEventPlanePhi->Fill(gammaToEPAngle);
4000          return kTRUE;
4001       }
4002       else return kFALSE;
4003    }
4004    else if(fInPlaneOutOfPlane == 2){
4005       if(abs(gammaToEPAngle)>0.25*TMath::Pi() && abs(gammaToEPAngle)<0.75*TMath::Pi()){
4006          if(fill&&hEventPlanePhi)hEventPlanePhi->Fill(gammaToEPAngle);
4007          return kTRUE;
4008       }
4009       else return kFALSE;
4010    }
4011    return kFALSE;
4012
4013 }
4014