]> git.uio.no Git - u/mrichter/AliRoot.git/blob - PWGGA/GammaConv/AliConversionCuts.cxx
fix for reader
[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 "AliStack.h"
36 #include "AliAODConversionMother.h"
37 #include "TObjString.h"
38 #include "AliAODEvent.h"
39 #include "AliESDEvent.h"
40 #include "AliCentrality.h"
41 #include "TList.h"
42 #include "AliLog.h"
43 #include "AliGenCocktailEventHeader.h"
44
45 class iostream;
46
47 using namespace std;
48
49 ClassImp(AliConversionCuts)
50
51
52 const char* AliConversionCuts::fgkCutNames[AliConversionCuts::kNCuts] = {
53    "HeavyIon",
54    "CentralityMin",
55    "CentralityMax",
56    "SelectV0AND",
57    "MultiplicityMethod",
58    "RemovePileUp",
59    "RejectExtraSignals",
60    "V0FinderType",
61    "EtaCut",
62    "MinRCut",
63    "SinglePtCut",
64    "ClsTPCCut",
65    "ededxSigmaCut",
66    "pidedxSigmaCut",
67    "piMomdedxSigmaCut",
68    "piMaxMomdedxSigmaCut",
69    "LowPRejectionSigmaCut",
70    "TOFelectronPID",
71    "QtMaxCut",
72    "Chi2GammaCut",
73    "PsiPair",
74    "DoPhotonAsymmetryCut",
75    "CosinePointingAngle",
76    "SharedElectronCuts",
77    "RejectToCloseV0s",
78 };
79
80
81 //________________________________________________________________________
82 AliConversionCuts::AliConversionCuts(const char *name,const char *title) :
83    AliAnalysisCuts(name,title),
84    fHistograms(NULL),
85    fHeaderList(NULL),                                                     
86    fPIDResponse(NULL),
87    fEventQuality(-1),
88    fMaxR(200),
89    fMinR(0),
90    fEtaCut(0.9),
91    fEtaCutMin(-0.1),
92    fPtCut(0.02),
93    fSinglePtCut(0),
94    fMaxZ(1000),
95    fMinClsTPC(0.),
96    fMinClsTPCToF(0.),
97    fLineCutZRSlope(0.),
98    fLineCutZValue(0),
99    fLineCutZRSlopeMin(0.),
100    fLineCutZValueMin(0),
101    fChi2CutConversion(1000),
102    fPIDProbabilityCutNegativeParticle(0),
103    fPIDProbabilityCutPositiveParticle(0),
104    fDodEdxSigmaCut(kTRUE),
105    fDoTOFsigmaCut(kFALSE), 
106    fPIDTRDEfficiency(1),
107    fDoTRDPID(kFALSE),
108    fPIDnSigmaAboveElectronLine(100),
109    fPIDnSigmaBelowElectronLine(-100),
110    fTofPIDnSigmaAboveElectronLine(100), 
111    fTofPIDnSigmaBelowElectronLine(-100), 
112    fPIDnSigmaAbovePionLine(0),
113    fPIDnSigmaAbovePionLineHighPt(-100),
114    fPIDMinPnSigmaAbovePionLine(0),
115    fPIDMaxPnSigmaAbovePionLine(0),
116    fDoKaonRejectionLowP(kFALSE),
117    fDoProtonRejectionLowP(kFALSE),
118    fDoPionRejectionLowP(kFALSE),
119    fPIDnSigmaAtLowPAroundKaonLine(0),
120    fPIDnSigmaAtLowPAroundProtonLine(0),
121    fPIDnSigmaAtLowPAroundPionLine(0),
122    fPIDMinPKaonRejectionLowP(0),
123    fPIDMinPProtonRejectionLowP(0),
124    fPIDMinPPionRejectionLowP(0),
125    fDoQtGammaSelection(kTRUE),
126    fDoHighPtQtGammaSelection(kFALSE), // RRnew
127    fQtMax(100),
128    fHighPtQtMax(0.), // RRnew
129    fPtBorderForQt(0), // RRnew
130    fXVertexCut(0.),
131    fYVertexCut(0.),
132    fZVertexCut(0.),
133    fNSigmaMass(0.),
134    fUseEtaMinCut(kFALSE),
135    fUseOnFlyV0Finder(kTRUE),
136    fDoPhotonAsymmetryCut(kTRUE),
137    fMinPPhotonAsymmetryCut(100.),
138    fMinPhotonAsymmetry(0.),
139    fIsHeavyIon(kFALSE),
140    fDetectorCentrality(0),
141    fModCentralityClass(0),
142    fMaxVertexZ(10),
143    fCentralityMin(0),
144    fCentralityMax(0),
145    fUseCorrectedTPCClsInfo(kFALSE),
146    fUseTOFpid(kFALSE),
147    fMultiplicityMethod(0),
148    fSelectV0AND(kFALSE),
149    fRemovePileUp(kFALSE),
150    fOpeningAngle(0.005),
151    fPsiPairCut(10000),
152    fCosPAngleCut(10000),
153    fDoToCloseV0sCut(kFALSE),
154    fRejectExtraSignals(0),
155    fminV0Dist(200.),
156    fDoSharedElecCut(kFALSE),
157    fOfflineTriggerMask(0),
158    fRandom(0),
159    fSizeElectronArray(500),
160    fElectronLabelArray(NULL),
161    fConversionPointXArray(0),
162    fConversionPointYArray(0),
163    fConversionPointZArray(0),
164    fnHeaders(100),
165    fNotRejectedStart(NULL),
166    fNotRejectedEnd(NULL),
167    fGeneratorNames(NULL),
168    fCutString(NULL),
169    hdEdxCuts(NULL),
170    hTPCdEdxbefore(NULL),
171    hTPCdEdxafter(NULL),
172    hTOFbefore(NULL),
173    hTOFafter(NULL),
174    hTrackCuts(NULL),
175    hPhotonCuts(NULL),
176    hInvMassbefore(NULL),
177    hArmenterosbefore(NULL),
178    hInvMassafter(NULL),
179    hArmenterosafter(NULL),
180    hAcceptanceCuts(NULL),
181    hCutIndex(NULL),
182    hV0EventCuts(NULL),
183    hCentrality(NULL),
184    hVertexZ(NULL),
185    hTriggerClass(NULL)
186 {
187    InitPIDResponse();
188    for(Int_t jj=0;jj<kNCuts;jj++){fCuts[jj]=0;}
189    fCutString=new TObjString((GetCutNumber()).Data());
190
191    fElectronLabelArray = new Int_t[500];
192
193 }
194
195 //________________________________________________________________________
196 AliConversionCuts::~AliConversionCuts() {
197    // Destructor
198    //Deleting fHistograms leads to seg fault it it's added to output collection of a task
199    // if(fHistograms)
200    //   delete fHistograms;
201    // fHistograms = NULL;
202    if(fCutString != NULL){
203       delete fCutString;
204       fCutString = NULL;
205    }
206    if(fElectronLabelArray){
207       delete fElectronLabelArray;
208       fElectronLabelArray = NULL;
209    }
210    if(fNotRejectedStart){
211       delete[] fNotRejectedStart;
212       fNotRejectedStart = NULL;
213    }
214    if(fNotRejectedEnd){
215       delete[] fNotRejectedEnd;
216       fNotRejectedEnd = NULL;
217    }
218    if(fGeneratorNames){
219       delete[] fGeneratorNames;
220       fGeneratorNames = NULL;
221    }
222 }
223
224 //________________________________________________________________________
225 void AliConversionCuts::InitCutHistograms(TString name, Bool_t preCut){
226
227    // Initialize Cut Histograms for QA (only initialized and filled if function is called)
228
229    if(fHistograms != NULL){
230       delete fHistograms;
231       fHistograms=NULL;
232    }
233    if(fHistograms==NULL){
234       fHistograms=new TList();
235       if(name=="")fHistograms->SetName(Form("ConvCuts_%s",GetCutNumber().Data()));
236       else fHistograms->SetName(Form("%s_%s",name.Data(),GetCutNumber().Data()));
237    }
238
239    // IsPhotonSelected
240    hCutIndex=new TH1F(Form("IsPhotonSelected %s",GetCutNumber().Data()),"IsPhotonSelected",10,-0.5,9.5);
241    hCutIndex->GetXaxis()->SetBinLabel(kPhotonIn+1,"in");
242    hCutIndex->GetXaxis()->SetBinLabel(kOnFly+1,"onfly");
243    hCutIndex->GetXaxis()->SetBinLabel(kNoTracks+1,"no tracks");
244    hCutIndex->GetXaxis()->SetBinLabel(kdEdxCuts+1,"dEdx");
245    hCutIndex->GetXaxis()->SetBinLabel(kTrackCuts+1,"Track cuts");
246    hCutIndex->GetXaxis()->SetBinLabel(kConvPointFail+1,"ConvPoint fail");
247    hCutIndex->GetXaxis()->SetBinLabel(kPhotonCuts+1,"PhotonCuts");
248    hCutIndex->GetXaxis()->SetBinLabel(kPhotonOut+1,"out");
249    fHistograms->Add(hCutIndex);
250
251    // Track Cuts
252    hTrackCuts=new TH1F(Form("TrackCuts %s",GetCutNumber().Data()),"TrackCuts",10,-0.5,9.5);
253    hTrackCuts->GetXaxis()->SetBinLabel(1,"in");
254    hTrackCuts->GetXaxis()->SetBinLabel(2,"likesign");
255    hTrackCuts->GetXaxis()->SetBinLabel(3,"ntpccl");
256    hTrackCuts->GetXaxis()->SetBinLabel(4,"acceptance");
257    hTrackCuts->GetXaxis()->SetBinLabel(5,"singlept");
258    hTrackCuts->GetXaxis()->SetBinLabel(6,"TPCrefit");
259    hTrackCuts->GetXaxis()->SetBinLabel(7,"kink");
260    hTrackCuts->GetXaxis()->SetBinLabel(8,"out");
261    fHistograms->Add(hTrackCuts);
262
263    // Photon Cuts
264    hPhotonCuts=new TH1F(Form("PhotonCuts %s",GetCutNumber().Data()),"PhotonCuts",12,-0.5,11.5);
265    hPhotonCuts->GetXaxis()->SetBinLabel(1,"in");
266    hPhotonCuts->GetXaxis()->SetBinLabel(2,"qtcut");
267    hPhotonCuts->GetXaxis()->SetBinLabel(3,"chi2");
268    hPhotonCuts->GetXaxis()->SetBinLabel(4,"acceptance");
269    hPhotonCuts->GetXaxis()->SetBinLabel(5,"asymmetry");
270    hPhotonCuts->GetXaxis()->SetBinLabel(6,"pidprob");
271    hPhotonCuts->GetXaxis()->SetBinLabel(7,"cortpcclinfo");
272    hPhotonCuts->GetXaxis()->SetBinLabel(8,"PsiPair");
273    hPhotonCuts->GetXaxis()->SetBinLabel(9,"CosPAngle");
274    hPhotonCuts->GetXaxis()->SetBinLabel(10,"out");
275    fHistograms->Add(hPhotonCuts);
276
277    if(preCut){
278       hInvMassbefore=new TH1F(Form("InvMass_before %s",GetCutNumber().Data()),"InvMass_before",100,0,0.3);
279       fHistograms->Add(hInvMassbefore);
280       hArmenterosbefore=new TH2F(Form("Armenteros_before %s",GetCutNumber().Data()),"Armenteros_before",200,-1,1,250,0,0.25);
281       fHistograms->Add(hArmenterosbefore);
282    }
283    hInvMassafter=new TH1F(Form("InvMass_after %s",GetCutNumber().Data()),"InvMass_after",100,0,0.3);
284    fHistograms->Add(hInvMassafter);
285    hArmenterosafter=new TH2F(Form("Armenteros_after %s",GetCutNumber().Data()),"Armenteros_after",200,-1,1,250,0,0.25);
286    fHistograms->Add(hArmenterosafter);
287
288    hAcceptanceCuts=new TH1F(Form("PhotonAcceptanceCuts %s",GetCutNumber().Data()),"PhotonAcceptanceCuts",10,-0.5,9.5);
289    hAcceptanceCuts->GetXaxis()->SetBinLabel(1,"in");
290    hAcceptanceCuts->GetXaxis()->SetBinLabel(2,"maxR");
291    hAcceptanceCuts->GetXaxis()->SetBinLabel(3,"minR");
292    hAcceptanceCuts->GetXaxis()->SetBinLabel(4,"line");
293    hAcceptanceCuts->GetXaxis()->SetBinLabel(5,"maxZ");
294    hAcceptanceCuts->GetXaxis()->SetBinLabel(6,"eta");
295    hAcceptanceCuts->GetXaxis()->SetBinLabel(7,"minpt");
296    hAcceptanceCuts->GetXaxis()->SetBinLabel(8,"out");
297    fHistograms->Add(hAcceptanceCuts);
298
299    // dEdx Cuts
300    hdEdxCuts=new TH1F(Form("dEdxCuts %s",GetCutNumber().Data()),"dEdxCuts",10,-0.5,9.5);
301    hdEdxCuts->GetXaxis()->SetBinLabel(1,"in");
302    hdEdxCuts->GetXaxis()->SetBinLabel(2,"TPCelectron");
303    hdEdxCuts->GetXaxis()->SetBinLabel(3,"TPCpion");
304    hdEdxCuts->GetXaxis()->SetBinLabel(4,"TPCpionhighp");
305    hdEdxCuts->GetXaxis()->SetBinLabel(5,"TPCkaonlowprej");
306    hdEdxCuts->GetXaxis()->SetBinLabel(6,"TPCprotonlowprej");
307    hdEdxCuts->GetXaxis()->SetBinLabel(7,"TPCpionlowprej");
308    hdEdxCuts->GetXaxis()->SetBinLabel(8,"TOFelectron");
309    hdEdxCuts->GetXaxis()->SetBinLabel(9,"TRDelectron");
310    hdEdxCuts->GetXaxis()->SetBinLabel(10,"out");
311    fHistograms->Add(hdEdxCuts);
312     
313    TAxis *AxisBeforedEdx = NULL;
314    TAxis *AxisBeforeTOF = NULL;
315    if(preCut){
316       hTPCdEdxbefore=new TH2F(Form("Gamma_dEdx_before %s",GetCutNumber().Data()),"dEdx Gamma before" ,150,0.05,20,400,-10,10);
317       fHistograms->Add(hTPCdEdxbefore);
318       AxisBeforedEdx = hTPCdEdxbefore->GetXaxis(); 
319
320       hTOFbefore=new TH2F(Form("Gamma_TOF_before %s",GetCutNumber().Data()),"TOF Gamma before" ,150,0.05,20,400,-6,10);
321       fHistograms->Add(hTOFbefore);
322       AxisBeforeTOF = hTOFbefore->GetXaxis(); 
323    }
324    hTPCdEdxafter=new TH2F(Form("Gamma_dEdx_after %s",GetCutNumber().Data()),"dEdx Gamma after" ,150,0.05,20,400, -10,10);
325    fHistograms->Add(hTPCdEdxafter);
326
327    hTOFafter=new TH2F(Form("Gamma_TOF_after %s",GetCutNumber().Data()),"TOF Gamma after" ,150,0.05,20,400,-6,10);
328    fHistograms->Add(hTOFafter);
329
330    TAxis *AxisAfter = hTPCdEdxafter->GetXaxis(); 
331    Int_t bins = AxisAfter->GetNbins();
332    Double_t from = AxisAfter->GetXmin();
333    Double_t to = AxisAfter->GetXmax();
334    Double_t *newBins = new Double_t[bins+1];
335    newBins[0] = from;
336    Double_t factor = TMath::Power(to/from, 1./bins);
337    for(Int_t i=1; i<=bins; ++i) newBins[i] = factor * newBins[i-1];
338    AxisAfter->Set(bins, newBins);
339    AxisAfter = hTOFafter->GetXaxis(); 
340    AxisAfter->Set(bins, newBins);
341    if(preCut){
342       AxisBeforedEdx->Set(bins, newBins);
343       AxisBeforeTOF->Set(bins, newBins);
344    }
345    delete [] newBins;
346         
347    // Event Cuts and Info
348    if(preCut){
349       hV0EventCuts=new TH1F(Form("ESD_EventCuts %s",GetCutNumber().Data()),"Event Cuts",10,-0.5,9.5);
350       hV0EventCuts->GetXaxis()->SetBinLabel(1,"in");
351       hV0EventCuts->GetXaxis()->SetBinLabel(2,"OfflineTrigger");
352       hV0EventCuts->GetXaxis()->SetBinLabel(3,"VertexZ");
353       hV0EventCuts->GetXaxis()->SetBinLabel(4,"nvtxcontr");
354       hV0EventCuts->GetXaxis()->SetBinLabel(5,"pileup");
355       hV0EventCuts->GetXaxis()->SetBinLabel(6,"centrsel");
356       hV0EventCuts->GetXaxis()->SetBinLabel(7,"out");
357       fHistograms->Add(hV0EventCuts);
358        
359       hCentrality=new TH1F(Form("Centrality %s",GetCutNumber().Data()),"Centrality",100,0,100);
360       fHistograms->Add(hCentrality);
361       hVertexZ=new TH1F(Form("VertexZ %s",GetCutNumber().Data()),"VertexZ",1000,-50,50);
362       fHistograms->Add(hVertexZ);
363        
364       hTriggerClass= new TH1F(Form("OfflineTrigger %s",GetCutNumber().Data()),"OfflineTrigger",4,-0.5,5.5);
365       hTriggerClass->GetXaxis()->SetBinLabel(1,"kAny");
366       hTriggerClass->GetXaxis()->SetBinLabel(2,"kMB");
367       hTriggerClass->GetXaxis()->SetBinLabel(3,"kCentral");
368       hTriggerClass->GetXaxis()->SetBinLabel(4,"kSemiCentral");
369       hTriggerClass->GetXaxis()->SetBinLabel(4,"kCINT5");
370       hTriggerClass->GetXaxis()->SetBinLabel(5,"kINT7");
371       fHistograms->Add(hTriggerClass);
372    }
373 }
374
375 //________________________________________________________________________
376 Bool_t AliConversionCuts::InitPIDResponse(){
377    // Set Pointer to AliPIDResponse
378
379    AliAnalysisManager *man=AliAnalysisManager::GetAnalysisManager();
380    if(man) { 
381       AliInputEventHandler* inputHandler = (AliInputEventHandler*) (man->GetInputEventHandler());
382       fPIDResponse = (AliPIDResponse*)inputHandler->GetPIDResponse();
383       if(fPIDResponse)return kTRUE;
384     
385    }
386   
387
388    return kFALSE;
389 }
390 ///________________________________________________________________________
391 Bool_t AliConversionCuts::EventIsSelected(AliVEvent *fInputEvent, AliVEvent *fMCEvent){
392    // Process Event Selection
393
394    Int_t cutindex=0;
395    if(hV0EventCuts)hV0EventCuts->Fill(cutindex);
396    cutindex++;
397     
398    // Check for MC event
399    if(fMCEvent){
400       // Check if MC event is correctly loaded
401       AliMCEventHandler* mcHandler = dynamic_cast<AliMCEventHandler*> (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
402       if (!mcHandler){
403          fEventQuality = 2;
404          return kFALSE;
405       }
406       if (!mcHandler->InitOk() ){
407          fEventQuality = 2;
408          return kFALSE;
409       }
410       if (!mcHandler->TreeK() ){
411          fEventQuality = 2;
412          return kFALSE;
413       }
414       if (!mcHandler->TreeTR() ) {
415          fEventQuality = 2;
416          return kFALSE;
417       }
418    }
419     
420    // Event Trigger
421
422    if(!IsTriggerSelected()){
423       if(hV0EventCuts)hV0EventCuts->Fill(cutindex);
424       fEventQuality = 3;
425       return kFALSE;
426    }
427    cutindex++;
428
429    // Z Vertex Position Cut
430    if(!VertexZCut(fInputEvent)){
431       if(hV0EventCuts)hV0EventCuts->Fill(cutindex);
432       fEventQuality = 4;
433       return kFALSE;
434    }
435    cutindex++;
436
437    // Number of Contributors Cut
438    if(GetNumberOfContributorsVtx(fInputEvent)<=0) {
439       if(hV0EventCuts)hV0EventCuts->Fill(cutindex);
440       fEventQuality = 5;
441       return kFALSE;
442    }
443    cutindex++;
444
445    // Pile Up Rejection
446
447    if(fRemovePileUp){
448       if(fInputEvent->IsPileupFromSPD(3,0.8,3.,2.,5.)){
449          if(hV0EventCuts)hV0EventCuts->Fill(cutindex);
450          fEventQuality = 6;
451          return kFALSE;
452       }
453    }
454    cutindex++;
455
456    // Centrality Selection
457    if(!IsCentralitySelected(fInputEvent)){
458       if(hV0EventCuts)hV0EventCuts->Fill(cutindex);
459       fEventQuality = 1;
460       return kFALSE;
461    }
462    cutindex++;
463
464    // Fill Event Histograms
465    if(hV0EventCuts)hV0EventCuts->Fill(cutindex);
466    if(hVertexZ)hVertexZ->Fill(fInputEvent->GetPrimaryVertex()->GetZ());
467    if(hCentrality)hCentrality->Fill(GetCentrality(fInputEvent));
468
469    fEventQuality = 0;
470    return kTRUE;
471 }
472
473 ///________________________________________________________________________
474 Bool_t AliConversionCuts::PhotonIsSelectedMC(TParticle *particle,AliStack *fMCStack,Bool_t checkForConvertedGamma){
475    // MonteCarlo Photon Selection
476
477    if(!fMCStack)return kFALSE;
478
479    if (particle->GetPdgCode() == 22){
480               
481       if(particle->R() > fMaxR) return kFALSE;
482       if(TMath::Abs(particle->Eta())> fEtaCut || TMath::Abs(particle->Eta())< fEtaCutMin)       return kFALSE;
483
484       if(particle->GetMother(0) >-1 && fMCStack->Particle(particle->GetMother(0))->GetPdgCode() == 22){
485          return kFALSE; // no photon as mothers!
486       }
487
488       if(particle->GetMother(0) >= fMCStack->GetNprimary()){
489          return kFALSE; // the gamma has a mother, and it is not a primary particle
490       }
491
492       if(!checkForConvertedGamma) return kTRUE; // return in case of accepted gamma
493
494       // looking for conversion gammas (electron + positron from pairbuilding (= 5) )
495       TParticle* ePos = NULL;
496       TParticle* eNeg = NULL;
497
498       if(particle->GetNDaughters() >= 2){
499          for(Int_t daughterIndex=particle->GetFirstDaughter();daughterIndex<=particle->GetLastDaughter();daughterIndex++){
500             TParticle *tmpDaughter = fMCStack->Particle(daughterIndex);
501             if(tmpDaughter->GetUniqueID() == 5){
502                if(tmpDaughter->GetPdgCode() == 11){
503                   eNeg = tmpDaughter;
504                } else if(tmpDaughter->GetPdgCode() == -11){
505                   ePos = tmpDaughter;
506                }
507             }
508          }
509       }
510
511       if(ePos == NULL || eNeg == NULL){ // means we do not have two daughters from pair production
512          return kFALSE;
513       }
514         
515       if(ePos->Pt()<fSinglePtCut || eNeg->Pt()<fSinglePtCut){
516          return kFALSE; // no reconstruction below the Pt cut
517       }
518         
519       if( TMath::Abs(ePos->Eta())> fEtaCut || TMath::Abs(ePos->Eta())< fEtaCutMin || 
520           TMath::Abs(eNeg->Eta())> fEtaCut || TMath::Abs(eNeg->Eta())< fEtaCutMin ) {
521          return kFALSE;
522       }
523         
524       if(ePos->R()>fMaxR){
525          return kFALSE; // cuts on distance from collision point
526       }
527
528       if(TMath::Abs(ePos->Vz()) > fMaxZ){
529          return kFALSE;  // outside material
530       }
531       if(TMath::Abs(eNeg->Vz()) > fMaxZ){
532          return kFALSE;  // outside material
533       }
534
535       if( ePos->R() <= ((TMath::Abs(ePos->Vz()) * fLineCutZRSlope) - fLineCutZValue)){
536          return kFALSE;  // line cut to exclude regions where we do not reconstruct
537       } else if ( fEtaCutMin != -0.1 &&   ePos->R() >= ((TMath::Abs(ePos->Vz()) * fLineCutZRSlopeMin) - fLineCutZValueMin)){
538          return kFALSE;
539       }
540         
541       if( eNeg->R() <= ((TMath::Abs(eNeg->Vz()) * fLineCutZRSlope) - fLineCutZValue)){
542          return kFALSE; // line cut to exclude regions where we do not reconstruct
543       } else if ( fEtaCutMin != -0.1 &&   eNeg->R() >= ((TMath::Abs(eNeg->Vz()) * fLineCutZRSlopeMin) - fLineCutZValueMin)){
544          return kFALSE;
545       }
546         
547       return kTRUE;
548       //if(AcceptanceCut(particle,ePos,eNeg))return kTRUE;
549    }
550    return kFALSE;
551 }
552
553
554 ///________________________________________________________________________
555 Bool_t AliConversionCuts::PhotonCuts(AliConversionPhotonBase *photon,AliVEvent *event)
556 {   // Specific Photon Cuts
557
558    Int_t cutIndex = 0;
559    if(hPhotonCuts)hPhotonCuts->Fill(cutIndex);
560    cutIndex++;
561
562    // Fill Histos before Cuts
563    if(hInvMassbefore)hInvMassbefore->Fill(photon->GetMass());
564    if(hArmenterosbefore)hArmenterosbefore->Fill(photon->GetArmenterosAlpha(),photon->GetArmenterosQt());
565
566    // Gamma selection based on QT from Armenteros
567    if(fDoQtGammaSelection == kTRUE){
568       if(!ArmenterosQtCut(photon)){
569          if(hPhotonCuts)hPhotonCuts->Fill(cutIndex); //1
570          return kFALSE;
571       }
572    }
573    cutIndex++; //2
574
575    // Chi Cut
576    if(photon->GetChi2perNDF() > fChi2CutConversion || photon->GetChi2perNDF() <=0){
577       {
578          if(hPhotonCuts)hPhotonCuts->Fill(cutIndex); //2
579          return kFALSE;
580       }
581    }
582    cutIndex++;//3
583
584    // Reconstruction Acceptance Cuts
585    if(!AcceptanceCuts(photon)){
586       if(hPhotonCuts)hPhotonCuts->Fill(cutIndex); //3
587       return kFALSE;
588    }
589
590    cutIndex++; //4
591    // Asymmetry Cut
592    if(fDoPhotonAsymmetryCut == kTRUE){
593       if(!AsymmetryCut(photon,event)){
594          if(hPhotonCuts)hPhotonCuts->Fill(cutIndex); //4
595          return kFALSE;
596       }
597    }
598
599    //Check the pid probability
600    cutIndex++; //5
601    if(!PIDProbabilityCut(photon, event)) {
602       if(hPhotonCuts)hPhotonCuts->Fill(cutIndex); //5
603       return kFALSE;
604    }
605
606    cutIndex++; //6
607    if(!CorrectedTPCClusterCut(photon, event)) {
608       if(hPhotonCuts)hPhotonCuts->Fill(cutIndex); //6
609       return kFALSE;
610    }
611
612
613    cutIndex++; //7
614    if(!PsiPairCut(photon)) {
615       if(hPhotonCuts)hPhotonCuts->Fill(cutIndex); //7
616       return kFALSE;
617    }
618
619    cutIndex++; //8
620    if(!CosinePAngleCut(photon, event)) {
621       if(hPhotonCuts)hPhotonCuts->Fill(cutIndex); //8
622       return kFALSE;
623    }
624
625    cutIndex++; //9
626    if(hPhotonCuts)hPhotonCuts->Fill(cutIndex); //9
627
628    // Histos after Cuts
629    if(hInvMassafter)hInvMassafter->Fill(photon->GetMass());
630    if(hArmenterosafter)hArmenterosafter->Fill(photon->GetArmenterosAlpha(),photon->GetArmenterosQt());
631
632
633    return kTRUE;
634
635 }
636
637 ///________________________________________________________________________
638 Bool_t AliConversionCuts::CorrectedTPCClusterCut(AliConversionPhotonBase *photon, AliVEvent * event)
639 {   //Cut on corrected TPC Cluster Info
640
641    AliVTrack * negTrack = GetTrack(event, photon->GetTrackLabelNegative());
642    AliVTrack * posTrack = GetTrack(event, photon->GetTrackLabelPositive());
643
644    if(!negTrack||!posTrack)return kFALSE;
645
646    Double_t negclsToF=0;
647
648    if (!fUseCorrectedTPCClsInfo ){
649       if(negTrack->GetTPCNclsF()!=0){
650          negclsToF = (Double_t)negTrack->GetNcls(1)/(Double_t)negTrack->GetTPCNclsF();}// Ncluster/Nfindablecluster
651    }
652    else {
653       negclsToF = negTrack->GetTPCClusterInfo(2,0,GetFirstTPCRow(photon->GetConversionRadius()));
654    }
655
656    Double_t posclsToF = 0.;
657    if (!fUseCorrectedTPCClsInfo ){
658       if(posTrack->GetTPCNclsF()!=0     ){
659          posclsToF = (Double_t)posTrack->GetNcls(1)/(Double_t)posTrack->GetTPCNclsF();
660       }
661    }else{
662       posclsToF = posTrack->GetTPCClusterInfo(2,0,GetFirstTPCRow(photon->GetConversionRadius()));
663    }
664
665    if( negclsToF < fMinClsTPCToF ||     posclsToF < fMinClsTPCToF ){
666       return kFALSE;
667    }
668
669    return kTRUE;
670 }
671
672 ///________________________________________________________________________
673 Bool_t AliConversionCuts::PhotonIsSelected(AliConversionPhotonBase *photon, AliVEvent * event)
674 {
675    //Selection of Reconstructed Photons
676
677    FillPhotonCutIndex(kPhotonIn);
678
679    // Get Tracks
680    AliVTrack * negTrack = GetTrack(event, photon->GetTrackLabelNegative());
681    AliVTrack * posTrack = GetTrack(event, photon->GetTrackLabelPositive());
682
683    if(!negTrack || !posTrack) {
684       FillPhotonCutIndex(kNoTracks);
685       return kFALSE;
686    }
687
688    // dEdx Cuts
689    if(!dEdxCuts(negTrack) || !dEdxCuts(posTrack)) {
690       FillPhotonCutIndex(kdEdxCuts);
691       return kFALSE;
692    }
693
694    // Track Cuts
695    if(!TracksAreSelected(negTrack, posTrack)){
696       FillPhotonCutIndex(kTrackCuts);
697       return kFALSE;
698    }
699
700    // Photon Cuts
701    if(!PhotonCuts(photon,event)){
702       FillPhotonCutIndex(kPhotonCuts);
703       return kFALSE;
704    }
705
706    // Photon passed cuts
707    FillPhotonCutIndex(kPhotonOut);
708    return kTRUE;
709 }
710
711 ///________________________________________________________________________
712 Bool_t AliConversionCuts::ArmenterosQtCut(AliConversionPhotonBase *photon)
713 {   // Armenteros Qt Cut
714
715    if(fDoHighPtQtGammaSelection){
716       if(photon->GetPhotonPt() < fPtBorderForQt){
717          if(photon->GetArmenterosQt()>fQtMax){
718             return kFALSE;
719          }
720       } else {
721          if(photon->GetArmenterosQt()>fHighPtQtMax){
722             return kFALSE;
723          }
724       }
725    } else {
726
727       if(photon->GetArmenterosQt()>fQtMax){
728          return kFALSE;
729       }
730    }
731    return kTRUE;
732 }
733
734
735 ///________________________________________________________________________
736 Bool_t AliConversionCuts::AcceptanceCuts(AliConversionPhotonBase *photon) {
737    // Exclude certain areas for photon reconstruction
738
739    Int_t cutIndex=0;
740    if(hAcceptanceCuts)hAcceptanceCuts->Fill(cutIndex);
741    cutIndex++;
742
743    if(photon->GetConversionRadius()>fMaxR){ // cuts on distance from collision point
744       if(hAcceptanceCuts)hAcceptanceCuts->Fill(cutIndex);
745       return kFALSE;
746    }
747    cutIndex++;
748
749    if(photon->GetConversionRadius()<fMinR){ // cuts on distance from collision point
750       if(hAcceptanceCuts)hAcceptanceCuts->Fill(cutIndex);
751       return kFALSE;
752    }
753    cutIndex++;
754
755    if(photon->GetConversionRadius() <= ((TMath::Abs(photon->GetConversionZ())*fLineCutZRSlope)-fLineCutZValue)){
756       if(hAcceptanceCuts)hAcceptanceCuts->Fill(cutIndex);
757       return kFALSE;
758    }
759    else if (fUseEtaMinCut &&  photon->GetConversionRadius() >= ((TMath::Abs(photon->GetConversionZ())*fLineCutZRSlopeMin)-fLineCutZValueMin )){
760       if(hAcceptanceCuts)hAcceptanceCuts->Fill(cutIndex);
761       return kFALSE;
762    }
763    cutIndex++;
764
765    if(TMath::Abs(photon->GetConversionZ()) > fMaxZ ){ // cuts out regions where we do not reconstruct
766       if(hAcceptanceCuts)hAcceptanceCuts->Fill(cutIndex);
767       return kFALSE;
768    }
769    cutIndex++;
770
771
772    if(TMath::Abs(photon->GetPhotonEta())> fEtaCut || TMath::Abs(photon->GetPhotonEta())< fEtaCutMin){
773       if(hAcceptanceCuts)hAcceptanceCuts->Fill(cutIndex);
774       return kFALSE;
775    }
776    cutIndex++;
777
778
779    if(photon->GetPhotonPt()<fPtCut){
780       if(hAcceptanceCuts)hAcceptanceCuts->Fill(cutIndex);
781       return kFALSE;
782    }
783    cutIndex++;
784
785    if(hAcceptanceCuts)hAcceptanceCuts->Fill(cutIndex);
786  
787    return kTRUE;
788 }
789
790
791 ///________________________________________________________________________
792 Bool_t AliConversionCuts::SpecificTrackCuts(AliAODTrack * negTrack, AliAODTrack * posTrack,Int_t &cutIndex) {
793    // Track Cuts which require AOD/ESD specific implementation
794
795    if( !negTrack->IsOn(AliESDtrack::kTPCrefit)  || !posTrack->IsOn(AliESDtrack::kTPCrefit)   )  {
796       if(hTrackCuts)hTrackCuts->Fill(cutIndex);
797       return kFALSE;
798    }
799    cutIndex++;
800
801    AliAODVertex * NegVtxType=negTrack->GetProdVertex();
802    AliAODVertex * PosVtxType=posTrack->GetProdVertex();
803    if((NegVtxType->GetType())==AliAODVertex::kKink  || (PosVtxType->GetType())==AliAODVertex::kKink) {
804       if(hTrackCuts)hTrackCuts->Fill(cutIndex);
805       return kFALSE;
806    }
807    return kTRUE;
808
809 }
810
811
812 ///________________________________________________________________________
813 Bool_t AliConversionCuts::SpecificTrackCuts(AliESDtrack * negTrack, AliESDtrack * posTrack,Int_t &cutIndex) {
814    // Track Cuts which require AOD/ESD specific implementation
815
816    if( !negTrack->IsOn(AliESDtrack::kTPCrefit)  || !posTrack->IsOn(AliESDtrack::kTPCrefit)   )  {
817       if(hTrackCuts)hTrackCuts->Fill(cutIndex);
818       return kFALSE;
819    }
820    cutIndex++;
821
822    if(negTrack->GetKinkIndex(0) > 0  || posTrack->GetKinkIndex(0) > 0 ) {
823       if(hTrackCuts)hTrackCuts->Fill(cutIndex);
824       return kFALSE;
825    }
826    return kTRUE;
827 }
828
829
830
831 ///________________________________________________________________________
832 Bool_t AliConversionCuts::TracksAreSelected(AliVTrack * negTrack, AliVTrack * posTrack) {
833    // Track Selection for Photon Reconstruction
834
835    Int_t cutIndex=0;
836    if(hTrackCuts)hTrackCuts->Fill(cutIndex);
837    cutIndex++;
838
839    // avoid like sign
840    if(negTrack->Charge() == posTrack->Charge()) {
841       if(hTrackCuts)hTrackCuts->Fill(cutIndex);
842       return kFALSE;
843    }
844    cutIndex++;
845
846    // Number of TPC Clusters
847    if( negTrack->GetNcls(1) < fMinClsTPC || posTrack->GetNcls(1) < fMinClsTPC ) {
848       if(hTrackCuts)hTrackCuts->Fill(cutIndex);
849       return kFALSE;
850    }
851    cutIndex++;
852
853    // Acceptance
854
855    if(TMath::Abs(negTrack->Eta()) > fEtaCut || TMath::Abs(negTrack->Eta()) < fEtaCutMin ||
856       TMath::Abs(posTrack->Eta())> fEtaCut || TMath::Abs(posTrack->Eta())< fEtaCutMin) {
857       if(hTrackCuts)hTrackCuts->Fill(cutIndex);
858       return kFALSE;
859    }
860    cutIndex++;
861
862    // Single Pt Cut
863    if( negTrack->Pt()< fSinglePtCut ||  posTrack->Pt()< fSinglePtCut){
864       if(hTrackCuts)hTrackCuts->Fill(cutIndex);
865       return kFALSE;
866    }
867    cutIndex++;
868
869    // AOD ESD specific cuts
870    Bool_t passCuts = kTRUE;
871
872    if(negTrack->IsA()==AliAODTrack::Class()) {
873       passCuts = passCuts * SpecificTrackCuts(static_cast<AliAODTrack*>(negTrack), static_cast<AliAODTrack*>(posTrack),cutIndex);
874    } else { 
875       passCuts = passCuts * SpecificTrackCuts(static_cast<AliESDtrack*>(negTrack), static_cast<AliESDtrack*>(posTrack),cutIndex);
876    }    
877
878    if(!passCuts){
879       if(hTrackCuts)hTrackCuts->Fill(cutIndex);
880       return kFALSE;
881    }
882    cutIndex++;
883
884    if(hTrackCuts)hTrackCuts->Fill(cutIndex);
885
886    return kTRUE;
887                     
888 }
889
890 ///________________________________________________________________________
891 Bool_t AliConversionCuts::dEdxCuts(AliVTrack *fCurrentTrack){
892    // Electron Identification Cuts for Photon reconstruction
893
894    if(!fPIDResponse){InitPIDResponse();}// Try to reinitialize PID Response
895    if(!fPIDResponse){AliError("No PID Response"); return kTRUE;}// if still missing fatal error
896
897    Int_t cutIndex=0;
898    if(hdEdxCuts)hdEdxCuts->Fill(cutIndex);
899    if(hTPCdEdxbefore)hTPCdEdxbefore->Fill(fCurrentTrack->P(),fPIDResponse->NumberOfSigmasTPC(fCurrentTrack, AliPID::kElectron));
900    cutIndex++;
901     
902
903    if(fDodEdxSigmaCut == kTRUE){
904       // TPC Electron Line
905       if( fPIDResponse->NumberOfSigmasTPC(fCurrentTrack,AliPID::kElectron)<fPIDnSigmaBelowElectronLine ||
906           fPIDResponse->NumberOfSigmasTPC(fCurrentTrack,AliPID::kElectron)>fPIDnSigmaAboveElectronLine){
907
908          if(hdEdxCuts)hdEdxCuts->Fill(cutIndex);
909          return kFALSE;
910       }
911       cutIndex++;
912
913       // TPC Pion Line
914       if( fCurrentTrack->P()>fPIDMinPnSigmaAbovePionLine && fCurrentTrack->P()<fPIDMaxPnSigmaAbovePionLine ){
915          if(fPIDResponse->NumberOfSigmasTPC(fCurrentTrack,AliPID::kElectron)>fPIDnSigmaBelowElectronLine &&
916             fPIDResponse->NumberOfSigmasTPC(fCurrentTrack,AliPID::kElectron)<fPIDnSigmaAboveElectronLine&&
917             fPIDResponse->NumberOfSigmasTPC(fCurrentTrack,AliPID::kPion)<fPIDnSigmaAbovePionLine){
918
919             if(hdEdxCuts)hdEdxCuts->Fill(cutIndex);
920             return kFALSE;
921          }
922       }
923       cutIndex++;
924    
925       // High Pt Pion rej
926       if( fCurrentTrack->P()>fPIDMaxPnSigmaAbovePionLine ){
927          if(fPIDResponse->NumberOfSigmasTPC(fCurrentTrack,AliPID::kElectron)>fPIDnSigmaBelowElectronLine &&
928             fPIDResponse->NumberOfSigmasTPC(fCurrentTrack,AliPID::kElectron)<fPIDnSigmaAboveElectronLine&&
929             fPIDResponse->NumberOfSigmasTPC(fCurrentTrack,AliPID::kPion)<fPIDnSigmaAbovePionLineHighPt){
930
931             if(hdEdxCuts)hdEdxCuts->Fill(cutIndex);
932             return kFALSE;
933          }
934       }
935       cutIndex++;
936    }
937    else{cutIndex+=3;}
938
939    if(fDoKaonRejectionLowP == kTRUE){
940       if(fCurrentTrack->P()<fPIDMinPKaonRejectionLowP ){
941          if( TMath::Abs(fPIDResponse->NumberOfSigmasTPC(fCurrentTrack,AliPID::kKaon))<fPIDnSigmaAtLowPAroundKaonLine){
942
943             if(hdEdxCuts)hdEdxCuts->Fill(cutIndex);
944             return kFALSE;
945          }
946       }
947    }
948    cutIndex++;
949    
950    if(fDoProtonRejectionLowP == kTRUE){
951       if( fCurrentTrack->P()<fPIDMinPProtonRejectionLowP ){
952          if( TMath::Abs(fPIDResponse->NumberOfSigmasTPC(fCurrentTrack,AliPID::kProton))<fPIDnSigmaAtLowPAroundProtonLine){
953
954             if(hdEdxCuts)hdEdxCuts->Fill(cutIndex);
955             return kFALSE;
956          }
957       }
958    }
959    cutIndex++;
960    
961    if(fDoPionRejectionLowP == kTRUE){
962       if( fCurrentTrack->P()<fPIDMinPPionRejectionLowP ){
963          if( TMath::Abs(fPIDResponse->NumberOfSigmasTPC(fCurrentTrack,AliPID::kPion))<fPIDnSigmaAtLowPAroundPionLine){
964
965             if(hdEdxCuts)hdEdxCuts->Fill(cutIndex);
966             return kFALSE;
967          }
968       }
969    }
970    cutIndex++;
971    
972    if((fCurrentTrack->GetStatus() & AliESDtrack::kTOFpid) && !(fCurrentTrack->GetStatus() & AliESDtrack::kTOFmismatch)){
973       if(hTOFbefore) hTOFbefore->Fill(fCurrentTrack->P(),fPIDResponse->NumberOfSigmasTOF(fCurrentTrack, AliPID::kElectron));
974       if(fUseTOFpid){
975          if(fPIDResponse->NumberOfSigmasTOF(fCurrentTrack, AliPID::kElectron)>fTofPIDnSigmaAboveElectronLine ||
976             fPIDResponse->NumberOfSigmasTOF(fCurrentTrack, AliPID::kElectron)<fTofPIDnSigmaBelowElectronLine ){
977             if(hdEdxCuts)hdEdxCuts->Fill(cutIndex);
978             return kFALSE;
979          }
980       }
981       if(hTOFafter)hTOFafter->Fill(fCurrentTrack->P(),fPIDResponse->NumberOfSigmasTOF(fCurrentTrack, AliPID::kElectron));
982    }
983    cutIndex++;
984    
985    // Apply TRD PID
986    if(fDoTRDPID){
987       if(!fPIDResponse->IdentifiedAsElectronTRD(fCurrentTrack,fPIDTRDEfficiency)){
988
989          if(hdEdxCuts)hdEdxCuts->Fill(cutIndex);
990          return kFALSE;
991       }
992    }
993    cutIndex++;
994
995    if(hdEdxCuts)hdEdxCuts->Fill(cutIndex);
996    if(hTPCdEdxafter)hTPCdEdxafter->Fill(fCurrentTrack->P(),fPIDResponse->NumberOfSigmasTPC(fCurrentTrack, AliPID::kElectron));
997
998    return kTRUE;
999 }
1000
1001 ///________________________________________________________________________
1002 Bool_t AliConversionCuts::AsymmetryCut(AliConversionPhotonBase * photon,AliVEvent *event) {
1003    // Cut on Energy Assymetry
1004
1005    for(Int_t ii=0;ii<2;ii++){
1006
1007       AliVTrack *track=GetTrack(event,photon->GetTrackLabel(ii));
1008
1009       if( track->P() > fMinPPhotonAsymmetryCut ){
1010          Double_t trackNegAsy=0;
1011          if (photon->GetPhotonP()!=0.){
1012             trackNegAsy= track->P()/photon->GetPhotonP();
1013          }
1014
1015          if( trackNegAsy<fMinPhotonAsymmetry ||trackNegAsy>(1.- fMinPhotonAsymmetry)){
1016             return kFALSE;
1017          }
1018       }
1019
1020    }
1021    return kTRUE;
1022 }
1023
1024 ///________________________________________________________________________
1025 AliVTrack *AliConversionCuts::GetTrack(AliVEvent * event, Int_t label){
1026    //Returns pointer to the track with given ESD label
1027    //(Important for AOD implementation, since Track array in AOD data is different
1028    //from ESD array, but ESD tracklabels are stored in AOD Tracks)
1029
1030    AliESDEvent * esdEvent = dynamic_cast<AliESDEvent*>(event);
1031    if(esdEvent) {
1032       if(label > event->GetNumberOfTracks() ) return NULL;
1033       AliESDtrack * track = esdEvent->GetTrack(label);
1034       return track;
1035         
1036    } else { 
1037       for(Int_t ii=0; ii<event->GetNumberOfTracks(); ii++) {
1038          AliVTrack * track = dynamic_cast<AliVTrack*>(event->GetTrack(ii));
1039           
1040          if(track) { 
1041             if(track->GetID() == label) {
1042                return track;
1043             }
1044          }
1045       }
1046    }
1047   
1048    //AliDebug(5,(Form("track not found %d %d",label,event->GetNumberOfTracks()));
1049    return NULL;
1050 }
1051
1052
1053 ///________________________________________________________________________
1054 Bool_t AliConversionCuts::PIDProbabilityCut(AliConversionPhotonBase *photon, AliVEvent * event){
1055    // Cut on Electron Probability for Photon Reconstruction
1056
1057    AliESDEvent * esdEvent = dynamic_cast<AliESDEvent*>(event);
1058
1059    if(esdEvent){
1060         
1061       Bool_t iResult=kFALSE;
1062         
1063       Double_t *posProbArray = new Double_t[AliPID::kSPECIES];
1064       Double_t *negProbArray = new Double_t[AliPID::kSPECIES];
1065         
1066       AliESDtrack* negTrack     = esdEvent->GetTrack(photon->GetTrackLabelNegative());
1067       AliESDtrack* posTrack     = esdEvent->GetTrack(photon->GetTrackLabelPositive());
1068         
1069       if(negProbArray && posProbArray){
1070           
1071          negTrack->GetTPCpid(negProbArray);
1072          posTrack->GetTPCpid(posProbArray);
1073
1074          if(negProbArray[AliPID::kElectron]>=fPIDProbabilityCutNegativeParticle && posProbArray[AliPID::kElectron]>=fPIDProbabilityCutPositiveParticle){
1075             iResult=kTRUE;
1076          }
1077       }
1078         
1079       delete [] posProbArray;
1080       delete [] negProbArray;
1081       return iResult;
1082
1083    } else {
1084       ///Not possible for AODs
1085       return kTRUE;
1086    }
1087
1088
1089
1090   
1091 }
1092
1093
1094 ///________________________________________________________________________
1095 Bool_t AliConversionCuts::AcceptanceCut(TParticle *particle, TParticle * ePos,TParticle* eNeg){
1096    // MC Acceptance Cuts
1097    //(Certain areas were excluded for photon reconstruction)
1098
1099    if(particle->R()>fMaxR){
1100       return kFALSE;}
1101
1102    if(ePos->R()>fMaxR){
1103       return kFALSE;
1104    }
1105
1106    if(ePos->R()<fMinR){
1107       return kFALSE;
1108    }
1109
1110    if( ePos->R() <= ((TMath::Abs(ePos->Vz())*fLineCutZRSlope)-fLineCutZValue)){
1111       return kFALSE;
1112    }
1113    else if (fUseEtaMinCut &&  ePos->R() >= ((TMath::Abs(ePos->Vz())*fLineCutZRSlopeMin)-fLineCutZValueMin )){
1114       return kFALSE;
1115    }
1116
1117    if(TMath::Abs(eNeg->Vz()) > fMaxZ ){ // cuts out regions where we do not reconstruct
1118       return kFALSE;
1119    }
1120
1121    if(eNeg->Vz()!=ePos->Vz()||eNeg->R()!=ePos->R()){
1122       return kFALSE;
1123    }
1124
1125    if(TMath::Abs(ePos->Vz()) > fMaxZ ){ // cuts out regions where we do not reconstruct
1126       return kFALSE;
1127    }
1128
1129    if(TMath::Abs(particle->Eta())> fEtaCut || TMath::Abs(particle->Eta())< fEtaCutMin){
1130       return kFALSE;
1131    }
1132
1133    if(TMath::Abs(ePos->Eta())> fEtaCut || TMath::Abs(ePos->Eta())< fEtaCutMin){
1134       return kFALSE;
1135    }
1136
1137    if(TMath::Abs(eNeg->Eta())> fEtaCut || TMath::Abs(eNeg->Eta())< fEtaCutMin){
1138       return kFALSE;
1139    }
1140
1141    if( ePos->Pt()< fSinglePtCut ||  eNeg->Pt()< fSinglePtCut){
1142       return kFALSE;
1143    }
1144
1145    if(particle->Pt()<fPtCut){
1146       return kFALSE;
1147    }
1148
1149    return kTRUE;
1150 }
1151 ///________________________________________________________________________
1152 Bool_t AliConversionCuts::UpdateCutString(cutIds cutID, Int_t value) {
1153    ///Update the cut string (if it has been created yet)
1154
1155    if(fCutString && fCutString->GetString().Length() == kNCuts) {
1156       fCutString->SetString(GetCutNumber());
1157    } else {
1158       return kFALSE;
1159    }
1160    return kTRUE;
1161 }
1162
1163 ///________________________________________________________________________
1164 Bool_t AliConversionCuts::InitializeCutsFromCutString(const TString analysisCutSelection ) {
1165    // Initialize Cuts from a given Cut string
1166
1167    AliInfo(Form("Set Photoncut Number: %s",analysisCutSelection.Data()));
1168    if(analysisCutSelection.Length()!=kNCuts) {
1169       AliError(Form("Cut selection has the wrong length! size is %d, number of cuts is %d", analysisCutSelection.Length(), kNCuts));
1170       return kFALSE;
1171    }
1172    if(!analysisCutSelection.IsDigit()){
1173       AliError("Cut selection contains characters");
1174       return kFALSE;
1175    }
1176   
1177    const char *cutSelection = analysisCutSelection.Data();
1178 #define ASSIGNARRAY(i)  fCuts[i] = cutSelection[i] - '0'
1179    for(Int_t ii=0;ii<kNCuts;ii++){
1180       ASSIGNARRAY(ii);
1181    }
1182
1183    // Set Individual Cuts
1184    for(Int_t ii=0;ii<kNCuts;ii++){
1185       if(!SetCut(cutIds(ii),fCuts[ii]))return kFALSE;
1186    }
1187
1188    //PrintCuts();
1189
1190    // Set StandardTriggers
1191
1192    if(fIsHeavyIon)SelectCollisionCandidates(AliVEvent::kMB | AliVEvent::kCentral | AliVEvent::kSemiCentral);
1193    else SelectCollisionCandidates(AliVEvent::kMB);
1194
1195    return kTRUE;
1196 }
1197 ///________________________________________________________________________
1198 Bool_t AliConversionCuts::SetCut(cutIds cutID, const Int_t value) {
1199    ///Set individual cut ID
1200
1201
1202    switch (cutID) {
1203  
1204    case kv0FinderType:
1205       if( SetV0Finder(value)) {
1206          fCuts[kv0FinderType] = value;
1207          UpdateCutString(cutID, value);
1208          return kTRUE;
1209       } else return kFALSE;
1210
1211    case kededxSigmaCut:
1212       if( SetTPCdEdxCutElectronLine(value)) {
1213          fCuts[kededxSigmaCut] = value;
1214          UpdateCutString(cutID, value);
1215          return kTRUE;
1216       } else return kFALSE;
1217
1218    case kpidedxSigmaCut:
1219       if( SetTPCdEdxCutPionLine(value)) {
1220          fCuts[kpidedxSigmaCut] = value;
1221          UpdateCutString(cutID, value);
1222          return kTRUE;
1223       } else return kFALSE;
1224
1225    case kpiMomdedxSigmaCut:
1226       if( SetMinMomPiondEdxCut(value)) {
1227          fCuts[kpiMomdedxSigmaCut] = value;
1228          UpdateCutString(cutID, value);
1229          return kTRUE;
1230       } else return kFALSE;
1231
1232    case kchi2GammaCut:
1233       if( SetChi2GammaCut(value)) {
1234          fCuts[kchi2GammaCut] = value;
1235          UpdateCutString(cutID, value);
1236          return kTRUE;
1237       } else return kFALSE;
1238
1239    case ksinglePtCut:
1240       if( SetSinglePtCut(value)) {
1241          fCuts[ksinglePtCut] = value;
1242          UpdateCutString(cutID, value);
1243          return kTRUE;
1244       } else return kFALSE;
1245
1246    case kclsTPCCut:
1247       if( SetTPCClusterCut(value)) {
1248          fCuts[kclsTPCCut] = value;
1249          UpdateCutString(cutID, value);
1250          return kTRUE;
1251       } else return kFALSE;
1252
1253    case ketaCut:
1254       if( SetEtaCut(value)) {
1255          fCuts[ketaCut] = value;
1256          UpdateCutString(cutID, value);
1257          return kTRUE;
1258       } else return kFALSE;
1259
1260    case kLowPRejectionSigmaCut:
1261       if( SetLowPRejectionCuts(value)) {
1262          fCuts[kLowPRejectionSigmaCut] = value;
1263          UpdateCutString(cutID, value);
1264          return kTRUE;
1265       } else return kFALSE;
1266
1267    case kQtMaxCut:
1268       if( SetQtMaxCut(value)) {
1269          fCuts[kQtMaxCut] = value;
1270          UpdateCutString(cutID, value);
1271          return kTRUE;
1272       } else return kFALSE;
1273
1274    case kpiMaxMomdedxSigmaCut:
1275       if( SetMaxMomPiondEdxCut(value)) {
1276          fCuts[kpiMaxMomdedxSigmaCut] = value;
1277          UpdateCutString(cutID, value);
1278          return kTRUE;
1279       } else return kFALSE;
1280
1281    case kRCut:
1282       if( SetRCut(value)) {
1283          fCuts[kRCut] = value;
1284          UpdateCutString(cutID, value);
1285          return kTRUE;
1286       } else return kFALSE;
1287
1288    case kremovePileUp:
1289       if( SetRemovePileUp(value)) {
1290          fCuts[kremovePileUp] = value;
1291          UpdateCutString(cutID, value);
1292          return kTRUE;
1293       } else return kFALSE;
1294
1295    case kselectV0AND:
1296       if( SetSelectV0AND(value)) {
1297          fCuts[kselectV0AND] = value;
1298          UpdateCutString(cutID, value);
1299          return kTRUE;
1300       } else return kFALSE;
1301
1302    case kmultiplicityMethod:
1303       if( SetMultiplicityMethod(value)) {
1304          fCuts[kmultiplicityMethod] = value;
1305          UpdateCutString(cutID, value);
1306          return kTRUE;
1307       } else return kFALSE;
1308
1309    case kisHeavyIon:
1310       if( SetIsHeavyIon(value)) {
1311          fCuts[kisHeavyIon] = value;
1312          UpdateCutString(cutID, value);
1313          return kTRUE;
1314       } else return kFALSE;
1315
1316    case kCentralityMin:
1317       if( SetCentralityMin(value)) {
1318          fCuts[kCentralityMin] = value;
1319          UpdateCutString(cutID, value);
1320          return kTRUE;
1321       } else return kFALSE;
1322
1323    case kCentralityMax:
1324       if( SetCentralityMax(value)) {
1325          fCuts[kCentralityMax] = value;
1326          UpdateCutString(cutID, value);
1327          return kTRUE;
1328       } else return kFALSE;
1329
1330    case kTOFelectronPID:
1331       if( SetTOFElectronPIDCut(value)) {
1332          fCuts[kTOFelectronPID] = value;
1333          UpdateCutString(cutID, value);
1334          return kTRUE;
1335       } else return kFALSE;
1336
1337    case kdoPhotonAsymmetryCut:
1338       if( SetPhotonAsymmetryCut(value)) {
1339          fCuts[kdoPhotonAsymmetryCut] = value;
1340          UpdateCutString(cutID, value);
1341          return kTRUE;
1342       } else return kFALSE;
1343
1344    case kPsiPair:
1345       if( SetPsiPairCut(value)) {
1346          fCuts[kPsiPair] = value;
1347          UpdateCutString(cutID, value);
1348          return kTRUE;
1349       } else return kFALSE;
1350
1351    case kCosPAngle:
1352       if( SetCosPAngleCut(value)) {
1353          fCuts[kCosPAngle] = value;
1354          UpdateCutString(cutID, value);
1355          return kTRUE;
1356       } else return kFALSE;
1357
1358
1359    case kElecShare:
1360       if( SetSharedElectronCut(value)) {
1361          fCuts[kElecShare] = value;
1362          UpdateCutString(cutID, value);
1363          return kTRUE;
1364       } else return kFALSE;
1365
1366
1367    case kToCloseV0s:
1368       if( SetToCloseV0sCut(value)) {
1369          fCuts[kToCloseV0s] = value;
1370          UpdateCutString(cutID, value);
1371          return kTRUE;
1372       } else return kFALSE;
1373
1374    case kExtraSignals:
1375       if( SetRejectExtraSignalsCut(value)) {
1376          fCuts[kExtraSignals] = value;
1377          UpdateCutString(cutID, value);
1378          return kTRUE;
1379       } else return kFALSE;
1380
1381
1382    case kNCuts:
1383       AliError("Cut id out of range");
1384       return kFALSE;
1385    }
1386
1387    AliError("Cut id %d not recognized");
1388    return kFALSE;
1389
1390 }
1391 ///________________________________________________________________________
1392 Bool_t AliConversionCuts::SetRemovePileUp(Int_t removePileUp)
1393 {// Set Cut
1394    switch(removePileUp){
1395    case 0:
1396       fRemovePileUp=kFALSE;
1397       break;
1398    case 1:
1399       fRemovePileUp=kTRUE;
1400       break;
1401    default:
1402       AliError("RemovePileUpCut not defined");
1403       return kFALSE;
1404    }
1405    return kTRUE;
1406 }
1407
1408 ///________________________________________________________________________
1409 Bool_t AliConversionCuts::SetSelectV0AND(Int_t selectV0AND)
1410 {// Set Cut
1411    switch(selectV0AND){
1412    case 0:
1413       fSelectV0AND=kFALSE;
1414       break;
1415    case 1:
1416       fSelectV0AND=kTRUE;
1417       break;
1418    default:
1419       AliError("Warning: V0ANDCut not defined");
1420       return kFALSE;
1421    }
1422    return kTRUE;
1423 }
1424
1425 ///________________________________________________________________________
1426 Bool_t AliConversionCuts::SetMultiplicityMethod(Int_t multiplicityMethod)
1427 {
1428    // Set Cut
1429    fMultiplicityMethod=multiplicityMethod;
1430
1431    // 0 Photon Multiplicity
1432    // 1 TPC Track multiplicity
1433    // 2 V0 Mult
1434    // 3 SPD Mult
1435
1436    return kTRUE;
1437 }
1438
1439
1440 ///________________________________________________________________________
1441 void AliConversionCuts::PrintCuts() {
1442    // Print out current Cut Selection
1443    for(Int_t ic = 0; ic < kNCuts; ic++) {
1444       printf("%-30s : %d \n", fgkCutNames[ic], fCuts[ic]);
1445    }
1446
1447 }
1448
1449 ///________________________________________________________________________
1450 Bool_t AliConversionCuts::SetRCut(Int_t RCut){
1451    // Set Cut
1452    switch(RCut){
1453    case 0:
1454       fMinR=0;
1455       fMaxR = 180.;
1456       break;
1457    case 1:
1458       fMinR=2.8;
1459       fMaxR = 180.;
1460       break;
1461    case 2:
1462       fMinR=5.;
1463       fMaxR = 180.;
1464       break;
1465    case 3:
1466       fMaxR = 70.;
1467       fMinR = 10.;
1468       break;
1469    case 4:
1470       fMaxR = 70.;
1471       fMinR = 5.;
1472       break;
1473       // High purity cuts for PbPb (remove first layers of material)
1474    case 5:
1475       fMaxR = 180.;
1476       fMinR = 10.;  
1477       break;
1478    case 6:
1479       fMaxR = 180.;
1480       fMinR = 20.;
1481       break;
1482    case 7:
1483       fMaxR = 180.;
1484       fMinR = 26.;
1485       break;
1486    case 8:
1487       fMaxR = 180.;
1488       fMinR = 35.;
1489       break;
1490    case 9:
1491       fMaxR = 180.;
1492       fMinR = 60.;
1493       break;
1494
1495    default:
1496       AliError("RCut not defined");
1497       return kFALSE;
1498    }
1499    return kTRUE;
1500 }
1501
1502 ///________________________________________________________________________
1503 Bool_t AliConversionCuts::SetTPCdEdxCutElectronLine(Int_t ededxSigmaCut)
1504 {   // Set Cut
1505    switch(ededxSigmaCut){
1506    case 0: // -10,10
1507       fPIDnSigmaBelowElectronLine=-10;
1508       fPIDnSigmaAboveElectronLine=10;
1509       break;
1510    case 1: // -5,5
1511       fPIDnSigmaBelowElectronLine=-5;
1512       fPIDnSigmaAboveElectronLine=5;
1513       break;
1514    case 2: // -3,5
1515       fPIDnSigmaBelowElectronLine=-3;
1516       fPIDnSigmaAboveElectronLine=5;
1517       break;
1518    case 3: // -4,5
1519       fPIDnSigmaBelowElectronLine=-4;
1520       fPIDnSigmaAboveElectronLine=5;
1521       break;
1522    case 4: // -6,7
1523       fPIDnSigmaBelowElectronLine=-6;
1524       fPIDnSigmaAboveElectronLine=7;
1525       break;
1526    case 5: // -4,4
1527       fPIDnSigmaBelowElectronLine=-4;
1528       fPIDnSigmaAboveElectronLine=4;
1529       break;
1530    case 6: // -2.5,4
1531       fPIDnSigmaBelowElectronLine=-2.5;
1532       fPIDnSigmaAboveElectronLine=4;
1533       break;
1534    case 7: // -2,3.5
1535       fPIDnSigmaBelowElectronLine=-2;
1536       fPIDnSigmaAboveElectronLine=3.5;
1537       break;
1538    default:
1539       AliError("TPCdEdxCutElectronLine not defined");
1540       return kFALSE;
1541         
1542    }
1543    return kTRUE;
1544 }
1545
1546 ///________________________________________________________________________
1547 Bool_t AliConversionCuts::SetTPCdEdxCutPionLine(Int_t pidedxSigmaCut)
1548 {   // Set Cut
1549
1550    switch(pidedxSigmaCut){
1551    case 0:  // -10
1552       fPIDnSigmaAbovePionLine=-10;
1553       fPIDnSigmaAbovePionLineHighPt=-10;
1554       break;
1555    case 1:   // 0
1556       fPIDnSigmaAbovePionLine=0;
1557       fPIDnSigmaAbovePionLineHighPt=-10;
1558       break;
1559    case 2:  // 1
1560       fPIDnSigmaAbovePionLine=1;
1561       fPIDnSigmaAbovePionLineHighPt=-10;
1562       break;
1563    case 3:  // 1
1564       fPIDnSigmaAbovePionLine=-1;
1565       fPIDnSigmaAbovePionLineHighPt=-10;
1566       break;
1567    case 4:  // 1
1568       fPIDnSigmaAbovePionLine=2.5;
1569       fPIDnSigmaAbovePionLineHighPt=-10;
1570       break;
1571    case 5:  // 1
1572       fPIDnSigmaAbovePionLine=2.;
1573       fPIDnSigmaAbovePionLineHighPt=-10;
1574       break;
1575    case 6:  // 1
1576       fPIDnSigmaAbovePionLine=2.;
1577       fPIDnSigmaAbovePionLineHighPt=0.5;
1578       break;
1579    case 7:  // 1
1580       fPIDnSigmaAbovePionLine=3.5;
1581       fPIDnSigmaAbovePionLineHighPt=-10;
1582       break;
1583    case 8:  // 1
1584       fPIDnSigmaAbovePionLine=2.;
1585       fPIDnSigmaAbovePionLineHighPt=1.;
1586       break;
1587    case 9:
1588       fPIDnSigmaAbovePionLine=3.0; // We need a bit less tight cut on dE/dx
1589       fPIDnSigmaAbovePionLineHighPt=-10;
1590       break;
1591    default:
1592       AliError(Form("Warning: pidedxSigmaCut not defined %d",pidedxSigmaCut));
1593       return kFALSE;
1594    }
1595    return kTRUE;
1596 }
1597
1598 ///________________________________________________________________________
1599 Bool_t AliConversionCuts::SetMinMomPiondEdxCut(Int_t piMomdedxSigmaCut)
1600 {   // Set Cut
1601    switch(piMomdedxSigmaCut){
1602    case 0:  // 0.5 GeV
1603       fPIDMinPnSigmaAbovePionLine=0.5;
1604       break;
1605    case 1:  // 1. GeV
1606       fPIDMinPnSigmaAbovePionLine=1.;
1607       break;
1608    case 2:  // 1.5 GeV
1609       fPIDMinPnSigmaAbovePionLine=1.5;
1610       break;
1611    case 3:  // 20.0 GeV
1612       fPIDMinPnSigmaAbovePionLine=20.;
1613       break;
1614    case 4:  // 50.0 GeV
1615       fPIDMinPnSigmaAbovePionLine=50.;
1616       break;
1617    case 5:  // 0.3 GeV
1618       fPIDMinPnSigmaAbovePionLine=0.3;
1619       break;
1620    case 6:  // 0.25 GeV
1621       fPIDMinPnSigmaAbovePionLine=0.25;
1622       break;
1623    case 7:  // 0.4 GeV
1624       fPIDMinPnSigmaAbovePionLine=0.4;
1625       break;
1626    default:
1627       AliError(Form("piMomdedxSigmaCut not defined %d",piMomdedxSigmaCut));
1628       return kFALSE;
1629    }
1630    return kTRUE;
1631 }
1632
1633 ///________________________________________________________________________
1634 Bool_t AliConversionCuts::SetChi2GammaCut(Int_t chi2GammaCut)
1635 {   // Set Cut
1636
1637    switch(chi2GammaCut){
1638    case 0: // 100
1639       fChi2CutConversion = 100.;
1640       break;
1641    case 1:  // 50
1642       fChi2CutConversion = 50.;
1643       break;
1644    case 2:  // 30
1645       fChi2CutConversion = 30.;
1646       break;
1647    case 3:
1648       fChi2CutConversion = 200.;
1649       break;
1650    case 4:
1651       fChi2CutConversion = 500.;
1652       break;
1653    case 5:
1654       fChi2CutConversion = 100000.;
1655       break;
1656    case 6:
1657       fChi2CutConversion = 5.;
1658       break;
1659    case 7:
1660       fChi2CutConversion = 10.;
1661       break;
1662    case 8:
1663       fChi2CutConversion = 20.;
1664       break;
1665    case 9:
1666       fChi2CutConversion = 15.;
1667       break;
1668    default:
1669       AliError(Form("Warning: Chi2GammaCut not defined %d",chi2GammaCut));
1670       return kFALSE;
1671    }
1672    return kTRUE;
1673 }
1674
1675 ///________________________________________________________________________
1676 Bool_t AliConversionCuts::SetV0Finder(Int_t v0FinderType)
1677 {   // Set Cut
1678    switch (v0FinderType){
1679    case 0:  // on fly V0 finder
1680       fUseOnFlyV0Finder=kTRUE;
1681       break;
1682    case 1:  // offline V0 finder
1683       fUseOnFlyV0Finder=kFALSE;
1684       break;
1685    default:
1686       AliError(Form(" v0FinderType not defined %d",v0FinderType));
1687       return kFALSE;
1688    }
1689    return kTRUE;
1690 }
1691
1692
1693 ///________________________________________________________________________
1694 Bool_t AliConversionCuts::SetSinglePtCut(Int_t singlePtCut)
1695 {   // Set Cut
1696    switch(singlePtCut){
1697    case 0: // 0.050 GeV
1698       fSinglePtCut = 0.050;
1699       break;
1700    case 1:  // 0.100 GeV
1701       fSinglePtCut = 0.100;
1702       break;
1703    case 2:  // 0.150 GeV
1704       fSinglePtCut = 0.150;
1705       break;
1706    case 3:  // 0.200 GeV
1707       fSinglePtCut = 0.200;
1708       break;
1709    case 4:  // 0.075 GeV
1710       fSinglePtCut = 0.075;
1711       break;
1712    case 5:  // 0.125 GeV
1713       fSinglePtCut = 0.125;
1714       break;
1715    default:
1716       AliError(Form("singlePtCut not defined %d",singlePtCut));
1717       return kFALSE;
1718    }
1719    return kTRUE;
1720 }
1721 ///________________________________________________________________________
1722 Bool_t AliConversionCuts::SetTPCClusterCut(Int_t clsTPCCut)
1723 {   // Set Cut
1724    switch(clsTPCCut){
1725    case 0: // 0
1726       fMinClsTPC= 0.;
1727       break;
1728    case 1:  // 70
1729       fMinClsTPC= 70.;
1730       break;
1731    case 2:  // 80
1732       fMinClsTPC= 80.;
1733       break;
1734    case 3:  // 100
1735       fMinClsTPC= 100.;
1736       break;
1737    case 4:  // 60% of findable clusters
1738       fMinClsTPCToF= 0.6;
1739       fUseCorrectedTPCClsInfo=0;
1740       break;
1741    case 5:  // 0% of findable clusters
1742       fMinClsTPCToF= 0.0;
1743       fUseCorrectedTPCClsInfo=1;
1744       break;
1745    case 6:  // 0% of findable clusters
1746       fMinClsTPCToF= 0.7;
1747       fUseCorrectedTPCClsInfo=1;
1748       break;
1749    case 7:  // 0% of findable clusters
1750       fMinClsTPCToF= 0.35;
1751       fUseCorrectedTPCClsInfo=0;
1752       break;
1753    case 8:
1754       fMinClsTPCToF= 0.35;
1755       fUseCorrectedTPCClsInfo=1;
1756       break;
1757    case 9:
1758       fMinClsTPCToF= 0.6;
1759       fUseCorrectedTPCClsInfo=1;
1760       break;
1761    default:
1762       AliError(Form("Warning: clsTPCCut not defined %d",clsTPCCut));
1763       return kFALSE;
1764    }
1765    return kTRUE;
1766 }
1767
1768 ///________________________________________________________________________
1769 Bool_t AliConversionCuts::SetEtaCut(Int_t etaCut)
1770 {   // Set Cut
1771
1772    //Set Standard LineCutZValues
1773    fLineCutZValueMin = -2;
1774    fLineCutZValue = 7.;
1775    
1776    switch(etaCut){
1777    case 0: // 0.9
1778       fEtaCut           = 0.9;
1779       fLineCutZRSlope = tan(2*atan(exp(-fEtaCut)));
1780       fEtaCutMin                = -0.1;
1781       fLineCutZRSlopeMin = 0.;
1782       break;
1783    case 1:      // 1.2
1784       fEtaCut           = 1.2;
1785       fLineCutZRSlope = tan(2*atan(exp(-fEtaCut)));
1786       fEtaCutMin                = -0.1;
1787       fLineCutZRSlopeMin = 0.;
1788       break;
1789    case 2:      // 1.4
1790       fEtaCut           = 1.4;
1791       fLineCutZRSlope = tan(2*atan(exp(-fEtaCut)));
1792       fEtaCutMin                = -0.1;
1793       fLineCutZRSlopeMin = tan(2*atan(exp(-fEtaCut)));
1794       break;
1795    case 3: // 0.8
1796       fEtaCut           = 0.8;
1797       fLineCutZRSlope = tan(2*atan(exp(-fEtaCut)));
1798       fEtaCutMin                = -0.1;
1799       fLineCutZRSlopeMin = 0.;
1800       break;
1801    case 4: // 0.75
1802       fEtaCut           = 0.75;
1803       fLineCutZRSlope = tan(2*atan(exp(-fEtaCut)));
1804       fEtaCutMin                = -0.1;
1805       fLineCutZRSlopeMin = 0.;
1806       break;
1807    case 5: // 0.9 - 1.4
1808       fEtaCut           = 1.4;
1809       fLineCutZRSlope = tan(2*atan(exp(-fEtaCut)));
1810       fEtaCutMin                = 0.9;
1811       fLineCutZRSlopeMin = tan(2*atan(exp(-fEtaCutMin)));
1812       break;
1813    case 6: // 0.9 - 1.2
1814       fEtaCut           = 1.2;
1815       fLineCutZRSlope = tan(2*atan(exp(-fEtaCut)));
1816       fEtaCutMin                = 0.9;
1817       fLineCutZRSlopeMin = tan(2*atan(exp(-fEtaCutMin)));
1818       break;
1819    case 7: // 0.1 - 0.8
1820       fEtaCut           = 0.8;
1821       fLineCutZRSlope = tan(2*atan(exp(-fEtaCut)));
1822       fEtaCutMin                = 0.1;
1823       fLineCutZRSlopeMin = tan(2*atan(exp(-fEtaCutMin)));
1824       break;
1825    case 8: // 0.1 - 0.8
1826       fEtaCut           = 0.9;
1827       fLineCutZRSlope = tan(2*atan(exp(-fEtaCut)));
1828       fEtaCutMin                = 0.1;
1829       fLineCutZRSlopeMin = tan(2*atan(exp(-fEtaCutMin)));
1830       break;
1831    default:
1832       AliError(Form(" EtaCut not defined %d",etaCut));
1833       return kFALSE;
1834    }
1835    return kTRUE;
1836 }
1837
1838 ///________________________________________________________________________
1839 Bool_t AliConversionCuts::SetMaxMomPiondEdxCut(Int_t piMaxMomdedxSigmaCut)
1840 {   // Set Cut
1841    switch(piMaxMomdedxSigmaCut){
1842    case 0:  // 100. GeV
1843       fPIDMaxPnSigmaAbovePionLine=100.;
1844       break;
1845    case 1:  // 5. GeV
1846       fPIDMaxPnSigmaAbovePionLine=5.;
1847       break;
1848    case 2:  // 4. GeV
1849       fPIDMaxPnSigmaAbovePionLine=4.;
1850       break;
1851    case 3:  // 3.5 GeV
1852       fPIDMaxPnSigmaAbovePionLine=3.5;
1853       break;
1854    case 4:  // 3. GeV
1855       fPIDMaxPnSigmaAbovePionLine=3.;
1856       break;
1857    default:
1858       AliError(Form("piMaxMomdedxSigmaCut not defined %d",piMaxMomdedxSigmaCut));
1859       return kFALSE;
1860    }
1861    return kTRUE;
1862 }
1863
1864 ///________________________________________________________________________
1865 Bool_t AliConversionCuts::SetIsHeavyIon(Int_t isHeavyIon)
1866 {   // Set Cut
1867    switch(isHeavyIon){
1868    case 0:
1869       fIsHeavyIon=0;
1870       break;
1871    case 1:
1872       fIsHeavyIon=1;
1873       fDetectorCentrality=0;
1874       break;
1875    case 2:
1876       fIsHeavyIon=1;
1877       fDetectorCentrality=1;
1878       break;
1879    case 3: //allows to select centrality 0-45% in steps of 5% for V0 Multiplicity
1880       fIsHeavyIon=1;
1881       fDetectorCentrality=0;      
1882       fModCentralityClass=1;
1883       break;
1884    case 4: //allows to select centrality 45-90% in steps of 5% for V0 Multiplicity
1885       fIsHeavyIon=1;
1886       fDetectorCentrality=0;
1887       fModCentralityClass=2;
1888       break;
1889    default:
1890       AliError(Form("SetHeavyIon not defined %d",isHeavyIon));
1891       return kFALSE;
1892    }
1893    return kTRUE;
1894 }
1895
1896
1897 ///________________________________________________________________________
1898 Bool_t AliConversionCuts::SetLowPRejectionCuts(Int_t LowPRejectionSigmaCut)
1899 {   // Set Cut
1900    switch(LowPRejectionSigmaCut){
1901    case 0:  //
1902       fPIDnSigmaAtLowPAroundKaonLine=0;
1903       fPIDnSigmaAtLowPAroundProtonLine=0;
1904       fPIDnSigmaAtLowPAroundPionLine=0;
1905       break;
1906    case 1:  //
1907       fPIDnSigmaAtLowPAroundKaonLine=0.5;
1908       fPIDnSigmaAtLowPAroundProtonLine=0.5;
1909       fPIDnSigmaAtLowPAroundPionLine=0.5;
1910       break;
1911    case 2:  //
1912       fPIDnSigmaAtLowPAroundKaonLine=1;
1913       fPIDnSigmaAtLowPAroundProtonLine=1;
1914       fPIDnSigmaAtLowPAroundPionLine=1;
1915       break;
1916    case 3:  //
1917       fPIDnSigmaAtLowPAroundKaonLine=2.;
1918       fPIDnSigmaAtLowPAroundProtonLine=2.;
1919       fPIDnSigmaAtLowPAroundPionLine=2.;
1920       break;
1921    case 4:  //
1922       fPIDnSigmaAtLowPAroundKaonLine=0.;
1923       fPIDnSigmaAtLowPAroundProtonLine=0.;
1924       fPIDnSigmaAtLowPAroundPionLine=1;
1925       break;
1926    case 5:  //
1927       fPIDnSigmaAtLowPAroundKaonLine=0.;
1928       fPIDnSigmaAtLowPAroundProtonLine=0.;
1929       fPIDnSigmaAtLowPAroundPionLine=1.5;
1930       break;
1931    case 6:  //
1932       fPIDnSigmaAtLowPAroundKaonLine=0.;
1933       fPIDnSigmaAtLowPAroundProtonLine=0.;
1934       fPIDnSigmaAtLowPAroundPionLine=2.;
1935       break;
1936    default:
1937       AliError(Form("LowPRejectionSigmaCut not defined %d",LowPRejectionSigmaCut));
1938       return kFALSE;
1939    }
1940    return kTRUE;
1941 }
1942
1943 ///________________________________________________________________________
1944 Bool_t AliConversionCuts::SetTOFElectronPIDCut(Int_t TOFelectronPID){
1945    // Set Cut
1946    switch(TOFelectronPID){ 
1947    case 0: // no cut
1948       fUseTOFpid = kFALSE;
1949       fTofPIDnSigmaBelowElectronLine=-100;
1950       fTofPIDnSigmaAboveElectronLine=100;
1951       break;
1952    case 1: // -7,7
1953       fUseTOFpid = kTRUE;
1954       fTofPIDnSigmaBelowElectronLine=-7;
1955       fTofPIDnSigmaAboveElectronLine=7;
1956       break;
1957    case 2: // -5,5
1958       fUseTOFpid = kTRUE;
1959       fTofPIDnSigmaBelowElectronLine=-5;
1960       fTofPIDnSigmaAboveElectronLine=5;
1961       break;
1962    case 3: // -3,5
1963       fUseTOFpid = kTRUE;
1964       fTofPIDnSigmaBelowElectronLine=-3;
1965       fTofPIDnSigmaAboveElectronLine=5;
1966       break;
1967    case 4: // -2,3
1968       fUseTOFpid = kTRUE;
1969       fTofPIDnSigmaBelowElectronLine=-2;
1970       fTofPIDnSigmaAboveElectronLine=3;
1971       break;
1972    default:
1973       AliError(Form("TOFElectronCut not defined %d",TOFelectronPID));
1974       return kFALSE;
1975    } 
1976    return kTRUE;
1977 }
1978
1979 ///________________________________________________________________________
1980 Bool_t AliConversionCuts::SetTRDElectronCut(Int_t TRDElectronCut)
1981 {   // Set Cut
1982    switch(TRDElectronCut){
1983    case 0:
1984       fDoTRDPID=kFALSE;
1985       break;
1986    case 1:
1987       fDoTRDPID=kTRUE;
1988       fPIDTRDEfficiency=0.1;
1989       break;
1990    case 8:
1991       fDoTRDPID=kTRUE;
1992       fPIDTRDEfficiency=0.8;
1993       break;
1994    case 9:
1995       fDoTRDPID=kTRUE;
1996       fPIDTRDEfficiency=0.9;
1997       break;
1998    default:
1999       AliError(Form("TRDElectronCut not defined %d",TRDElectronCut));
2000       return kFALSE;
2001    }
2002
2003    return kTRUE;
2004 }
2005 ///________________________________________________________________________
2006 Bool_t AliConversionCuts::SetQtMaxCut(Int_t QtMaxCut)
2007 {   // Set Cut
2008    switch(QtMaxCut){
2009    case 0: //
2010       fQtMax=1.;
2011       fDoQtGammaSelection=kFALSE;      //No Qt selection (true by default)
2012       fDoHighPtQtGammaSelection=kFALSE; // RRnew
2013       fHighPtQtMax=100.;                // RRnew
2014       fPtBorderForQt=100.;              // RRnew
2015       break;
2016    case 1:
2017       fQtMax=0.1;
2018       fDoHighPtQtGammaSelection=kFALSE; // RRnew
2019       fHighPtQtMax=100.;                // RRnew
2020       fPtBorderForQt=100.;              // RRnew
2021       break;
2022    case 2:
2023       fQtMax=0.07;
2024       fDoHighPtQtGammaSelection=kFALSE; // RRnew
2025       fHighPtQtMax=100.;                // RRnew
2026       fPtBorderForQt=100.;              // RRnew
2027       break;
2028    case 3:
2029       fQtMax=0.05;
2030       fDoHighPtQtGammaSelection=kFALSE; // RRnew
2031       fHighPtQtMax=100.;                // RRnew
2032       fPtBorderForQt=100.;              // RRnew
2033       break;
2034    case 4:
2035       fQtMax=0.03;
2036       fDoHighPtQtGammaSelection=kFALSE; // RRnew
2037       fHighPtQtMax=100.;                // RRnew
2038       fPtBorderForQt=100.;              // RRnew
2039       break;
2040    case 5: // RR try to improve (get rid of) low InvMass peak in PbPb
2041       fQtMax=0.02;
2042       fDoHighPtQtGammaSelection=kFALSE; // RRnew
2043       fHighPtQtMax=100.;                // RRnew
2044       fPtBorderForQt=100.;              // RRnew
2045       break; // end RR ///////////////////////////////////////////////
2046    case 6:  // RRnew start: pT dependent qT cut
2047       fQtMax=0.02;
2048       fDoHighPtQtGammaSelection=kTRUE;
2049       fHighPtQtMax=0.06;
2050       fPtBorderForQt=2.5;
2051       break; // RRnew end ////////////////////////////////////////////
2052    case 7:
2053       fQtMax=0.15;
2054       fDoHighPtQtGammaSelection=kFALSE; // RRnew
2055       fHighPtQtMax=100.;                // RRnew
2056       fPtBorderForQt=100.;              // RRnew
2057       break;
2058    default:
2059       AliError(Form("Warning: QtMaxCut not defined %d",QtMaxCut));
2060       return kFALSE;
2061    }
2062    return kTRUE;
2063 }
2064
2065 //-------------------------------------------------------------
2066 Double_t AliConversionCuts::GetCentrality(AliVEvent *event)
2067 {   // Get Event Centrality
2068
2069    AliESDEvent *esdEvent=dynamic_cast<AliESDEvent*>(event);
2070    if(esdEvent){
2071       AliCentrality *fESDCentrality=(AliCentrality*)esdEvent->GetCentrality();
2072
2073       if(fDetectorCentrality==0){
2074          return fESDCentrality->GetCentralityPercentile("V0M"); // default
2075       }
2076       if(fDetectorCentrality==1){
2077          return fESDCentrality->GetCentralityPercentile("CL1");
2078       }
2079    }
2080
2081    AliAODEvent *aodEvent=dynamic_cast<AliAODEvent*>(event);
2082    if(aodEvent){
2083       if(aodEvent->GetHeader()){return aodEvent->GetHeader()->GetCentrality();}
2084    }
2085
2086    return -1;
2087 }
2088
2089 //-------------------------------------------------------------
2090 Bool_t AliConversionCuts::IsCentralitySelected(AliVEvent *event)
2091 {   // Centrality Selection
2092    if(!fIsHeavyIon)return kTRUE;
2093     
2094    if(fCentralityMin == 0 && fCentralityMax == 0) return kTRUE;//0-100%
2095    if(fCentralityMin >= fCentralityMax) return kTRUE;//0-100%
2096    
2097    Double_t centrality=GetCentrality(event);
2098    if(centrality<0)return kFALSE;
2099         
2100    Int_t centralityC=0;
2101    if (fModCentralityClass == 0){
2102       centralityC= Int_t(centrality/10);
2103       if(centralityC >= fCentralityMin && centralityC < fCentralityMax)
2104          return kTRUE;
2105       else 
2106          return kFALSE;
2107    }    else if (fModCentralityClass ==1){
2108       centralityC= Int_t(centrality);
2109       if(centralityC >= fCentralityMin*5 && centralityC < fCentralityMax*5){
2110          return kTRUE;
2111       } else { 
2112          return kFALSE;
2113       }
2114    } else if (fModCentralityClass ==2){
2115       centralityC= Int_t(centrality+1);
2116       if(centralityC >= (fCentralityMin*5+45) && centralityC < (fCentralityMax*5+45))
2117          return kTRUE;
2118       else 
2119          return kFALSE;
2120    }
2121    return kFALSE;
2122 }
2123
2124 //-------------------------------------------------------------
2125 Bool_t AliConversionCuts::SetCentralityMin(Int_t minCentrality)
2126 {
2127    // Set Cut
2128    if(minCentrality<0||minCentrality>9){
2129       AliError(Form("minCentrality not defined %d",minCentrality));
2130       return kFALSE;
2131    }
2132
2133    fCentralityMin=minCentrality;
2134    return kTRUE;
2135 }
2136 //-------------------------------------------------------------
2137 Bool_t AliConversionCuts::SetCentralityMax(Int_t maxCentrality)
2138 {
2139    // Set Cut
2140    if(maxCentrality<0||maxCentrality>9){
2141       AliError(Form("maxCentrality not defined %d",maxCentrality));
2142       return kFALSE;
2143    }
2144
2145    fCentralityMax=maxCentrality;
2146    return kTRUE;
2147 }
2148
2149 ///________________________________________________________________________
2150 Bool_t AliConversionCuts::SetPhotonAsymmetryCut(Int_t doPhotonAsymmetryCut){
2151    // Set Cut
2152    switch(doPhotonAsymmetryCut){
2153    case 0:
2154       fDoPhotonAsymmetryCut=0;
2155       fMinPPhotonAsymmetryCut=100.;
2156       fMinPhotonAsymmetry=0.;
2157       break;
2158    case 1:
2159       fDoPhotonAsymmetryCut=1;
2160       fMinPPhotonAsymmetryCut=3.5;
2161       fMinPhotonAsymmetry=0.04;
2162       break;
2163    case 2:
2164       fDoPhotonAsymmetryCut=1;
2165       fMinPPhotonAsymmetryCut=3.5;
2166       fMinPhotonAsymmetry=0.06;
2167       break;
2168    default:
2169       AliError(Form("PhotonAsymmetryCut not defined %d",doPhotonAsymmetryCut));
2170       return kFALSE;
2171    }
2172    fCuts[kdoPhotonAsymmetryCut]=doPhotonAsymmetryCut;
2173    return kTRUE;
2174 }
2175
2176
2177 ///________________________________________________________________________
2178 Bool_t AliConversionCuts::SetPsiPairCut(Int_t psiCut) {
2179   
2180
2181    switch(psiCut) {
2182    case 0:
2183       fPsiPairCut = 10000; // 
2184       break;
2185    case 1:
2186       fPsiPairCut = 0.1; // 
2187       break;
2188    case 2:
2189       fPsiPairCut = 0.05; // Standard
2190       break;
2191    case 3:
2192       fPsiPairCut = 0.035; // 
2193       break;
2194    case 4:
2195       fPsiPairCut = 0.15; // 
2196       break;
2197    case 5:
2198       fPsiPairCut = 0.2; // 
2199       break;
2200    case 6:
2201       fPsiPairCut = 0.03; // 
2202       break;
2203    case 7:
2204       fPsiPairCut = 0.025; // 
2205       break;
2206    case 8:
2207       fPsiPairCut = 0.01; // 
2208       break;
2209    default:
2210       AliError(Form("PsiPairCut not defined %d",psiCut));
2211       return kFALSE;
2212    }
2213
2214    return kTRUE;
2215 }
2216
2217 ///________________________________________________________________________
2218 Bool_t AliConversionCuts::SetCosPAngleCut(Int_t cosCut) {
2219
2220    switch(cosCut){
2221    case 0:
2222       fCosPAngleCut = TMath::Pi(); //
2223       break;
2224    case 1:
2225       fCosPAngleCut = 0.1; //
2226       break;
2227    case 2:
2228       fCosPAngleCut = 0.05; //
2229       break;
2230    case 3:
2231       fCosPAngleCut = 0.025; // Standard
2232       break;
2233    case 4:
2234       fCosPAngleCut = 0.01; //
2235       break;
2236    default:
2237       AliError(Form("Cosine Pointing Angle cut not defined %d",cosCut));
2238       return kFALSE;
2239    }
2240         
2241    return kTRUE;
2242 }
2243
2244
2245 ///________________________________________________________________________
2246 Bool_t AliConversionCuts::VertexZCut(AliVEvent *event){
2247    // Cut on z position of primary vertex
2248    Double_t fVertexZ=event->GetPrimaryVertex()->GetZ();
2249
2250    if(TMath::Abs(fVertexZ)>fMaxVertexZ)return kFALSE;
2251    return kTRUE;
2252 }
2253
2254 ///________________________________________________________________________
2255 Bool_t AliConversionCuts::SetSharedElectronCut(Int_t sharedElec) {
2256
2257    switch(sharedElec){
2258    case 0:
2259       fDoSharedElecCut = kFALSE;
2260       break;
2261    case 1:
2262       fDoSharedElecCut = kTRUE;
2263       break;
2264    default:
2265       AliError(Form("Shared Electron Cut not defined %d",sharedElec));
2266       return kFALSE;
2267    }
2268         
2269    return kTRUE;
2270 }
2271
2272 ///________________________________________________________________________
2273 Bool_t AliConversionCuts::SetToCloseV0sCut(Int_t toClose) {
2274
2275    switch(toClose){
2276    case 0:
2277       fDoToCloseV0sCut = kFALSE;
2278       fminV0Dist = 250;
2279       break;
2280    case 1:
2281       fDoToCloseV0sCut = kTRUE;
2282       fminV0Dist = 1;
2283       break;
2284    case 2:
2285       fDoToCloseV0sCut = kTRUE;
2286       fminV0Dist = 2;
2287       break;
2288    case 3:
2289       fDoToCloseV0sCut = kTRUE;
2290       fminV0Dist = 3;
2291       break;
2292    default:
2293       AliError(Form("Shared Electron Cut not defined %d",toClose));
2294       return kFALSE;
2295    }
2296    return kTRUE;
2297 }
2298 ///________________________________________________________________________
2299 Bool_t AliConversionCuts::SetRejectExtraSignalsCut(Int_t extraSignal) {
2300
2301    switch(extraSignal){
2302    case 0:
2303       fRejectExtraSignals = 0;
2304       break; // No Rejection
2305    case 1:
2306       fRejectExtraSignals = 1;
2307       break; // MinBias Header
2308    case 2:
2309       fRejectExtraSignals = 2;
2310       break; // User String Array
2311    default:
2312       AliError(Form("Extra Signal Rejection not defined %d",extraSignal));
2313       return kFALSE;
2314    }
2315    return kTRUE;
2316 }
2317
2318 ///________________________________________________________________________
2319
2320 Int_t AliConversionCuts::GetNumberOfContributorsVtx(AliVEvent *event){
2321    // returns number of contributors to the vertex
2322     
2323    AliESDEvent *fESDEvent=dynamic_cast<AliESDEvent*>(event);
2324    if(fESDEvent){
2325       if(fESDEvent->GetPrimaryVertexTracks()->GetNContributors()>0) {
2326          return fESDEvent->GetPrimaryVertexTracks()->GetNContributors();
2327       }
2328      
2329       if(fESDEvent->GetPrimaryVertexTracks()->GetNContributors()<1) {
2330          //             return 0;
2331          //-AM test pi0s without SPD only vertex
2332          if(fESDEvent->GetPrimaryVertexSPD()->GetNContributors()>0) {
2333             return fESDEvent->GetPrimaryVertexSPD()->GetNContributors();
2334
2335          }
2336          if(fESDEvent->GetPrimaryVertexSPD()->GetNContributors()<1) {
2337             return 0;
2338          }
2339       }
2340    }
2341
2342    AliAODEvent *fAODEvent=dynamic_cast<AliAODEvent*>(event);
2343    if(fAODEvent){
2344       if(fAODEvent->GetPrimaryVertex()->GetNContributors()>0) {
2345          return fAODEvent->GetPrimaryVertex()->GetNContributors();
2346       }
2347       if(fAODEvent->GetPrimaryVertex()->GetNContributors()<1) {
2348          if(fAODEvent->GetPrimaryVertexSPD()->GetNContributors()>0) {
2349             return fAODEvent->GetPrimaryVertexSPD()->GetNContributors();
2350          }
2351          if(fAODEvent->GetPrimaryVertexSPD()->GetNContributors()<1) {
2352             AliWarning(Form("Number of contributors from bad vertex type:: %s",fAODEvent->GetPrimaryVertex()->GetName()));
2353             return 0;
2354          }
2355       }
2356    }
2357
2358
2359    return 0;
2360 }
2361
2362 ///________________________________________________________________________
2363
2364 Bool_t AliConversionCuts::IsTriggerSelected()
2365 {
2366    AliInputEventHandler *fInputHandler=(AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
2367
2368    UInt_t isSelected = AliVEvent::kAny;
2369    if( fInputHandler && fInputHandler->GetEventSelection()) {
2370       // Get the actual offline trigger mask for the event and AND it with the
2371       // requested mask. If no mask requested select by default the event.
2372       if (fOfflineTriggerMask)
2373          isSelected = fOfflineTriggerMask & fInputHandler->IsEventSelected();
2374    }
2375    
2376    if(!isSelected)return kFALSE;
2377
2378    // Fill Histogram
2379    if(hTriggerClass){
2380       if (fInputHandler->IsEventSelected() & AliVEvent::kAny)hTriggerClass->Fill(0);
2381       if (fInputHandler->IsEventSelected() & AliVEvent::kMB)hTriggerClass->Fill(1);
2382       if (fInputHandler->IsEventSelected() & AliVEvent::kCentral)hTriggerClass->Fill(2);
2383       if (fInputHandler->IsEventSelected() & AliVEvent::kSemiCentral)hTriggerClass->Fill(3);
2384       if (fInputHandler->IsEventSelected() & AliVEvent::kCINT5)hTriggerClass->Fill(4);
2385       if (fInputHandler->IsEventSelected() & AliVEvent::kCINT5)hTriggerClass->Fill(5);
2386    }
2387
2388     
2389    return kTRUE;
2390
2391 }
2392
2393 ///________________________________________________________________________
2394 Int_t AliConversionCuts::GetFirstTPCRow(Double_t radius){
2395    // Get first TPC row
2396    Int_t firstTPCRow=0;
2397    Double_t radiusI     =       84.8;
2398    Double_t radiusO     = 134.6;
2399    Double_t radiusOB = 198.;
2400    Double_t rSizeI       = 0.75;
2401    Double_t rSizeO       = 1.;
2402    Double_t rSizeOB     = 1.5;
2403    Int_t nClsI=63;
2404    Int_t nClsIO=127;
2405
2406    if(radius <= radiusI){
2407       return firstTPCRow;
2408    }
2409    if(radius>radiusI && radius<=radiusO){
2410       firstTPCRow = (Int_t)((radius-radiusI)/rSizeI);
2411    }
2412    if(radius>radiusO && radius<=radiusOB){
2413       firstTPCRow = (Int_t)(nClsI+(radius-radiusO)/rSizeO);
2414    }
2415
2416    if(radius>radiusOB){
2417       firstTPCRow =(Int_t)(nClsIO+(radius-radiusOB)/rSizeOB);
2418    }
2419
2420    return firstTPCRow;
2421 }
2422
2423 Bool_t AliConversionCuts::CosinePAngleCut(const AliConversionPhotonBase * photon, AliVEvent * event) const {
2424    ///Check if passes cosine of pointing angle cut
2425    if(GetCosineOfPointingAngle(photon, event) < (TMath::Cos(fCosPAngleCut))){
2426       return kFALSE;
2427    }
2428    return kTRUE;
2429 }
2430
2431 Double_t AliConversionCuts::GetCosineOfPointingAngle( const AliConversionPhotonBase * photon, AliVEvent * event) const{
2432    // calculates the pointing angle of the recalculated V0 
2433
2434    Double_t momV0[3] = {0,0,0};
2435    if(event->IsA()==AliESDEvent::Class()){
2436       AliESDEvent *esdEvent = dynamic_cast<AliESDEvent*>(event);
2437       AliESDv0 *v0 = esdEvent->GetV0(photon->GetV0Index());
2438       v0->GetPxPyPz(momV0[0],momV0[1],momV0[2]);
2439    }
2440    if(event->IsA()==AliAODEvent::Class()){
2441       momV0[0] = photon->GetPx();
2442       momV0[1] = photon->GetPy();
2443       momV0[2] = photon->GetPz();
2444    }
2445    
2446    //Double_t momV0[3] = { photon->GetPx(), photon->GetPy(), photon->GetPz() }; //momentum of the V0
2447    Double_t PosV0[3] = { photon->GetConversionX() - event->GetPrimaryVertex()->GetX(), 
2448                          photon->GetConversionY() - event->GetPrimaryVertex()->GetY(), 
2449                          photon->GetConversionZ() - event->GetPrimaryVertex()->GetZ() }; //Recalculated V0 Position vector
2450    
2451    Double_t momV02 = momV0[0]*momV0[0] + momV0[1]*momV0[1] + momV0[2]*momV0[2];
2452    Double_t PosV02 = PosV0[0]*PosV0[0] + PosV0[1]*PosV0[1] + PosV0[2]*PosV0[2];
2453
2454    Double_t cosinePointingAngle = (PosV0[0]*momV0[0] +  PosV0[1]*momV0[1] + PosV0[2]*momV0[2] ) / TMath::Sqrt(momV02 * PosV02);
2455   
2456    return cosinePointingAngle;
2457 }
2458
2459
2460 Bool_t AliConversionCuts::PsiPairCut(const AliConversionPhotonBase * photon) const {
2461
2462    if(photon->GetPsiPair() > fPsiPairCut){
2463       return kFALSE;}
2464    else{return kTRUE;}
2465
2466 }
2467
2468 ///________________________________________________________________________
2469 TString AliConversionCuts::GetCutNumber(){
2470    // returns TString with current cut number
2471    TString a(kNCuts);
2472    for(Int_t ii=0;ii<kNCuts;ii++){
2473       a.Append(Form("%d",fCuts[ii]));
2474    }
2475    return a;
2476 }
2477
2478 ///________________________________________________________________________
2479 void AliConversionCuts::FillElectonLabelArray(AliAODConversionPhoton* photon, Int_t nV0){
2480    
2481    Int_t posLabel = photon->GetTrackLabelPositive();
2482    Int_t negLabel = photon->GetTrackLabelNegative();
2483    
2484    fElectronLabelArray[nV0*2] = posLabel;
2485    fElectronLabelArray[(nV0*2)+1] = negLabel;
2486 }
2487 ///________________________________________________________________________
2488 Bool_t AliConversionCuts::RejectSharedElectronV0s(AliAODConversionPhoton* photon, Int_t nV0, Int_t nV0s){
2489
2490    Int_t posLabel = photon->GetTrackLabelPositive();
2491    Int_t negLabel = photon->GetTrackLabelNegative();
2492    
2493    for(Int_t i = 0; i<nV0s*2;i++){
2494       if(i==nV0*2)     continue;
2495       if(i==(nV0*2)+1) continue;
2496       if(fElectronLabelArray[i] == posLabel){
2497          return kFALSE;}
2498       if(fElectronLabelArray[i] == negLabel){
2499          return kFALSE;}
2500    }
2501
2502    return kTRUE;
2503 }
2504 ///________________________________________________________________________
2505 Bool_t AliConversionCuts::RejectToCloseV0s(AliAODConversionPhoton* photon, TList *photons, Int_t nV0){
2506
2507
2508    Double_t posX = photon->GetConversionX();
2509    Double_t posY = photon->GetConversionY();
2510    Double_t posZ = photon->GetConversionZ();
2511
2512    for(Int_t i = 0;i<photons->GetEntries();i++){
2513       if(nV0 == i) continue;
2514       AliAODConversionPhoton *photonComp = (AliAODConversionPhoton*) photons->At(i);
2515       Double_t posCompX = photonComp->GetConversionX();
2516       Double_t posCompY = photonComp->GetConversionY();
2517       Double_t posCompZ = photonComp->GetConversionZ();
2518       
2519       Double_t dist = pow((posX - posCompX),2)+pow((posY - posCompY),2)+pow((posZ - posCompZ),2);
2520
2521       if(dist < fminV0Dist*fminV0Dist){
2522          if(photon->GetChi2perNDF() < photonComp->GetChi2perNDF()) return kTRUE;
2523          else {
2524             return kFALSE;}
2525       }
2526       
2527    }
2528    return kTRUE;
2529 }
2530 ///________________________________________________________________________
2531 void AliConversionCuts::GetNotRejectedParticles(Int_t rejection, TList *HeaderList, AliMCEvent *MCEvent){
2532    
2533    if(fNotRejectedStart){
2534       delete[] fNotRejectedStart;
2535       fNotRejectedStart = NULL;
2536    }
2537    if(fNotRejectedEnd){
2538       delete[] fNotRejectedEnd;
2539       fNotRejectedEnd = NULL;
2540    }
2541    if(fGeneratorNames){
2542       delete[] fGeneratorNames;
2543       fGeneratorNames = NULL;
2544    }
2545    
2546    if(rejection == 0) return; // No Rejection
2547    AliGenCocktailEventHeader *cHeader = dynamic_cast<AliGenCocktailEventHeader*>(MCEvent->GenEventHeader());
2548    if(cHeader){
2549       TList *genHeaders = cHeader->GetHeaders();
2550       AliGenEventHeader* gh = 0;
2551       fnHeaders = 0;
2552       if(rejection == 1) fnHeaders = 1; // MinBiasHeader
2553       if(rejection == 2){ // TList of Headers Names
2554          for(Int_t i = 0; i<genHeaders->GetEntries();i++){
2555             gh = (AliGenEventHeader*)genHeaders->At(i);
2556             TString GeneratorName = gh->GetName();
2557             for(Int_t j = 0; j<HeaderList->GetEntries();j++){
2558                TString GeneratorInList = ((TObjString*)HeaderList->At(j))->GetString();
2559                if(GeneratorName.CompareTo(GeneratorInList) == 0){
2560                   fnHeaders++;
2561                   continue;
2562                }
2563             }
2564          }
2565       }
2566
2567       fNotRejectedStart = new Int_t[fnHeaders];
2568       fNotRejectedEnd = new Int_t[fnHeaders];
2569       fGeneratorNames = new TString[fnHeaders];
2570
2571       if(rejection == 1){
2572          fNotRejectedStart[0] = 0;
2573          fNotRejectedEnd[0] = ((AliGenEventHeader*)genHeaders->At(0))->NProduced();
2574          fGeneratorNames[0] = ((AliGenEventHeader*)genHeaders->At(0))->GetName();
2575
2576          return;
2577       }
2578       
2579       Int_t firstindex = 0;
2580       Int_t lastindex = -1;
2581       Int_t nummer = 0;
2582       for(Int_t i = 0; i<genHeaders->GetEntries();i++){
2583          gh = (AliGenEventHeader*)genHeaders->At(i);
2584          TString GeneratorName = gh->GetName();
2585          lastindex = lastindex + gh->NProduced();
2586          for(Int_t j = 0; j<HeaderList->GetEntries();j++){
2587             TString GeneratorInList = ((TObjString*)HeaderList->At(j))->GetString();
2588             if(GeneratorName.CompareTo(GeneratorInList) == 0){
2589                fNotRejectedStart[nummer] = firstindex;
2590                fNotRejectedEnd[nummer] = lastindex;
2591                fGeneratorNames[nummer] = GeneratorName;
2592                nummer++;
2593                continue;
2594             }
2595          }
2596          firstindex = firstindex + gh->NProduced();
2597       }
2598    }
2599    else{
2600       fNotRejectedStart = new Int_t[1];
2601       fNotRejectedEnd = new Int_t[1];
2602       fGeneratorNames = new TString[1];
2603       if(rejection == 1 || rejection == 0){  fnHeaders = 1; } else {fnHeaders = 0;}
2604       fNotRejectedStart[0] = 0;
2605       fNotRejectedEnd[0] = MCEvent->Stack()->GetNprimary();
2606       fGeneratorNames[0] = "NoGenerator";
2607       //cout << "rejection: " << rejection << " start: " << fNotRejectedStart[0] << " end: " << fNotRejectedEnd[0] << endl;
2608    }
2609    
2610 }
2611 //_________________________________________________________________________
2612 Bool_t AliConversionCuts::IsParticleFromBGEvent(Int_t index, AliStack *MCStack){
2613    
2614    if(index < 0) return kFALSE; // No Particle
2615    
2616    Bool_t accepted = kFALSE; 
2617    if( index >= MCStack->GetNprimary()){ // Secondary Particle
2618       if( ((TParticle*)MCStack->Particle(index))->GetMother(0) < 0) return kTRUE; // Secondary Particle without Mother??
2619       //cout<<"Secondary "<<index<<"   "<<((TParticle*)MCStack->Particle(index))->GetPdgCode()<<"  "<<((TParticle*)MCStack->Particle(index))->GetMother(0)<<endl;
2620       return IsParticleFromBGEvent(((TParticle*)MCStack->Particle(index))->GetMother(0),MCStack);
2621    }
2622    //   cout<<index<<"   "<<((TParticle*)MCStack->Particle(index))->GetPdgCode()<<"  "<<fNotRejectedStart[0]<<"  "<<fNotRejectedEnd[0]<<"  "<<MCStack->GetNprimary()<<endl;
2623    for(Int_t i = 0;i<fnHeaders;i++){
2624       if(index >= fNotRejectedStart[i] && index <= fNotRejectedEnd[i])
2625          accepted = kTRUE;
2626    }
2627    
2628    return accepted;
2629 }