PWG4/GammaConv/AliV0ReaderV1.cxx

   1 #include "AliV0ReaderV1.h"
   2 #include "AliKFParticle.h"
   3 #include "AliAODv0.h"
   4 #include "AliESDv0.h"
   5 #include "AliAODEvent.h"
   6 #include "AliESDEvent.h"
   7 #include "AliKFParticle.h"
   8 #include "AliKFConversionPhoton.h"
   9 #include "AliAODConversionPhoton.h"
  10 #include "AliConversionPhotonBase.h"
  11 #include "TVector.h"
  12 #include "AliKFVertex.h"
  13 #include "AliAODTrack.h"
  14 #include "AliESDtrack.h"
  15 #include "AliAnalysisManager.h"
  16 #include "AliInputEventHandler.h"
  17 #include "AliAODHandler.h"
  18 #include "AliPIDResponse.h"
  19 #include "TH1.h"
  20 #include "TH2.h"
  21 #include "TChain.h"
  22 #include "AliStack.h"
  23
  24 class iostream;
  25
  26
  27 using namespace std;
  28
  29 ClassImp(AliV0ReaderV1)
  30
  31 //________________________________________________________________________
  32 AliV0ReaderV1::AliV0ReaderV1(const char *name) : AliAnalysisTaskSE(name),
  33 fConversionGammas(NULL),
  34 fESDEvent(NULL),
  35 fAODEvent(NULL),
  36 fMCStack(NULL),
  37 fOutputList(NULL),
  38 fCurrentMotherKFCandidate(NULL),
  39 fCurrentPositiveKFParticle(NULL),
  40 fCurrentNegativeKFParticle(NULL),
  41 fCurrentTrackLabels(NULL),
  42 fCurrentV0Index(-1),
  43 fNCentralityBins(5),
  44     fCentralityBin(0),
  45     fBGHandler(NULL),
  46     fVertexZ(-999),
  47     fCentrality(-1),
  48     fEPAngle(-1),
  49 fMaxVertexZ(10),
  50 fMaxR(180),// 100 meter(outside of ALICE)
  51 fMinR(5),// 100 meter(outside of ALICE)
  52 fEtaCut(0.9),
  53 fEtaCutMin(-0.1),
  54 fPtCut(0.),
  55 fSinglePtCut(0.),
  56 fMaxZ(240),
  57 fMinClsTPC(0.),
  58 fMinClsTPCToF(0.),
  59 fLineCutZRSlope(0.),
  60 fLineCutZValue(7),
  61 fLineCutZRSlopeMin(0.),
  62 fLineCutZValueMin(-2),
  63 fChi2CutConversion(30),
  64 fPIDProbabilityCutNegativeParticle(0),
  65 fPIDProbabilityCutPositiveParticle(0),
  66 fDodEdxSigmaCut(kTRUE),
  67 fDoTOFsigmaCut(kFALSE), // RRnewTOF
  68 fPIDTRDEfficiency(0.95),
  69 fDoTRDPID(kFALSE),
  70 fPIDnSigmaAboveElectronLine(5),
  71 fPIDnSigmaBelowElectronLine(-3),
  72 fTofPIDnSigmaAboveElectronLine(100), // RRnewTOF
  73 fTofPIDnSigmaBelowElectronLine(-100), // RRnewTOF
  74 fPIDnSigmaAbovePionLine(0),
  75 fPIDnSigmaAbovePionLineHighPt(-100),
  76 fPIDMinPnSigmaAbovePionLine(1),
  77 fPIDMaxPnSigmaAbovePionLine(3),
  78 fDoKaonRejectionLowP(kTRUE),
  79     fDoProtonRejectionLowP(kTRUE),
  80     fDoPionRejectionLowP(kTRUE),
  81     fPIDnSigmaAtLowPAroundKaonLine(0),
  82     fPIDnSigmaAtLowPAroundProtonLine(0),
  83     fPIDnSigmaAtLowPAroundPionLine(0),
  84     fPIDMinPKaonRejectionLowP(1.5),
  85     fPIDMinPProtonRejectionLowP(2),
  86     fPIDMinPPionRejectionLowP(0.5),
  87     fDoQtGammaSelection(kTRUE),
  88     fDoHighPtQtGammaSelection(kTRUE), // RRnew
  89     fQtMax(0.05),
  90     fHighPtQtMax(0.06), // RRnew
  91     fPtBorderForQt(2.5), // RRnew
  92         fNSigmaMass(0.),
  93         fUseImprovedVertex(kTRUE),
  94         fUseOwnXYZCalculation(kTRUE),
  95         fUseConstructGamma(kFALSE),
  96         fUseEtaMinCut(kFALSE),
  97         fUseOnFlyV0Finder(kTRUE),
  98         fDoPhotonAsymmetryCut(kTRUE),
  99         fMinPPhotonAsymmetryCut(100.),
 100 fMinPhotonAsymmetry(0.),
 101     kUseAODConversionPhoton(kFALSE),
 102     fIsHeavyIon(kTRUE),
 103     fCreateAOD(kFALSE),
 104     fDeltaAODFilename("AliAODGammaConversion.root")
 105 {
 106
 107     fLineCutZRSlope=tan(2*atan(exp(-fEtaCut)));
 108         // Input slot #0 works with a TChain
 109   DefineInput(0, TChain::Class());
 110   // Output slot #0 id reserved by the base class for AOD
 111   // Output slot #1 writes into a TH1 container
 112   DefineOutput(1, TList::Class());
 113
 114   fCurrentTrackLabels=new Int_t[2];
 115 }
 116
 117
 118 //________________________________________________________________________
 119 AliV0ReaderV1::~AliV0ReaderV1()
 120 {
 121     if(fConversionGammas){
 122         fConversionGammas->Delete();// Clear Objects
 123         delete fConversionGammas;
 124         fConversionGammas=0x0;
 125     }
 126
 127     if(fCurrentTrackLabels){
 128         delete[] fCurrentTrackLabels;
 129         fCurrentTrackLabels=NULL;}
 130
 131 }
 132
 133 //________________________________________________________________________
 134 void AliV0ReaderV1::UserCreateOutputObjects()
 135 {
 136
 137 if(fOutputList != NULL){
 138     delete fOutputList;
 139     fOutputList = NULL;
 140   }
 141   if(fOutputList == NULL){
 142     fOutputList = new TList();
 143     fOutputList->SetOwner(kTRUE);
 144   }
 145
 146   TList *fCutList=new TList();
 147   fCutList->SetName("GammaReconstruction");
 148   fCutList->SetOwner(kTRUE);
 149   fOutputList->Add(fCutList);
 150
 151 //GammaMass-plots
 152 Int_t kGCnXBinsGammaMass                = 4000;
 153 Double_t kGCfirstXBinGammaMass= 0.;
 154 Double_t kGClastXBinGammaMass = 1.;
 155
 156   Int_t kGCnYBinsSpectra = 250;
 157   Double_t kGCfirstYBinSpectra = 0.;
 158   Double_t kGClastYBinSpectra = 25.;
 159
 160 // Process Gammas Histograms
 161
 162 hV0CurrentFinder=new TH1F("ESD_V0sCurrentFinder_InvMass","V0sCurrentFinder",kGCnXBinsGammaMass,kGCfirstXBinGammaMass,kGClastXBinGammaMass);
 163 fCutList->Add(hV0CurrentFinder);
 164 hV0AllArmenteros=new TH2F("ESD_V0sCurrentFinder_Armenteros","Armenteros Alpha Qt",200,-1,1,250,0,0.25);
 165 fCutList->Add(hV0AllArmenteros);
 166 hV0Good=new TH1F("ESD_GoodV0s_InvMass","GoodV0s",kGCnXBinsGammaMass,kGCfirstXBinGammaMass,kGClastXBinGammaMass);
 167 fCutList->Add(hV0Good);
 168 hV0GoodArmenteros=new TH2F("ESD_GoodV0s_Armenteros","Armenteros Alpha Qt",200,-1,1,250,0,0.25);
 169 fCutList->Add(hV0GoodArmenteros);
 170
 171
 172 // Track Cuts
 173 hV0CutLikeSign=new TH1F("ESD_CutLikeSign_InvMass","LikeSign",kGCnXBinsGammaMass,kGCfirstXBinGammaMass,kGClastXBinGammaMass);
 174 fCutList->Add(hV0CutLikeSign);
 175 hV0CutRefit=new TH1F("ESD_CutRefit_InvMass","No TPC refit",kGCnXBinsGammaMass,kGCfirstXBinGammaMass,kGClastXBinGammaMass);
 176 fCutList->Add(hV0CutRefit);
 177 hV0CutKinks=new TH1F("ESD_CutKink_InvMass","Kinks",kGCnXBinsGammaMass,kGCfirstXBinGammaMass,kGClastXBinGammaMass);
 178 fCutList->Add(hV0CutKinks);
 179 hV0CutMinNclsTPCToF=new TH1F("ESD_CutMinNClsTPCToF_InvMass","Min Ncls TPC ToF",kGCnXBinsGammaMass,kGCfirstXBinGammaMass,kGClastXBinGammaMass);
 180 fCutList->Add(hV0CutMinNclsTPCToF);
 181
 182 // Event Cuts
 183 /*hV0CutNContributors=new TH1F("ESD_CutNContributors_InvMass","NContributors<=0",kGCnXBinsGammaMass,kGCfirstXBinGammaMass,kGClastXBinGammaMass);
 184 fCutList->Add(hV0CutNContributors);
 185 hV0CutVertexZ=new TH1F("ESD_CutVertexZ_InvMass","VertexZ",kGCnXBinsGammaMass,kGCfirstXBinGammaMass,kGClastXBinGammaMass);
 186 fCutList->Add(hV0CutVertexZ);
 187 */
 188
 189
 190 // dEdx Cuts
 191
 192 hV0CutdEdxElectron=new TH1F("ESD_CutdEdxSigmaElectronLine_InvMass" ,"dedx ElectronLine" , kGCnXBinsGammaMass, kGCfirstXBinGammaMass, kGClastXBinGammaMass);
 193 fCutList->Add(hV0CutdEdxElectron);
 194 hV0CutdEdxPion=new TH1F("ESD_CutdEdxSigmaPionLine_InvMass" ,"dedx PionLine" , kGCnXBinsGammaMass, kGCfirstXBinGammaMass, kGClastXBinGammaMass);
 195 fCutList->Add(hV0CutdEdxPion);
 196 hV0CutdEdxKaonLowP=new TH1F("ESD_CutKaonRejectionLowP_InvMass" ,"dedx KaonRejection LowP" , kGCnXBinsGammaMass, kGCfirstXBinGammaMass, kGClastXBinGammaMass);
 197 fCutList->Add(hV0CutdEdxKaonLowP);
 198 hV0CutdEdxProtonLowP=new TH1F("ESD_CutProtonRejectionLowP_InvMass" ,"dedx ProtonRejection LowP" , kGCnXBinsGammaMass, kGCfirstXBinGammaMass, kGClastXBinGammaMass);
 199 fCutList->Add(hV0CutdEdxProtonLowP);
 200 hV0CutdEdxPionLowP=new TH1F("ESD_CutPionRejectionLowP_InvMass" ,"dedx PionRejection LowP" , kGCnXBinsGammaMass, kGCfirstXBinGammaMass, kGClastXBinGammaMass);
 201 fCutList->Add(hV0CutdEdxPionLowP);
 202 hV0CutdEdxTOFElectron=new TH1F("ESD_CutTOFsigmaElec_InvMass", "ESD_CutTOFsigmaElec_InvMass",kGCnXBinsGammaMass, kGCfirstXBinGammaMass, kGClastXBinGammaMass);
 203 fCutList->Add(hV0CutdEdxTOFElectron);
 204 hV0CutdEdxTRD=new TH1F("ESD_CutTRD_InvMass", "ESD_CutTRD_InvMass",kGCnXBinsGammaMass, kGCfirstXBinGammaMass, kGClastXBinGammaMass);
 205 fCutList->Add(hV0CutdEdxTRD);
 206 hGammadEdxbefore=new TH2F("Gamma_dEdx_before","dEdx Gamma before" ,kGCnYBinsSpectra, kGCfirstYBinSpectra, kGClastYBinSpectra,400, 0,200);
 207 fCutList->Add(hGammadEdxbefore);
 208 hGammadEdxafter=new TH2F("Gamma_dEdx_after","dEdx Gamma after" ,kGCnYBinsSpectra, kGCfirstYBinSpectra, kGClastYBinSpectra,400, 0,200);
 209 fCutList->Add(hGammadEdxafter);
 210
 211
 212 // Armenteros
 213
 214 hV0CutQt=new TH1F("ESD_CutQt_InvMass","ESD_CutQt_InvMass",kGCnXBinsGammaMass, kGCfirstXBinGammaMass, kGClastXBinGammaMass);
 215 fCutList->Add(hV0CutQt);
 216
 217 // Kinematic Cuts
 218
 219 hV0CutR=new TH1F("ESD_CutR_InvMass" ,"Above RMax" , kGCnXBinsGammaMass, kGCfirstXBinGammaMass, kGClastXBinGammaMass);
 220 fCutList->Add(hV0CutR);
 221 hV0CutMinR=new TH1F("ESD_CutMinR_InvMass" ,"Above RMax" , kGCnXBinsGammaMass, kGCfirstXBinGammaMass, kGClastXBinGammaMass);
 222 fCutList->Add(hV0CutMinR);
 223 hV0CutLine=new TH1F("ESD_CutLine_InvMass" ,"Out of reconstruction area" , kGCnXBinsGammaMass, kGCfirstXBinGammaMass, kGClastXBinGammaMass);
 224 fCutList->Add(hV0CutLine);
 225 hV0CutZ=new TH1F("ESD_CutZ_InvMass" ,"Out of reconstruction area" , kGCnXBinsGammaMass, kGCfirstXBinGammaMass, kGClastXBinGammaMass);
 226 fCutList->Add(hV0CutZ);
 227 hV0CutEta=new TH1F("ESD_CutEta_InvMass" ,"Above #eta max" , kGCnXBinsGammaMass, kGCfirstXBinGammaMass, kGClastXBinGammaMass);
 228 fCutList->Add(hV0CutEta);
 229
 230
 231 hV0CutSinglePt=new TH1F("ESD_CutSinglePt_InvMass" ,"Below p_{t} min" , kGCnXBinsGammaMass, kGCfirstXBinGammaMass, kGClastXBinGammaMass);
 232 fCutList->Add(hV0CutSinglePt);
 233 hV0CutNDF=new TH1F("ESD_CutNDF_InvMass" ,"#chi^{2} > Max" , kGCnXBinsGammaMass, kGCfirstXBinGammaMass, kGClastXBinGammaMass);
 234 fCutList->Add(hV0CutNDF);
 235 hV0CutChi2=new TH1F("ESD_CutChi2_InvMass" ,"#chi^{2} > Max" , kGCnXBinsGammaMass, kGCfirstXBinGammaMass, kGClastXBinGammaMass);
 236 fCutList->Add(hV0CutChi2);
 237 hV0CutPt= new TH1F("ESD_CutPt_InvMass" ,"Below p_{t} min" , kGCnXBinsGammaMass, kGCfirstXBinGammaMass, kGClastXBinGammaMass);
 238 fCutList->Add(hV0CutPt);
 239
 240 // Asymmetry Cut
 241
 242 hV0CutAsymmetry=new TH1F("ESD_CutPhotonAsymmetry_InvMass" ,"Out of reconstruction area" , kGCnXBinsGammaMass, kGCfirstXBinGammaMass, kGClastXBinGammaMass);
 243 fCutList->Add(hV0CutAsymmetry);
 244
 245 // PID Prob
 246
 247 hV0CutPIDProb=new TH1F("ESD_CutPIDProb_InvMass" ,"wrong TPC PID" , kGCnXBinsGammaMass, kGCfirstXBinGammaMass, kGClastXBinGammaMass);
 248 fCutList->Add(hV0CutPIDProb);
 249
 250 // Event Info
 251
 252    // Other
 253   TList *fOtherList=new TList();
 254   fOtherList->SetName("EventInfo");
 255   fOtherList->SetOwner(kTRUE);
 256   fOutputList->Add(fOtherList);
 257
 258
 259   hV0EventCuts=new TH1F("ESD_EventCuts","Event Cuts",10,-0.5,9.5);
 260   fOtherList->Add(hV0EventCuts);
 261   hNEvents=new TH1F("NEvents_vs_Centrality","NEvents vs Centrality",fNCentralityBins,-0.5,fNCentralityBins-0.5);
 262   fOtherList->Add(hNEvents);
 263   hCentrality=new TH1F("Centrality","Centrality",100,0,100);
 264   fOtherList->Add(hCentrality);
 265   hVertexZ=new TH1F("VertexZ","VertexZ",1000,-50,50);
 266   fOtherList->Add(hVertexZ);
 267
 268 // QA
 269
 270   TList *fGammaList=new TList();
 271   fGammaList->SetName("GammaInfo");
 272   fGammaList->SetOwner(kTRUE);
 273   fOutputList->Add(fGammaList);
 274
 275   hGammaPt=new TH1F*[fNCentralityBins];
 276
 277   for(int i=0;i<fNCentralityBins;i++){
 278       hGammaPt[i]=new TH1F(Form("GammaSpectrum_RECO_Pt_%d",i),"Reco Gamma Pt",kGCnYBinsSpectra, kGCfirstYBinSpectra, kGClastYBinSpectra);
 279       fGammaList->Add(hGammaPt[i]);
 280   }
 281
 282   hGammaPhi=new TH1F("Phi_Gamma","Phi Gamma" ,36, 0, 2*TMath::Pi());
 283   fGammaList->Add(hGammaPhi);
 284   hGammaConversionMapXY=new TH2F("Gamma_ConversionMap_XY","Conversion Point xy",400,-200,200,400,-200,200);
 285   fGammaList->Add(hGammaConversionMapXY);
 286   hGammaConversionMapZR=new TH2F("Gamma_ConversionMap_ZR","Conversion Point zr",500,-250,250,180,0,180);
 287   fGammaList->Add(hGammaConversionMapZR);
 288
 289
 290   // FILL MC PART only if MC is available
 291   if(AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler()){
 292
 293     hMCPtTRUE=new TH1F*[fNCentralityBins];
 294     hMCPtRECOTRUE=new TH1F*[fNCentralityBins];
 295
 296     for(int i=0;i<fNCentralityBins;i++){
 297         hMCPtTRUE[i]=new TH1F(Form("GammaSpectrum_MC_Pt_%d",i),"TRUE Gamma Pt",kGCnYBinsSpectra, kGCfirstYBinSpectra, kGClastYBinSpectra);
 298         fGammaList->Add(hMCPtTRUE[i]);
 299         hMCPtRECOTRUE[i]=new TH1F(Form("GammaSpectrum_RECOTRUE_Pt_%d",i),"True Reco Gamma Pt",kGCnYBinsSpectra, kGCfirstYBinSpectra, kGClastYBinSpectra);
 300         fGammaList->Add(hMCPtRECOTRUE[i]);
 301     }
 302
 303     // Call Sumw2 Option
 304     for (Int_t i=0; i<fGammaList->GetEntries(); i++) {
 305         TH1 *h1 = dynamic_cast<TH1*>(fGammaList->At(i));
 306         if (h1){h1->Sumw2();}
 307     }
 308
 309     hMCPtResolution=new TH2F("Resolution_Gamma_dPt_Pt","dPt vs Pt", kGCnYBinsSpectra, kGCfirstYBinSpectra, kGClastYBinSpectra,200,-10,10);
 310     fGammaList->Add(hMCPtResolution);
 311     hMCPtResolutionPhi=new TH2F("Resolution_Gamma_dPt_Phi","dPt vs Phi",180,0,2*TMath::Pi(),200,-10,10);
 312     fGammaList->Add(hMCPtResolutionPhi);
 313     hMCRResolutionvsR=new TH2F("Resolution_dRAbs_VS_R","dR vs R", 720,0,360,100,-5,5);
 314     fGammaList->Add(hMCRResolutionvsR);
 315     hMCZResolutionvsZ=new TH2F("Resolution_dZAbs_VS_Z","dZ vs Z", 200,-50,50,100,-5,5);
 316     fGammaList->Add(hMCZResolutionvsZ);
 317
 318 }
 319 // Gamma Output
 320
 321 if(fCreateAOD){kUseAODConversionPhoton=kTRUE;}
 322
 323 if(fConversionGammas == NULL){
 324     if(kUseAODConversionPhoton){
 325         fConversionGammas = new TClonesArray("AliAODConversionPhoton",100);}
 326     else{
 327         fConversionGammas = new TClonesArray("AliKFConversionPhoton",100);}
 328 }
 329 fConversionGammas->Delete();//Reset the TClonesArray
 330
 331 // Create AODs
 332
 333 if(fCreateAOD){
 334     fConversionGammas->SetName(Form("GammaConv_gamma"));
 335
 336     AddAODBranch("TClonesArray", &fConversionGammas, fDeltaAODFilename.Data());
 337     AliAnalysisManager::GetAnalysisManager()->RegisterExtraFile(fDeltaAODFilename.Data());
 338 }
 339
 340   PostData(1, fOutputList);
 341
 342 }
 343 //________________________________________________________________________
 344 void AliV0ReaderV1::UserExec(Option_t *){
 345
 346     fAODEvent = dynamic_cast<AliAODEvent*>(fInputEvent);
 347     fESDEvent = dynamic_cast<AliESDEvent*>(fInputEvent);
 348     if(!fAODEvent&&!fESDEvent) {
 349         AliError("No Input event");
 350         return;
 351     }
 352
 353     fMCStack=NULL;
 354     if(fMCEvent){
 355         fMCStack = fMCEvent->Stack();}
 356
 357
 358     if(fESDEvent)AliKFParticle::SetField(fESDEvent->GetMagneticField());
 359     if(fAODEvent)AliKFParticle::SetField(fAODEvent->GetMagneticField());
 360
 361     fConversionGammas->Delete();//Reset the TClonesArray
 362
 363     // Event Cuts
 364
 365     EventCuts();
 366
 367     if(EventIsSelected()){
 368
 369         // Process V0s
 370         for(fCurrentV0Index=0;fCurrentV0Index<fInputEvent->GetNumberOfV0s();fCurrentV0Index++){
 371             if(CheckV0Status()){
 372
 373                 ProcessV0();
 374             }
 375         }
 376         ProcessMCGammasForEfficiency();
 377
 378         // Set AOD Output
 379
 380         ///Make sure delta aod is filled if standard aod is filled (for synchronization when reading aod with standard aod)
 381         if(fCreateAOD) {
 382             AliAODHandler * aodhandler = dynamic_cast<AliAODHandler*>(AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler());
 383             if (aodhandler && aodhandler->GetFillAOD()) {
 384                 AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler()->SetFillExtension(kTRUE);
 385             }}
 386     }
 387
 388   PostData(1, fOutputList);
 389 }
 390  ///________________________________________________________________________
 391 Bool_t AliV0ReaderV1::IsGammaCandidate(AliConversionPhotonBase *fPhotonCandidate)
 392 {
 393
 394     // Fill Histos before Cuts
 395     hV0CurrentFinder->Fill(fPhotonCandidate->GetPhotonMass());
 396     hV0AllArmenteros->Fill(fPhotonCandidate->GetArmenterosAlpha(),fPhotonCandidate->GetArmenterosQt());
 397
 398     Bool_t passcuts=kTRUE;
 399
 400     // Gamma selection based on QT from Armenteros
 401     if(fDoQtGammaSelection == kTRUE){
 402         if(!ArmenterosQtCut(fPhotonCandidate))return kFALSE;//passcuts=kFALSE;
 403     }
 404
 405     // Chi Cut
 406
 407     if(fPhotonCandidate->GetChi2perNDF() > fChi2CutConversion || fPhotonCandidate->GetChi2perNDF() <=0){
 408         hV0CutChi2->Fill(fPhotonCandidate->GetPhotonMass());
 409         return kFALSE;
 410     }
 411
 412     // Reconstruction Acceptance Cuts
 413     if(!AcceptanceCuts(fPhotonCandidate))return kFALSE;//passcuts=kFALSE;
 414
 415
 416     // Track Cuts
 417     if(!TrackCuts(fPhotonCandidate))return kFALSE;//passcuts=kFALSE;
 418
 419
 420     // PID Cuts
 421     if(!dEdxCuts(fPhotonCandidate))return kFALSE;//passcuts=kFALSE;
 422
 423     // Asymmetry Cut
 424     if(fDoPhotonAsymmetryCut == kTRUE){
 425         if(!AsymmetryCut(fPhotonCandidate))return kFALSE;//passcuts=kFALSE;
 426     }
 427
 428     //Check the pid probability
 429
 430    if(!PIDProbabilityCut(fPhotonCandidate))return kFALSE;//passcuts=kFALSE;
 431
 432 return passcuts;
 433
 434 }
 435
 436 ///________________________________________________________________________
 437 const AliExternalTrackParam *AliV0ReaderV1::GetExternalTrackParam(Int_t charge){
 438
 439
 440     if(!(charge==1||charge==-1)){AliError("Charge not defined");return 0x0;}
 441
 442     Int_t label;
 443     if(charge>0)label=0;
 444     else label=1;
 445     // Check for sign flip
 446
 447     if(fESDEvent){
 448         AliESDv0 *fCurrentV0=dynamic_cast<AliESDv0*>(fESDEvent->GetV0(fCurrentV0Index));
 449         if(fCurrentV0){
 450             if(!fCurrentV0->GetParamN()||!fCurrentV0->GetParamP())return 0x0;
 451             if(!fESDEvent->GetTrack(fCurrentV0->GetNindex())||!fESDEvent->GetTrack(fCurrentV0->GetPindex()))return 0x0;
 452             if((fESDEvent->GetTrack(fCurrentV0->GetPindex()))->Charge()==charge){
 453                 fCurrentTrackLabels[label]=fCurrentV0->GetPindex();
 454                 return fCurrentV0->GetParamP();}
 455             if((fESDEvent->GetTrack(fCurrentV0->GetNindex()))->Charge()==charge){
 456                 fCurrentTrackLabels[label]=fCurrentV0->GetNindex();
 457                 return fCurrentV0->GetParamN();}
 458         }
 459
 460     }
 461
 462     if(fAODEvent){
 463 /*      AliAODv0 *fCurrentV0=dynamic_cast<AliAODv0*>(fAODEvent->GetV0(fCurrentV0Index));
 464         if(fCurrentV0){
 465             if(!fCurrentV0->GetParamN()||!fCurrentV0->GetParamP())return 0x0;
 466             if(!fAODEvent->GetTrack(fCurrentV0->GetNegID())||!fAODEvent->GetTrack(fCurrentV0->GetPosID()))return 0x0;
 467             if((fAODEvent->GetTrack(fCurrentV0->GetPosID()))->Charge()==1){
 468                 fCurrentTrackLabels[label]=fCurrentV0->GetPosID();
 469                 return fCurrentV0->GetParamP();}
 470             if((fAODEvent->GetTrack(fCurrentV0->GetNegID()))->Charge()==1){
 471                 fCurrentTrackLabels[label]=fCurrentV0->GetNegID();
 472                 return fCurrentV0->GetParamN();}
 473         }
 474   */  }
 475     return 0x0;
 476 }
 477
 478 ///________________________________________________________________________
 479 Bool_t AliV0ReaderV1::CheckV0Status()
 480 {
 481     if(fESDEvent){
 482         AliESDv0 *fCurrentV0=(AliESDv0*)(fESDEvent->GetV0(fCurrentV0Index));
 483         if(!fCurrentV0){
 484             printf("Requested V0 does not exist");
 485             return kFALSE;}
 486         //checks if on the fly mode is set
 487         if(fCurrentV0->GetOnFlyStatus()==fUseOnFlyV0Finder)return kTRUE;
 488     }
 489
 490     if(fAODEvent){
 491              AliAODv0 *fCurrentV0=dynamic_cast<AliAODv0*>(fAODEvent->GetV0(fCurrentV0Index));
 492              if(!fCurrentV0){
 493                  AliWarning("Requested V0 does not exist");
 494                  return kFALSE;}
 495
 496              //checks if on the fly mode is set
 497              if(fCurrentV0->GetOnFlyStatus()==fUseOnFlyV0Finder)return kTRUE;
 498     }
 499     return kFALSE;
 500 }
 501
 502
 503
 504
 505 ///________________________________________________________________________
 506 void AliV0ReaderV1::ProcessV0(){
 507
 508     // Reset TrackLabels
 509     fCurrentTrackLabels[0]=-1;
 510     fCurrentTrackLabels[1]=-1;
 511
 512
 513    //  cout<<"V0ReaderV1 ProcessV0 "<<fCurrentV0Index<<endl;
 514
 515     // Get Daughter KF Particles
 516
 517     const AliExternalTrackParam *fCurrentExternalTrackParamPositive=GetExternalTrackParamP();
 518     const AliExternalTrackParam *fCurrentExternalTrackParamNegative=GetExternalTrackParamN();
 519
 520     if(fCurrentExternalTrackParamPositive&&fCurrentExternalTrackParamNegative){
 521
 522         fCurrentNegativeKFParticle=new AliKFParticle(*(fCurrentExternalTrackParamNegative),11);
 523         fCurrentPositiveKFParticle=new AliKFParticle(*(fCurrentExternalTrackParamPositive),-11);
 524
 525     }
 526     //else{hV0CutLikeSign->Fill(fCurrentMotherKFCandidate->M());}// Like Sign error is already catched here
 527
 528
 529     // Reconstruct Gamma
 530
 531     if(fCurrentNegativeKFParticle&&fCurrentPositiveKFParticle){
 532
 533         if(fUseConstructGamma==kTRUE){
 534
 535             fCurrentMotherKFCandidate = new AliKFConversionPhoton();
 536                 fCurrentMotherKFCandidate->ConstructGamma(*fCurrentNegativeKFParticle,*fCurrentPositiveKFParticle);
 537         }else{
 538                 fCurrentMotherKFCandidate = new AliKFConversionPhoton(*fCurrentNegativeKFParticle,*fCurrentPositiveKFParticle);
 539                 fCurrentMotherKFCandidate->SetMassConstraint(0,fNSigmaMass);
 540         }
 541
 542         if(fCurrentNegativeKFParticle){delete fCurrentNegativeKFParticle;
 543             fCurrentNegativeKFParticle=0x0;}
 544         if(fCurrentPositiveKFParticle){ delete fCurrentPositiveKFParticle;
 545             fCurrentPositiveKFParticle=0x0;}
 546
 547
 548         // Update Vertex
 549         if(fUseImprovedVertex == kTRUE){
 550                 AliKFVertex primaryVertexImproved(*GetPrimaryVertex());
 551                 primaryVertexImproved+=*fCurrentMotherKFCandidate;
 552                 fCurrentMotherKFCandidate->SetProductionVertex(primaryVertexImproved);
 553         }
 554
 555         // Set Track Labels
 556
 557         fCurrentMotherKFCandidate->SetV0Index(fCurrentV0Index);
 558         fCurrentMotherKFCandidate->SetTrackLabels(fCurrentTrackLabels[0],fCurrentTrackLabels[1]);
 559
 560         //Set MC Label
 561
 562         if(fMCStack){
 563             Int_t labeln=TMath::Abs(GetTrack(fCurrentMotherKFCandidate->GetTrackLabelPositive())->GetLabel());
 564             Int_t labelp=TMath::Abs(GetTrack(fCurrentMotherKFCandidate->GetTrackLabelNegative())->GetLabel());
 565
 566                TParticle *fNegativeMCParticle = fMCStack->Particle(labeln);
 567                TParticle *fPositiveMCParticle = fMCStack->Particle(labelp);
 568
 569                if(fPositiveMCParticle&&fNegativeMCParticle){
 570                    fCurrentMotherKFCandidate->SetMCLabelPositive(labelp);
 571                    fCurrentMotherKFCandidate->SetMCLabelNegative(labeln);
 572                }
 573         }
 574
 575
 576
 577         //Add PID information with ESD tender (AOD implementation is not complete)
 578
 579
 580         AliAnalysisManager *man=AliAnalysisManager::GetAnalysisManager();
 581         AliInputEventHandler* inputHandler = (AliInputEventHandler*) (man->GetInputEventHandler());
 582         AliPIDResponse *fPIDResponse = (AliPIDResponse*)inputHandler->GetPIDResponse();
 583
 584         if(fESDEvent){
 585             Int_t labelp=((AliESDv0*)fESDEvent->GetV0(fCurrentMotherKFCandidate->GetV0Index()))->GetNindex();
 586             Int_t labeln=((AliESDv0*)fESDEvent->GetV0(fCurrentMotherKFCandidate->GetV0Index()))->GetNindex();
 587
 588         AliESDtrack *trackpos=fESDEvent->GetTrack(labelp);
 589         AliESDtrack *trackneg=fESDEvent->GetTrack(labeln);
 590
 591         if(trackpos&&trackneg){
 592
 593             Float_t fNSigmadEdxPositive[5];
 594             Float_t fNSigmadEdxNegative[5];
 595
 596             fNSigmadEdxPositive[0]=fPIDResponse->NumberOfSigmasTPC(trackpos,AliPID::kElectron);
 597             fNSigmadEdxPositive[1]=fPIDResponse->NumberOfSigmasTPC(trackpos,AliPID::kMuon);
 598             fNSigmadEdxPositive[2]=fPIDResponse->NumberOfSigmasTPC(trackpos,AliPID::kPion);
 599             fNSigmadEdxPositive[3]=fPIDResponse->NumberOfSigmasTPC(trackpos,AliPID::kKaon);
 600             fNSigmadEdxPositive[4]=fPIDResponse->NumberOfSigmasTPC(trackpos,AliPID::kProton);
 601
 602             fNSigmadEdxNegative[0]=fPIDResponse->NumberOfSigmasTPC(trackneg,AliPID::kElectron);
 603             fNSigmadEdxNegative[1]=fPIDResponse->NumberOfSigmasTPC(trackneg,AliPID::kMuon);
 604             fNSigmadEdxNegative[2]=fPIDResponse->NumberOfSigmasTPC(trackneg,AliPID::kPion);
 605             fNSigmadEdxNegative[3]=fPIDResponse->NumberOfSigmasTPC(trackneg,AliPID::kKaon);
 606             fNSigmadEdxNegative[4]=fPIDResponse->NumberOfSigmasTPC(trackneg,AliPID::kProton);
 607
 608         fCurrentMotherKFCandidate->SetNSigmadEdx(fNSigmadEdxPositive,fNSigmadEdxNegative);
 609         }
 610     }
 611
 612         // Calculate ConversionPoint
 613
 614         if(fUseOwnXYZCalculation){
 615
 616             Double_t convpos[3]={0,0,0};
 617             GetConversionPoint(fCurrentExternalTrackParamPositive,fCurrentExternalTrackParamNegative,convpos);
 618             fCurrentMotherKFCandidate->SetConversionPoint(convpos);
 619
 620         }
 621
 622        // if(kTRUE){
 623        if(IsGammaCandidate(fCurrentMotherKFCandidate)){
 624            // Fill Histos after Cuts
 625
 626            // Process MC
 627
 628            ProcessMC(fCurrentMotherKFCandidate);
 629
 630
 631             hV0Good->Fill(fCurrentMotherKFCandidate->M());
 632             hV0GoodArmenteros->Fill(fCurrentMotherKFCandidate->GetArmenterosAlpha(),fCurrentMotherKFCandidate->GetArmenterosQt());
 633
 634
 635             // Set Mass Zero for Gammas
 636             SetGammaMassZero();
 637
 638             // Add Gamma to the TClonesArray
 639
 640             if(kUseAODConversionPhoton){
 641                 new((*fConversionGammas)[fConversionGammas->GetEntriesFast()]) AliAODConversionPhoton(fCurrentMotherKFCandidate);
 642             }
 643             else{
 644                 new((*fConversionGammas)[fConversionGammas->GetEntriesFast()]) AliKFConversionPhoton(*fCurrentMotherKFCandidate);
 645             }
 646             // Fill QA Histos
 647
 648             hGammaPhi->Fill(fCurrentMotherKFCandidate->Phi());
 649             hGammaPt[fCentralityBin]->Fill(fCurrentMotherKFCandidate->Pt());
 650             hGammaConversionMapXY->Fill(fCurrentMotherKFCandidate->GetConversionX(),fCurrentMotherKFCandidate->GetConversionY());
 651             hGammaConversionMapZR->Fill(fCurrentMotherKFCandidate->GetConversionZ(),fCurrentMotherKFCandidate->GetConversionRadius());
 652
 653        }
 654         delete fCurrentMotherKFCandidate;
 655         fCurrentMotherKFCandidate=NULL;
 656     }
 657
 658 }
 659
 660 ///________________________________________________________________________
 661 Bool_t AliV0ReaderV1::ArmenterosQtCut(AliConversionPhotonBase *fPhotonCandidate)
 662 {
 663     if(fDoHighPtQtGammaSelection){
 664         if(fPhotonCandidate->GetPhotonPt() < fPtBorderForQt){
 665             if(fPhotonCandidate->GetArmenterosQt()>fQtMax){
 666                 hV0CutQt->Fill(fPhotonCandidate->GetPhotonMass());
 667                 return kFALSE;
 668             }
 669         } else {
 670             if(fPhotonCandidate->GetArmenterosQt()>fHighPtQtMax){
 671                 hV0CutQt->Fill(fPhotonCandidate->GetPhotonMass());
 672                 return kFALSE;
 673             }
 674         }
 675     } else {
 676
 677         if(fPhotonCandidate->GetArmenterosQt()>fQtMax){
 678            hV0CutQt->Fill(fPhotonCandidate->GetPhotonMass());
 679            return kFALSE;
 680         }
 681     }
 682     return kTRUE;
 683 }
 684
 685 ///________________________________________________________________________
 686 Bool_t AliV0ReaderV1::AcceptanceCuts(AliConversionPhotonBase *fPhotonCandidate)
 687 {
 688     AliVTrack *fCurrentNegativeTrack=GetTrack(fPhotonCandidate->GetTrackLabelNegative());
 689     AliVTrack *fCurrentPositiveTrack=GetTrack(fPhotonCandidate->GetTrackLabelPositive());
 690
 691     if(fPhotonCandidate->GetConversionRadius()>fMaxR){ // cuts on distance from collision point
 692         hV0CutR->Fill(fPhotonCandidate->GetPhotonMass());
 693         return kFALSE;
 694              }
 695
 696     if(fPhotonCandidate->GetConversionRadius()<fMinR){ // cuts on distance from collision point
 697         hV0CutMinR->Fill(fPhotonCandidate->GetPhotonMass());
 698         return kFALSE;
 699     }
 700
 701    if(fPhotonCandidate->GetConversionRadius() <= ((TMath::Abs(fPhotonCandidate->GetConversionZ())*fLineCutZRSlope)-fLineCutZValue)){
 702         hV0CutLine->Fill(fPhotonCandidate->GetPhotonMass());
 703         return kFALSE;
 704    }
 705    /*else if (fUseEtaMinCut &&  fPhotonCandidate->GetConversionRadius() >= ((TMath::Abs(fPhotonCandidate->GetConversionZ())*fLineCutZRSlopeMin)-fLineCutZValueMin )){
 706         hV0CutLine->Fill(fPhotonCandidate->GetPhotonMass());
 707         return kFALSE;
 708     }*/
 709
 710
 711     if(TMath::Abs(fPhotonCandidate->GetConversionZ()) > fMaxZ ){ // cuts out regions where we do not reconstruct
 712         hV0CutZ->Fill(fPhotonCandidate->GetPhotonMass());
 713         return kFALSE;
 714     }
 715
 716     if(TMath::Abs(fPhotonCandidate->GetPhotonEta())> fEtaCut || TMath::Abs(fPhotonCandidate->GetPhotonEta())< fEtaCutMin){
 717         hV0CutEta->Fill(fPhotonCandidate->GetPhotonMass());
 718         return kFALSE;
 719     }
 720
 721     if(TMath::Abs(fCurrentNegativeTrack->Eta())> fEtaCut || TMath::Abs(fCurrentNegativeTrack->Eta())< fEtaCutMin){
 722         hV0CutEta->Fill(fPhotonCandidate->GetPhotonMass());
 723         return kFALSE;
 724     }
 725
 726     if(TMath::Abs(fCurrentPositiveTrack->Eta())> fEtaCut || TMath::Abs(fCurrentPositiveTrack->Eta())< fEtaCutMin){
 727         hV0CutEta->Fill(fPhotonCandidate->GetPhotonMass());
 728         return kFALSE;
 729     }
 730
 731     if( fCurrentNegativeTrack->Pt()< fSinglePtCut ||    fCurrentNegativeTrack->Pt()< fSinglePtCut){
 732         hV0CutSinglePt->Fill(fPhotonCandidate->GetPhotonMass());
 733         return kFALSE;
 734     }
 735
 736
 737     if(fPhotonCandidate->GetPhotonPt()<fPtCut){
 738         hV0CutPt->Fill(fPhotonCandidate->GetPhotonMass());
 739         return kFALSE;
 740     }
 741 return kTRUE;
 742 }
 743
 744
 745
 746 ///________________________________________________________________________
 747 Bool_t AliV0ReaderV1::TrackCuts(AliConversionPhotonBase *fPhotonCandidate){
 748
 749     Bool_t passtrackcuts=kTRUE;
 750
 751
 752
 753     if(fESDEvent){
 754
 755         AliESDtrack *fCurrentNegativeESDTrack=(AliESDtrack*)fESDEvent->GetTrack(fPhotonCandidate->GetTrackLabelNegative());
 756         AliESDtrack *fCurrentPositiveESDTrack=(AliESDtrack*)fESDEvent->GetTrack(fPhotonCandidate->GetTrackLabelPositive());
 757
 758           if(!fCurrentNegativeESDTrack||!fCurrentPositiveESDTrack)return kFALSE;
 759
 760         // avoid like sign
 761         if(fCurrentNegativeESDTrack->Charge() == fCurrentPositiveESDTrack->Charge()){
 762
 763             hV0CutLikeSign->Fill(fPhotonCandidate->GetPhotonMass());
 764             passtrackcuts=kFALSE;
 765         }
 766
 767         if( (!(fCurrentNegativeESDTrack->IsOn(AliESDtrack::kTPCrefit))||(!(fCurrentPositiveESDTrack->IsOn(AliESDtrack::kTPCrefit))))){
 768             hV0CutRefit->Fill(fPhotonCandidate->GetPhotonMass());
 769             passtrackcuts=kFALSE;
 770
 771         }
 772
 773         if( fCurrentNegativeESDTrack->GetKinkIndex(0) > 0 ||
 774            fCurrentPositiveESDTrack->GetKinkIndex(0) > 0) {
 775             passtrackcuts=kFALSE;
 776         }
 777
 778         if(fCurrentNegativeESDTrack->GetNcls(1) < fMinClsTPC || fCurrentPositiveESDTrack->GetNcls(1) < fMinClsTPC ){
 779             passtrackcuts=kFALSE;
 780             hV0CutKinks->Fill(fPhotonCandidate->GetPhotonMass());
 781         }
 782
 783              /* Double_t negclsToF = 0.;
 784                 if (!fUseCorrectedTPCClsInfo ){
 785                         if(fCurrentNegativeESDTrack->GetTPCNclsF()!=0   ){
 786                                 negclsToF = (Double_t)fCurrentNegativeESDTrack->GetNcls(1)/(Double_t)fCurrentNegativeESDTrack->GetTPCNclsF();
 787                         }
 788                 } else {
 789                         negclsToF = fCurrentNegativeESDTrack->GetTPCClusterInfo(2,0,GetFirstTPCRow(GetXYRadius()));
 790                 }
 791
 792                 Double_t posclsToF = 0.;
 793                 if (!fUseCorrectedTPCClsInfo ){
 794                         if(fCurrentTrack->GetTPCNclsF()!=0      ){
 795                                 posclsToF = (Double_t)fCurrentTrack->GetNcls(1)/(Double_t)fCurrentTrack->GetTPCNclsF();
 796                         }
 797                 }else{
 798                         posclsToF = fCurrentTrack->GetTPCClusterInfo(2,0,GetFirstTPCRow(GetXYRadius()));
 799                 }
 800
 801                 if( negclsToF < fMinClsTPCToF ||        posclsToF < fMinClsTPCToF ){
 802                     hV0CutMinNclsTPCToF->Fill(fPhotonCandidate->GetPhotonMass());
 803                     passtrackcuts=kFALSE; }
 804                */
 805
 806
 807
 808
 809     }
 810
 811     if(fAODEvent){
 812
 813         AliAODTrack *fCurrentNegativeESDTrack=(AliAODTrack*)fAODEvent->GetTrack(fPhotonCandidate->GetTrackLabelNegative());
 814         AliAODTrack *fCurrentPositiveESDTrack=(AliAODTrack*)fAODEvent->GetTrack(fPhotonCandidate->GetTrackLabelPositive());
 815
 816         if(!fCurrentNegativeESDTrack||!fCurrentPositiveESDTrack)return kFALSE;
 817
 818         // avoid like sign
 819         if(fCurrentNegativeESDTrack->Charge() == fCurrentPositiveESDTrack->Charge()){
 820
 821             hV0CutLikeSign->Fill(fPhotonCandidate->GetPhotonMass());
 822             passtrackcuts=kFALSE;
 823         }
 824
 825         if( !(fCurrentNegativeESDTrack->IsOn(AliESDtrack::kTPCrefit))){
 826             hV0CutRefit->Fill(fPhotonCandidate->GetPhotonMass());
 827             passtrackcuts=kFALSE;
 828         }
 829
 830         if( !(fCurrentPositiveESDTrack->IsOn(AliESDtrack::kTPCrefit))){
 831             hV0CutRefit->Fill(fPhotonCandidate->GetPhotonMass());
 832             passtrackcuts=kFALSE;
 833         }
 834
 835         // to be implemented
 836         /*
 837          if( fCurrentNegativeESDTrack->GetKinkIndex(0) > 0 ||
 838          fCurrentPositiveESDTrack->GetKinkIndex(0) > 0) {
 839          }*/
 840
 841         if(fCurrentNegativeESDTrack->GetNcls(1) < fMinClsTPC || fCurrentPositiveESDTrack->GetNcls(1) < fMinClsTPC ){
 842             passtrackcuts=kFALSE;}
 843
 844
 845     }
 846
 847     return passtrackcuts;
 848 }
 849
 850 ///________________________________________________________________________
 851 Bool_t AliV0ReaderV1::dEdxCuts(AliConversionPhotonBase *fPhotonCandidate){
 852
 853     AliAnalysisManager *man=AliAnalysisManager::GetAnalysisManager();
 854     AliInputEventHandler* inputHandler = (AliInputEventHandler*) (man->GetInputEventHandler());
 855     AliPIDResponse *fPIDResponse = (AliPIDResponse*)inputHandler->GetPIDResponse();
 856
 857     AliVTrack *fCurrentTrack=0x0;
 858
 859     for(Int_t ilabel=0;ilabel<2;ilabel++){
 860
 861         fCurrentTrack=GetTrack(fPhotonCandidate->GetTrackLabel(ilabel));
 862
 863         if(!fCurrentTrack)return kFALSE;
 864
 865         // Fill dEdx before cuts
 866
 867         hGammadEdxbefore->Fill(fCurrentTrack->P(),fCurrentTrack->GetTPCsignal());
 868
 869          if(fDodEdxSigmaCut == kTRUE){
 870              if( fPIDResponse->NumberOfSigmasTPC(fCurrentTrack,AliPID::kElectron)<fPIDnSigmaBelowElectronLine ||
 871                 fPIDResponse->NumberOfSigmasTPC(fCurrentTrack,AliPID::kElectron)>fPIDnSigmaAboveElectronLine){
 872
 873                 hV0CutdEdxElectron->Fill(fPhotonCandidate->GetPhotonMass());
 874                  return kFALSE;
 875              }
 876
 877              if( fCurrentTrack->P()>fPIDMinPnSigmaAbovePionLine && fCurrentTrack->P()<fPIDMaxPnSigmaAbovePionLine ){
 878                  if(fPIDResponse->NumberOfSigmasTPC(fCurrentTrack,AliPID::kElectron)>fPIDnSigmaBelowElectronLine &&
 879                     fPIDResponse->NumberOfSigmasTPC(fCurrentTrack,AliPID::kElectron)<fPIDnSigmaAboveElectronLine&&
 880                     fPIDResponse->NumberOfSigmasTPC(fCurrentTrack,AliPID::kPion)<fPIDnSigmaAbovePionLine){
 881
 882                     hV0CutdEdxPion->Fill(fPhotonCandidate->GetPhotonMass());
 883                      return kFALSE;
 884                  }
 885              }
 886
 887              // High Pt
 888              if( fCurrentTrack->P()>fPIDMaxPnSigmaAbovePionLine ){
 889                  if(fPIDResponse->NumberOfSigmasTPC(fCurrentTrack,AliPID::kElectron)>fPIDnSigmaBelowElectronLine &&
 890                     fPIDResponse->NumberOfSigmasTPC(fCurrentTrack,AliPID::kElectron)<fPIDnSigmaAboveElectronLine&&
 891                     fPIDResponse->NumberOfSigmasTPC(fCurrentTrack,AliPID::kPion)<fPIDnSigmaAbovePionLineHighPt){
 892
 893                     hV0CutdEdxPion->Fill(fPhotonCandidate->GetPhotonMass());
 894                      return kFALSE;
 895                  }
 896              }
 897          }
 898
 899          if(fDoKaonRejectionLowP == kTRUE){
 900              if(fCurrentTrack->P()<fPIDMinPKaonRejectionLowP ){
 901                  if( TMath::Abs(fPIDResponse->NumberOfSigmasTPC(fCurrentTrack,AliPID::kKaon))<fPIDnSigmaAtLowPAroundKaonLine){
 902                     hV0CutdEdxKaonLowP->Fill(fPhotonCandidate->GetPhotonMass());
 903                     return kFALSE;
 904                  }
 905              }
 906          }
 907
 908          if(fDoProtonRejectionLowP == kTRUE){
 909              if( fCurrentTrack->P()<fPIDMinPProtonRejectionLowP ){
 910                  if( TMath::Abs(fPIDResponse->NumberOfSigmasTPC(fCurrentTrack,AliPID::kProton))<fPIDnSigmaAtLowPAroundProtonLine){
 911                     hV0CutdEdxProtonLowP->Fill(fPhotonCandidate->GetPhotonMass());
 912                      return kFALSE;
 913                  }
 914              }
 915          }
 916
 917          if(fDoPionRejectionLowP == kTRUE){
 918              if( fCurrentTrack->P()<fPIDMinPPionRejectionLowP ){
 919                  if( TMath::Abs(fPIDResponse->NumberOfSigmasTPC(fCurrentTrack,AliPID::kPion))<fPIDnSigmaAtLowPAroundPionLine){
 920                     hV0CutdEdxPionLowP->Fill(fPhotonCandidate->GetPhotonMass());
 921                  return kFALSE;
 922                  }
 923              }
 924          }
 925
 926
 927          if( fDoTOFsigmaCut == kTRUE ){ // RRnewTOF start /////////////////////////////////////////////////////////////////////////////
 928
 929             if((fPIDResponse->NumberOfSigmasTOF(fCurrentTrack,AliPID::kElectron)>fTofPIDnSigmaAboveElectronLine) || (fPIDResponse->NumberOfSigmasTOF(fCurrentTrack,AliPID::kElectron)<fTofPIDnSigmaBelowElectronLine)){
 930               hV0CutdEdxTOFElectron->Fill(fPhotonCandidate->GetPhotonMass());
 931               return kFALSE;
 932             }
 933          } /////////////////////////////// RRnewTOF end ///////////////////////////////////////////////////////////////////////////////
 934
 935          // Apply TRD PID
 936          if(fDoTRDPID){
 937              if(!fPIDResponse->IdentifiedAsElectronTRD(fCurrentTrack,fPIDTRDEfficiency)){
 938                  hV0CutdEdxTRD->Fill(fPhotonCandidate->GetPhotonMass());
 939                  return kFALSE;
 940              }
 941          }
 942
 943          // Fill dEdx Histogram after Cuts
 944
 945          hGammadEdxafter->Fill(fCurrentTrack->P(),fCurrentTrack->GetTPCsignal());
 946
 947     }
 948
 949     return kTRUE;
 950 }
 951
 952 ///________________________________________________________________________
 953 Bool_t AliV0ReaderV1::AsymmetryCut(AliConversionPhotonBase *fPhotonCandidate)
 954 {
 955
 956     for(Int_t ii=0;ii<2;ii++){
 957         AliVTrack *fCurrentTrack=GetTrack(fPhotonCandidate->GetTrackLabel(ii));
 958
 959         if( fCurrentTrack->P()>fMinPPhotonAsymmetryCut ){
 960             Double_t trackNegAsy=0;
 961             if (fPhotonCandidate->GetPhotonP()!=0.){
 962                 trackNegAsy= fCurrentTrack->P()/fPhotonCandidate->GetPhotonP();
 963             }
 964             if( trackNegAsy<fMinPhotonAsymmetry ||trackNegAsy>(1.- fMinPhotonAsymmetry)){
 965                 hV0CutAsymmetry->Fill(fPhotonCandidate->GetPhotonMass());
 966                 return kFALSE;
 967             }
 968         }
 969     }
 970
 971     return kTRUE;
 972 }
 973
 974
 975 ///________________________________________________________________________
 976
 977 Int_t AliV0ReaderV1::GetNumberOfContributorsVtx(){
 978         // returns number of contributors to the vertex
 979     if(fESDEvent){
 980         if(fESDEvent->GetPrimaryVertexTracks()->GetNContributors()>0) {
 981             return fESDEvent->GetPrimaryVertexTracks()->GetNContributors();
 982         }
 983
 984         if(fESDEvent->GetPrimaryVertexTracks()->GetNContributors()<1) {
 985             //          return 0;
 986             //-AM test pi0s without SPD only vertex
 987             if(fESDEvent->GetPrimaryVertexSPD()->GetNContributors()>0) {
 988                 return fESDEvent->GetPrimaryVertexSPD()->GetNContributors();
 989
 990             }
 991             if(fESDEvent->GetPrimaryVertexSPD()->GetNContributors()<1) {
 992              //   cout<<"number of contributors from bad vertex type::"<< fESDEvent->GetPrimaryVertex()->GetName() << endl;
 993                 return 0;
 994             }
 995         }
 996     }
 997     if(fAODEvent){
 998           if(fAODEvent->GetPrimaryVertex()->GetNContributors()>0) {
 999             return fESDEvent->GetPrimaryVertex()->GetNContributors();
1000         }
1001         if(fAODEvent->GetPrimaryVertex()->GetNContributors()<1) {
1002             if(fAODEvent->GetPrimaryVertexSPD()->GetNContributors()>0) {
1003                 return fAODEvent->GetPrimaryVertexSPD()->GetNContributors();
1004
1005             }
1006             if(fAODEvent->GetPrimaryVertexSPD()->GetNContributors()<1) {
1007                 AliWarning(Form("Number of contributors from bad vertex type:: %s",fAODEvent->GetPrimaryVertex()->GetName()));
1008                 return 0;
1009             }
1010         }
1011     }
1012
1013     return 0;
1014
1015 }
1016
1017 ///________________________________________________________________________
1018 Double_t AliV0ReaderV1::GetCentrality(){
1019     if(fAODEvent){
1020         if(fAODEvent->GetHeader()){return fAODEvent->GetHeader()->GetCentrality();}
1021     }
1022
1023     if(fESDEvent){
1024         return ((AliCentrality *)fESDEvent->GetCentrality())->GetCentralityPercentile("V0M");
1025     }
1026
1027     return -1;
1028 }
1029
1030 ///________________________________________________________________________
1031
1032 AliEventplane *AliV0ReaderV1::GetEventPlane(){
1033   if(fESDEvent){
1034            return fESDEvent->GetEventplane();
1035   }
1036   if(fAODEvent){
1037       return fAODEvent->GetEventplane();
1038   }
1039 return 0x0;
1040 }
1041
1042 //________________________________________________________________________
1043 Bool_t AliV0ReaderV1::SetEventPlane()
1044 {
1045     if(fMCStack||!fIsHeavyIon){
1046         // NO EP in MC mode and pp mode
1047         fEPAngle=0;
1048        return kTRUE;
1049     }
1050
1051     if(fIsHeavyIon){
1052         AliEventplane *fEP=GetEventPlane();
1053
1054         if(!fEP)return kFALSE;
1055
1056         fEPAngle=fEP->GetEventplane("Q");
1057
1058         return kTRUE;
1059     }
1060     return kFALSE;
1061 }
1062
1063 ///________________________________________________________________________
1064 Bool_t AliV0ReaderV1::EventCuts(){
1065
1066     fEventIsSelected=kTRUE;
1067
1068     // Z Vertex Position Cut
1069
1070     if(!VertexZCut()){
1071         hV0EventCuts->Fill(0);
1072         fEventIsSelected=kFALSE;;
1073     }
1074     // Number of Contributors Cut
1075
1076     if(GetNumberOfContributorsVtx()<=0) {
1077         hV0EventCuts->Fill(1);
1078         fEventIsSelected=kFALSE;;
1079     }
1080     // Centrality Selection
1081
1082     if(!CentralitySelection()){
1083         hV0EventCuts->Fill(2);
1084         fEventIsSelected=kFALSE;;
1085     }
1086
1087     // Event Plane
1088     if(!SetEventPlane()){hV0EventCuts->Fill(4);}
1089
1090     // Fill Event Histograms
1091
1092     if(fEventIsSelected){
1093           // Fill Event Histograms
1094           hV0EventCuts->Fill(9);
1095           hVertexZ->Fill(fVertexZ);
1096           hCentrality->Fill(fCentrality);
1097           hNEvents->Fill(fCentralityBin);
1098     }
1099     else{hV0EventCuts->Fill(8);}
1100
1101     return fEventIsSelected;
1102 }
1103
1104 ///________________________________________________________________________
1105 Bool_t AliV0ReaderV1::VertexZCut(){
1106
1107     fVertexZ=GetPrimaryVertex()->GetZ();
1108
1109     if(fBGHandler){
1110         if(fBGHandler->GetZBinIndex(fVertexZ)<0)return kFALSE;
1111     }
1112     else{
1113         if(fVertexZ>fMaxVertexZ)return kFALSE;
1114     }
1115     return kTRUE;
1116 }
1117
1118 ///________________________________________________________________________
1119 AliVTrack *AliV0ReaderV1::GetTrack(Int_t label){
1120     if(fESDEvent){
1121         return (AliESDtrack*)fESDEvent->GetTrack(label);
1122     }
1123     if(fAODEvent)return (AliAODTrack*)fAODEvent->GetTrack(label);
1124     return 0x0;
1125 }
1126
1127 ///________________________________________________________________________
1128 Bool_t AliV0ReaderV1::PIDProbabilityCut(AliConversionPhotonBase *fPhotonCandidate){
1129
1130     if(fESDEvent){
1131
1132         Bool_t iResult=kFALSE;
1133
1134         Double_t *posProbArray = new Double_t[AliPID::kSPECIES];
1135         Double_t *negProbArray = new Double_t[AliPID::kSPECIES];
1136
1137         AliESDtrack* negTrack   = (AliESDtrack*)fESDEvent->GetTrack(fPhotonCandidate->GetTrackLabelNegative());
1138         AliESDtrack* posTrack   = (AliESDtrack*)fESDEvent->GetTrack(fPhotonCandidate->GetTrackLabelPositive());
1139
1140         if(negProbArray && posProbArray){
1141
1142             negTrack->GetTPCpid(negProbArray);
1143             posTrack->GetTPCpid(posProbArray);
1144
1145             if(negProbArray[AliPID::kElectron]>=fPIDProbabilityCutNegativeParticle && posProbArray[AliPID::kElectron]>=fPIDProbabilityCutPositiveParticle){
1146                 iResult=kTRUE;
1147             }
1148             else{hV0CutPIDProb->Fill(fPhotonCandidate->GetPhotonMass());}
1149         }
1150
1151         delete [] posProbArray;
1152         delete [] negProbArray;
1153         return iResult;
1154
1155     }
1156     if(fAODEvent){
1157         // not possible to implement
1158         return kTRUE;}
1159     return kFALSE;
1160 }
1161
1162 ///________________________________________________________________________
1163 const AliVertex *AliV0ReaderV1::GetPrimaryVertex() {
1164     if(fESDEvent)return fESDEvent->GetPrimaryVertex();
1165     if(fAODEvent)return const_cast<const AliVertex*>(dynamic_cast<AliVertex*>(fAODEvent->GetPrimaryVertex()));
1166 return 0x0;
1167 }
1168
1169 ///________________________________________________________________________
1170 Bool_t AliV0ReaderV1::GetHelixCenter(const AliExternalTrackParam *track, Double_t b,Int_t charge, Double_t center[2]){
1171         // see header file for documentation
1172
1173         Double_t        helix[6];
1174         track->GetHelixParameters(helix,b);
1175
1176         Double_t xpos = helix[5];
1177         Double_t ypos = helix[0];
1178         Double_t radius = TMath::Abs(1./helix[4]);
1179         Double_t phi = helix[2];
1180
1181         if(phi < 0){
1182                 phi = phi + 2*TMath::Pi();
1183         }
1184
1185         phi -= TMath::Pi()/2.;
1186         Double_t xpoint =       radius * TMath::Cos(phi);
1187         Double_t ypoint =       radius * TMath::Sin(phi);
1188
1189         if(b<0){
1190                 if(charge > 0){
1191                         xpoint = - xpoint;
1192                         ypoint = - ypoint;
1193                 }
1194
1195                 if(charge < 0){
1196                         xpoint =        xpoint;
1197                         ypoint =        ypoint;
1198                 }
1199         }
1200         if(b>0){
1201                 if(charge > 0){
1202                         xpoint =        xpoint;
1203                         ypoint =        ypoint;
1204                 }
1205
1206                 if(charge < 0){
1207                         xpoint = - xpoint;
1208                         ypoint = - ypoint;
1209                 }
1210         }
1211         center[0] =     xpos + xpoint;
1212         center[1] =     ypos + ypoint;
1213
1214         return 1;
1215 }
1216 ///________________________________________________________________________
1217 Bool_t AliV0ReaderV1::GetConversionPoint(const AliExternalTrackParam *pparam,const AliExternalTrackParam *nparam,Double_t convpos[3]){
1218
1219     if(!pparam||!nparam)return kFALSE;
1220
1221         Double_t helixcenterpos[2];
1222         GetHelixCenter(pparam,GetMagneticField(),pparam->Charge(),helixcenterpos);
1223
1224         Double_t helixcenterneg[2];
1225         GetHelixCenter(nparam,GetMagneticField(),nparam->Charge(),helixcenterneg);
1226
1227         Double_t helixpos[6];
1228         pparam->GetHelixParameters(helixpos,GetMagneticField());
1229         Double_t posradius = TMath::Abs(1./helixpos[4]);
1230
1231         Double_t helixneg[6];
1232         nparam->GetHelixParameters(helixneg,GetMagneticField());
1233         Double_t negradius = TMath::Abs(1./helixneg[4]);
1234
1235         // Calculate xy-position
1236
1237         Double_t xpos = helixcenterpos[0];
1238         Double_t ypos = helixcenterpos[1];
1239         Double_t xneg = helixcenterneg[0];
1240         Double_t yneg = helixcenterneg[1];
1241
1242         convpos[0] = (xpos*negradius + xneg*posradius)/(negradius+posradius);
1243         convpos[1] = (ypos*negradius+   yneg*posradius)/(negradius+posradius);
1244
1245
1246         // Calculate z-position
1247
1248         Double_t deltaXPos = convpos[0] -       xpos;
1249          Double_t deltaYPos = convpos[1] -      ypos;
1250
1251          Double_t deltaXNeg = convpos[0] -      xneg;
1252          Double_t deltaYNeg = convpos[1] -      yneg;
1253
1254          Double_t alphaPos =    TMath::Pi() + TMath::ATan2(-deltaYPos,-deltaXPos);
1255          Double_t alphaNeg =    TMath::Pi() + TMath::ATan2(-deltaYNeg,-deltaXNeg);
1256
1257          Double_t vertexXNeg =  xneg +  TMath::Abs(negradius)*
1258          TMath::Cos(alphaNeg);
1259          Double_t vertexYNeg =  yneg +  TMath::Abs(negradius)*
1260          TMath::Sin(alphaNeg);
1261
1262          Double_t vertexXPos =  xpos +  TMath::Abs(posradius)*
1263          TMath::Cos(alphaPos);
1264          Double_t vertexYPos =  ypos +  TMath::Abs(posradius)*
1265          TMath::Sin(alphaPos);
1266
1267          Double_t x0neg =        helixneg[5];
1268          Double_t y0neg =        helixneg[0];
1269
1270          Double_t x0pos =        helixpos[5];
1271          Double_t y0pos =        helixpos[0];
1272
1273          Double_t dNeg = TMath::Sqrt((vertexXNeg -      x0neg)*(vertexXNeg - x0neg)
1274                                                                                                                          +(vertexYNeg - y0neg)*(vertexYNeg - y0neg));
1275
1276          Double_t dPos = TMath::Sqrt((vertexXPos -      x0pos)*(vertexXPos - x0pos)
1277                                                                                                                          +(vertexYPos - y0pos)*(vertexYPos - y0pos));
1278
1279          Double_t rNeg =        TMath::Sqrt(negradius*negradius -
1280          dNeg*dNeg/4.);
1281
1282          Double_t rPos = TMath::Sqrt(posradius*posradius -
1283          dPos*dPos/4.);
1284
1285          Double_t deltabetaNeg =        2*(TMath::Pi() +         TMath::ATan2(-dNeg/2.,-rNeg));
1286          Double_t deltabetaPos = 2*(TMath::Pi() + TMath::ATan2(-dPos/2.,-rPos));
1287
1288          Double_t deltaUNeg = negradius*deltabetaNeg;
1289          Double_t deltaUPos = posradius*deltabetaPos;
1290
1291          Double_t zphaseNeg = nparam->GetZ() +  deltaUNeg * nparam->GetTgl();
1292          Double_t zphasePos = pparam->GetZ() +  deltaUPos * pparam->GetTgl();
1293
1294          convpos[2] = (zphasePos*negradius+zphaseNeg*posradius)/(negradius+posradius);
1295
1296 return kTRUE;
1297 }
1298 ///________________________________________________________________________
1299 void AliV0ReaderV1::ProcessMC(AliKFConversionPhoton *fCurrentReconstructedGamma){
1300
1301     if(!fMCStack)return;
1302
1303     TParticle *fMotherMCParticle=NULL;
1304     TParticle *fNegativeMCParticle=NULL;
1305     TParticle *fPositiveMCParticle=NULL;
1306
1307     // Get MC Particles
1308     fMotherMCParticle   = fCurrentReconstructedGamma->GetMCParticle(fMCStack);
1309     fNegativeMCParticle = fCurrentReconstructedGamma->GetNegativeMCDaughter(fMCStack);
1310     fPositiveMCParticle = fCurrentReconstructedGamma->GetPositiveMCDaughter(fMCStack);
1311
1312     if(fPositiveMCParticle&&fNegativeMCParticle&&fMotherMCParticle){
1313
1314         // Check if it is a true photon
1315
1316         if(fMotherMCParticle->GetPdgCode()==22){
1317
1318             hMCPtRECOTRUE[fCentralityBin]->Fill(fCurrentReconstructedGamma->GetPt());
1319
1320             // Pt Resolution
1321
1322             Double_t mcpt        = fMotherMCParticle->Pt();
1323             Double_t esdpt      = fCurrentReconstructedGamma->GetPt();
1324             Double_t resdPt = 0.;
1325             if(mcpt > 0){
1326                 resdPt = ((esdpt - mcpt)/mcpt)*100.;
1327             } else if(mcpt < 0){
1328                 AliWarning("Pt of MC particle is negative, this will cause wrong calculation of resPt");
1329             }
1330
1331
1332             hMCPtResolution->Fill(mcpt,resdPt);
1333             hMCPtResolutionPhi->Fill(fMotherMCParticle->Phi(),resdPt);
1334
1335             // Conversion Point Resolution
1336
1337             Double_t resdR = 0.;
1338             if(fNegativeMCParticle->R() != 0){
1339                 resdR = ((fCurrentReconstructedGamma->GetConversionRadius() - fNegativeMCParticle->R())/fNegativeMCParticle->R())*100.;
1340             }
1341             Double_t resdRAbs = 0.;
1342             resdRAbs = (fCurrentReconstructedGamma->GetConversionRadius() - fNegativeMCParticle->R());
1343
1344             hMCRResolutionvsR->Fill(fNegativeMCParticle->R(),resdRAbs);
1345
1346             //          fHistograms->FillHistogram("Resolution_dR", fV0Reader->GetNegativeMCParticle()->R(), resdR);
1347             //          fHistograms->FillHistogram("Resolution_MC_R", fV0Reader->GetNegativeMCParticle()->R());
1348             //          fHistograms->FillHistogram("Resolution_ESD_R", fV0Reader->GetXYRadius());
1349             //          fHistograms->FillHistogram("Resolution_R_dPt", fV0Reader->GetNegativeMCParticle()->R(), resdPt);
1350
1351             Double_t resdZ = 0.;
1352             if(fNegativeMCParticle->Vz() != 0){
1353                 resdZ = ((fCurrentReconstructedGamma->GetZ() -fNegativeMCParticle->Vz())/fNegativeMCParticle->Vz())*100.;
1354             }
1355             Double_t resdZAbs = 0.;
1356             resdZAbs = fCurrentReconstructedGamma->GetZ() -fNegativeMCParticle->Vz();
1357
1358             hMCZResolutionvsZ->Fill( fNegativeMCParticle->Vz(), resdZAbs);
1359         }
1360     }
1361
1362 }
1363
1364 ///________________________________________________________________________
1365 void AliV0ReaderV1::ProcessMCGammasForEfficiency(){
1366
1367     if(!fMCStack)return;
1368
1369     for (Int_t iTracks = 0; iTracks < fMCStack->GetNprimary(); iTracks++) {
1370         TParticle* particle = (TParticle *)fMCStack->Particle(iTracks);
1371
1372         //process the gammas
1373
1374         if(IsMCConversionGammaInAcceptance(particle)){
1375
1376             hMCPtTRUE[fCentralityBin]->Fill(particle->Pt());
1377
1378         }
1379     }
1380 }
1381
1382
1383 ///________________________________________________________________________
1384 Bool_t AliV0ReaderV1::IsMCConversionGammaInAcceptance(TParticle *particle){
1385     if(!fMCStack)return kFALSE;
1386
1387     if (particle->GetPdgCode() == 22){
1388         if(TMath::Abs(particle->Eta())> fEtaCut || TMath::Abs(particle->Eta())< fEtaCutMin)     return kFALSE;
1389
1390         if(particle->GetMother(0) >-1 && fMCStack->Particle(particle->GetMother(0))->GetPdgCode() == 22){
1391             return kFALSE; // no photon as mothers!
1392         }
1393
1394         if(particle->GetMother(0) >= fMCStack->GetNprimary()){
1395             return kFALSE; // the gamma has a mother, and it is not a primary particle
1396         }
1397
1398         // looking for conversion (electron + positron from pairbuilding (= 5) )
1399         TParticle* ePos = NULL;
1400         TParticle* eNeg = NULL;
1401
1402         if(particle->GetNDaughters() >= 2){
1403             for(Int_t daughterIndex=particle->GetFirstDaughter();daughterIndex<=particle->GetLastDaughter();daughterIndex++){
1404                 TParticle *tmpDaughter = fMCStack->Particle(daughterIndex);
1405                 if(tmpDaughter->GetUniqueID() == 5){
1406                     if(tmpDaughter->GetPdgCode() == 11){
1407                         eNeg = tmpDaughter;
1408                     } else if(tmpDaughter->GetPdgCode() == -11){
1409                         ePos = tmpDaughter;
1410                     }
1411                 }
1412             }
1413         }
1414
1415         if(ePos == NULL || eNeg == NULL){ // means we do not have two daughters from pair production
1416             return kFALSE;
1417         }
1418
1419         if(AcceptanceCut(particle,ePos,eNeg))return kTRUE;
1420     }
1421     return kFALSE;
1422 }
1423
1424 ///________________________________________________________________________
1425 Bool_t AliV0ReaderV1::AcceptanceCut(TParticle *particle, TParticle * ePos,TParticle* eNeg){
1426
1427     // cuts on distance from collision point
1428
1429     if(particle->R()>fMaxR){
1430         return kFALSE;}
1431
1432     if(ePos->R()>fMaxR){
1433         return kFALSE;
1434     }
1435
1436     if(ePos->R()<fMinR){
1437         return kFALSE;
1438     }
1439
1440     if( ePos->R() <= ((TMath::Abs(ePos->Vz())*fLineCutZRSlope)-fLineCutZValue)){
1441         return kFALSE;
1442     }
1443     /*else if (fUseEtaMinCut &&  ePos->R() >= ((TMath::Abs(ePos->Vz())*fLineCutZRSlopeMin)-fLineCutZValueMin )){
1444         return kFALSE;
1445     } */
1446
1447     if(TMath::Abs(eNeg->Vz()) > fMaxZ ){ // cuts out regions where we do not reconstruct
1448         return kFALSE;
1449     }
1450
1451     if(eNeg->Vz()!=ePos->Vz()||eNeg->R()!=ePos->R()){
1452         return kFALSE;
1453     }
1454
1455     if(TMath::Abs(ePos->Vz()) > fMaxZ ){ // cuts out regions where we do not reconstruct
1456         return kFALSE;
1457     }
1458
1459     if(TMath::Abs(particle->Eta())> fEtaCut || TMath::Abs(particle->Eta())< fEtaCutMin){
1460         return kFALSE;
1461     }
1462
1463     if(TMath::Abs(ePos->Eta())> fEtaCut || TMath::Abs(ePos->Eta())< fEtaCutMin){
1464         return kFALSE;
1465     }
1466
1467     if(TMath::Abs(eNeg->Eta())> fEtaCut || TMath::Abs(eNeg->Eta())< fEtaCutMin){
1468         return kFALSE;
1469     }
1470
1471     if( ePos->Pt()< fSinglePtCut ||  eNeg->Pt()< fSinglePtCut){
1472         return kFALSE;
1473     }
1474
1475     if(particle->Pt()<fPtCut){
1476         return kFALSE;
1477     }
1478
1479     return kTRUE;
1480 }
1481 ///________________________________________________________________________
1482 void AliV0ReaderV1::PrintCuts(){
1483
1484     cout<<"V0 Reader initialized with following settings"<<endl;
1485
1486
1487     cout<<"Acceptance Eta:"<<endl;
1488     cout<<fEtaCutMin<<" < eta < "<<fEtaCut<<endl;
1489     cout<<"Conversion Point"<<endl;
1490     cout<<"Z <"<<fMaxZ<<endl;
1491     cout<<fMinR<<" < R < "<<fMaxR<<endl;
1492     cout<<"Line Cut Slope"<<fLineCutZRSlope<<" ZValue "<<fLineCutZValue<<endl;
1493
1494     cout<<"Pt Gamma > "<<fPtCut<<endl;
1495     cout<<"Pt Daughters > "<<fSinglePtCut<<endl;
1496
1497     cout<<"Armenteros Qt Cut"<<endl;
1498
1499     if(fDoHighPtQtGammaSelection){
1500         cout<<" qt < "<<fQtMax<<" for pt < "<<fPtBorderForQt<<endl;
1501         cout<<" qt < "<<fHighPtQtMax<<" for pt > "<<fPtBorderForQt<<endl;
1502     }
1503     else{
1504          cout<<" qt < "<<fQtMax<<endl;
1505     }
1506     cout<<"Chi2perNDF > "<<fChi2CutConversion<<endl;
1507
1508
1509
1510 }
1511
1512 //_______________________________________________________________________
1513
1514 Bool_t AliV0ReaderV1::CentralitySelection(){
1515
1516     fCentralityBin=-1;
1517
1518     if(!fIsHeavyIon){fCentralityBin=0;fCentrality=0;return kTRUE;}
1519
1520     fCentrality=GetCentrality();
1521
1522     if(fIsHeavyIon){
1523
1524         if(fBGHandler){
1525             fCentralityBin=fBGHandler->GetCentralityBinIndex(Int_t(fCentrality));
1526         }
1527         else{
1528             Double_t fCentralityBins[fNCentralityBins+1];
1529             for(int i=0;i<fNCentralityBins;i++){
1530                 fCentralityBins[i]=i*100/Double_t(fNCentralityBins);
1531             }
1532
1533             for(int i=0;i<fNCentralityBins;i++){
1534                 if(fCentrality>fCentralityBins[i]&&fCentrality<fCentralityBins[i+1]){
1535                     fCentralityBin=i;
1536                     return kTRUE;}
1537             }
1538         }
1539
1540     }
1541     if(fCentralityBin>=0&&fCentrality>=0){
1542
1543         return kTRUE;
1544     }
1545
1546     AliWarning("Centrality not defined");
1547     return kFALSE;
1548
1549 }
1550
1551 //________________________________________________________________________
1552 void AliV0ReaderV1::Terminate(Option_t *)
1553 {
1554
1555     printf(Form("V0ReaderV1: V0s processed: %4.0f reconstructed photons: %4.0f \n",hV0CurrentFinder->GetEntries(),hV0Good->GetEntries()));
1556
1557 }