AliPhysics  64bcec2 (64bcec2)
AliV0ReaderV1.cxx
Go to the documentation of this file.
1 /**************************************************************************
2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3  * *
4  * Developers: Friederike Bock *
5  * Original Authors: Svein Lindal, Daniel Lohner *
6  * *
7  * Version 1.0 *
8  * *
9  * *
10  * based on: on older version (see aliroot up to v5-04-42-AN) *
11  * AliV0Reader.cxx *
12  * Authors: Kathrin Koch, Kenneth Aamodt, Ana Marin *
13  * *
14  * Permission to use, copy, modify and distribute this software and its *
15  * documentation strictly for non-commercial purposes is hereby granted *
16  * without fee, provided that the above copyright notice appears in all *
17  * copies and that both the copyright notice and this permission notic *
18  * appear in the supporting documentation. The authors make no claims *
19  * about the suitability of this software for any purpose. It is *
20  * provided "as is" without express or implied warranty. *
21  **************************************************************************/
22 
24 //---------------------------------------------
25 // Class reconstructing conversion photons from V0s
26 //---------------------------------------------
28 
29 // The AliAODConversionPhotons will return the Position (== ESDTrack->GetID())of the positive (negative) ESDTrack
30 // in the Array of ESDTracks stored in the ESDEvent via GetTrackLabelPositive() (GetTrackLabelNagative()).
31 // Since it is the Position in the Array it is always positive.
32 // After the conversion ESD --> AOD each AOD track will give you the former Position in the ESDArray via the
33 // Function AODTrack->GetID(). AODTracks are stored for different TrackParameter (e.g. HybridTracks, TPC only ...). No Standard
34 // AODTracks (not copies of AliExternalTrackParam) will be stored with negative values of AODTrack->GetID() (Formula: (ESDTrack->GetID()+1)*-1)
35 // This leads to the possibility of more than one AOD track with the same negative ID. For that reason all AOD tracks are additionally flaged.
36 // In this analysis we should therefore only use AODTracks with positive values of AODTrack->GetID().
37 
38 // If you want to find the AODTrack corresponding to the daugher track of a AliAODConversionPhoton you have to find the AODTrack with
39 // AODTrack->GetID() == GetTrackLabelPositive() (GetTrackLabelNagative()).
40 
41 #include <vector>
42 #include <TGeoGlobalMagField.h>
43 
44 #include "AliV0ReaderV1.h"
45 #include "AliKFParticle.h"
46 #include "AliAODv0.h"
47 #include "AliESDv0.h"
48 #include "AliAODEvent.h"
49 #include "AliESDEvent.h"
50 #include "AliPID.h"
51 #include "AliMCEvent.h"
52 #include "AliMCEventHandler.h"
53 #include "AliESDpid.h"
54 #include "AliESDtrackCuts.h"
55 #include "TRandom3.h"
56 #include "AliGenCocktailEventHeader.h"
57 #include "TList.h"
58 #include "AliKFConversionPhoton.h"
59 #include "AliAODConversionPhoton.h"
61 #include "TVector.h"
62 #include "AliKFVertex.h"
63 #include "AliAODTrack.h"
64 #include "AliESDtrack.h"
65 #include "AliAnalysisManager.h"
66 #include "AliInputEventHandler.h"
67 #include "AliAODHandler.h"
68 #include "AliPIDResponse.h"
69 #include "TChain.h"
70 #include "TFile.h"
71 #include "TString.h"
72 #include "TObjArray.h"
73 
74 class iostream;
75 
76 using namespace std;
77 
78 ClassImp(AliV0ReaderV1)
79 
80 //________________________________________________________________________
81 AliV0ReaderV1::AliV0ReaderV1(const char *name) : AliAnalysisTaskSE(name),
82  kAddv0sInESDFilter(kFALSE),
83  fPCMv0BitField(NULL),
84  fConversionCuts(NULL),
85  fEventCuts(NULL),
86  fConversionGammas(NULL),
87  fUseImprovedVertex(kTRUE),
88  fUseOwnXYZCalculation(kTRUE),
89  fUseConstructGamma(kFALSE),
90  kUseAODConversionPhoton(kTRUE),
91  fCreateAOD(kFALSE),
92  fDeltaAODBranchName("GammaConv"),
93  fDeltaAODFilename("AliAODGammaConversion.root"),
94  fRelabelAODs(kFALSE),
95  fPreviousV0ReaderPerformsAODRelabeling(0),
96  fEventIsSelected(kFALSE),
97  fNumberOfPrimaryTracks(0),
98  fNumberOfTPCoutTracks(0),
99  fPeriodName(""),
100  fPtHardBin(0),
101  fUseMassToZero(kTRUE),
102  fProduceV0findingEffi(kFALSE),
103  fProduceImpactParamHistograms(kFALSE),
104  fCurrentInvMassPair(0),
105  fImprovedPsiPair(3),
106  fHistograms(NULL),
107  fImpactParamHistograms(NULL),
108  fHistoMCGammaPtvsR(NULL),
109  fHistoMCGammaPtvsPhi(NULL),
110  fHistoMCGammaPtvsEta(NULL),
111  fHistoMCGammaRvsPhi(NULL),
112  fHistoMCGammaRvsEta(NULL),
113  fHistoMCGammaPhivsEta(NULL),
114  fHistoRecMCGammaPtvsR(NULL),
115  fHistoRecMCGammaPtvsPhi(NULL),
116  fHistoRecMCGammaPtvsEta(NULL),
117  fHistoRecMCGammaRvsPhi(NULL),
118  fHistoRecMCGammaRvsEta(NULL),
119  fHistoRecMCGammaPhivsEta(NULL),
120  fHistoRecMCGammaMultiPt(NULL),
121  fHistoRecMCGammaMultiPtvsEta(NULL),
122  fHistoRecMCGammaMultiR(NULL),
123  fHistoRecMCGammaMultiPhi(NULL),
124  fHistoPosTrackImpactParamZ(NULL),
125  fHistoPosTrackImpactParamY(NULL),
126  fHistoPosTrackImpactParamX(NULL),
127  fHistoPosTrackImpactParamZvsPt(NULL),
128  fHistoPosTrackImpactParamYvsPt(NULL),
129  fHistoPosTrackImpactParamXvsPt(NULL),
130  fHistoNegTrackImpactParamZ(NULL),
131  fHistoNegTrackImpactParamY(NULL),
132  fHistoNegTrackImpactParamX(NULL),
133  fHistoNegTrackImpactParamZvsPt(NULL),
134  fHistoNegTrackImpactParamYvsPt(NULL),
135  fHistoNegTrackImpactParamXvsPt(NULL),
136  fHistoImpactParamZvsR(NULL),
137  fHistoImpactParamZvsR2(NULL),
138  fHistoPt(NULL),
139  fHistoPt2(NULL),
140  fHistoDCAzPhoton(NULL),
141  fHistoDCAzPhoton2(NULL),
142  fHistoR(NULL),
143  fHistoRrecalc(NULL),
144  fHistoRviaAlpha(NULL),
145  fHistoRviaAlphaRecalc(NULL),
146  fHistoRdiff(NULL),
147  fHistoImpactParameterStudy(NULL),
148  fImpactParamTree(NULL),
149  fVectorFoundGammas(0),
150  fCurrentFileName(""),
151  fMCFileChecked(kFALSE)
152 {
153  // Default constructor
154 
155  DefineInput(0, TChain::Class());
156  DefineOutput(1,TBits::Class());
157 }
158 
159 //________________________________________________________________________
161 {
162  // default deconstructor
163 
164  if(fConversionGammas){
165  fConversionGammas->Delete();// Clear Objects
166  delete fConversionGammas;
167  fConversionGammas=0x0;
168  }
169 }
170 
180  if(fConversionGammas && (index < fConversionGammas->GetEntriesFast())) {
181  TObject *tmp = fConversionGammas->At(index);
182  if(tmp->IsA() == AliAODConversionPhoton::Class()){
183  AliAODConversionPhoton *photon = static_cast<AliAODConversionPhoton *>(tmp);
184  return photon;
185  } else {
186  AliKFConversionPhoton *photon = static_cast<AliKFConversionPhoton *>(tmp);
187  return photon;
188  }
189  } else {
190  return nullptr;
191  }
192 }
193 
194 /*
195 //________________________________________________________________________
196 AliV0ReaderV1::AliV0ReaderV1(AliV0ReaderV1 &original) : AliAnalysisTaskSE(original),
197 fConversionCuts(NULL),
198 fConversionGammas(NULL),
199 fUseImprovedVertex(original.fUseImprovedVertex),
200 fUseOwnXYZCalculation(original.fUseOwnXYZCalculation),
201 fUseConstructGamma(original.fUseConstructGamma),
202 kUseAODConversionPhoton(original.kUseAODConversionPhoton),
203 fCreateAOD(original.fCreateAOD),
204 fDeltaAODBranchName(original.fDeltaAODBranchName),
205 fDeltaAODFilename(original.fDeltaAODFilename),
206 fEventIsSelected(original.fEventIsSelected)
207 {
208 // Default constructor
209 
210 DefineInput(0, TChain::Class());
211 }
212 
213 //____________________________________________________________
214 AliV0ReaderV1 &AliV0ReaderV1::operator=(const AliV0ReaderV1 &ref){
215 //
216 // Assignment operator
217 // Only copies pointers, object is not the owner of the references
218 //
219 if(this != &ref){
220 AliAnalysisTaskSE::operator=(ref);
221 fUseImprovedVertex=ref.fUseImprovedVertex;
222 fUseOwnXYZCalculation=ref.fUseOwnXYZCalculation;
223 fUseConstructGamma=ref.fUseConstructGamma;
224 kUseAODConversionPhoton=ref.kUseAODConversionPhoton;
225 fCreateAOD=ref.fCreateAOD;
226 fDeltaAODBranchName=ref.fDeltaAODBranchName;
227 fDeltaAODFilename=ref.fDeltaAODFilename;
228 fEventIsSelected=ref.fEventIsSelected;
229 }
230 return *this;
231 }
232 */
233 
234 //________________________________________________________________________
236 {
237  // Initialize function to be called once before analysis
238  if(fConversionCuts==NULL){
239  if(fConversionCuts==NULL)AliError("No Conversion Cut Selection initialized");
240  }
241  if(fEventCuts==NULL){
242  if(fEventCuts==NULL)AliError("No Event Cut Selection initialized");
243  }
244 
245  if(fCreateAOD){kUseAODConversionPhoton=kTRUE;}
246 
247  if(fConversionGammas != NULL){
248  delete fConversionGammas;
249  fConversionGammas=NULL;
250  }
251 
252  if(fConversionGammas == NULL){
253  if(kUseAODConversionPhoton){
254  fConversionGammas = new TClonesArray("AliAODConversionPhoton",100);}
255  else{
256  fConversionGammas = new TClonesArray("AliKFConversionPhoton",100);}
257  }
258  fConversionGammas->Delete();//Reset the TClonesArray
259 }
260 
261 //________________________________________________________________________
263 {
264  // Create AODs
265 
266  if(fCreateAOD){
267  fPCMv0BitField = new TBits();
268 
269  if (fEventCuts){
270  fDeltaAODBranchName.Append("_");
271  fDeltaAODBranchName.Append(fEventCuts->GetCutNumber());
272  }
273  if(fConversionCuts){
274  fDeltaAODBranchName.Append("_");
275  fDeltaAODBranchName.Append(fConversionCuts->GetCutNumber());
276  fDeltaAODBranchName.Append("_gamma");
277  }
278  fConversionGammas->SetName(fDeltaAODBranchName.Data());
279 
280  AddAODBranch("TClonesArray", &fConversionGammas, fDeltaAODFilename.Data());
281  AliAnalysisManager::GetAnalysisManager()->RegisterExtraFile(fDeltaAODFilename.Data());
282  }
283 
284  if(fProduceImpactParamHistograms){
285  if(fImpactParamHistograms != NULL){
286  delete fImpactParamHistograms;
287  fImpactParamHistograms = NULL;
288  }
289  if(fImpactParamHistograms==NULL){
290  fImpactParamHistograms = new TList();
291  fImpactParamHistograms->SetOwner(kTRUE);
292  fImpactParamHistograms->SetName(Form("ImpactParamHistograms_%s_%s",fEventCuts->GetCutNumber().Data(),fConversionCuts->GetCutNumber().Data()));
293  }
294  fHistoPosTrackImpactParamZ = new TH1F("fHistoPosTrackImpactParamZ","",480,-80,80);
295  fHistoPosTrackImpactParamZ->SetXTitle("Z (cm)");
296  fImpactParamHistograms->Add(fHistoPosTrackImpactParamZ);
297 
298  fHistoPosTrackImpactParamY = new TH1F("fHistoPosTrackImpactParamY","",720,-120,120);
299  fHistoPosTrackImpactParamY->SetXTitle("Y (cm)");
300  fImpactParamHistograms->Add(fHistoPosTrackImpactParamY);
301 
302  fHistoPosTrackImpactParamX = new TH1F("fHistoPosTrackImpactParamX","",30,-3,3);
303  fHistoPosTrackImpactParamX->SetXTitle("X (cm)");
304  fImpactParamHistograms->Add(fHistoPosTrackImpactParamX);
305 
306  fHistoNegTrackImpactParamZ = new TH1F("fHistoNegTrackImpactParamZ","",480,-80,80);
307  fHistoNegTrackImpactParamZ->SetXTitle("Z (cm)");
308  fImpactParamHistograms->Add(fHistoNegTrackImpactParamZ);
309 
310  fHistoNegTrackImpactParamY = new TH1F("fHistoNegTrackImpactParamY","",720,-120,120);
311  fHistoNegTrackImpactParamY->SetXTitle("Y (cm)");
312  fImpactParamHistograms->Add(fHistoNegTrackImpactParamY);
313 
314  fHistoNegTrackImpactParamX = new TH1F("fHistoNegTrackImpactParamX","",30,-3,3);
315  fHistoNegTrackImpactParamX->SetXTitle("X (cm)");
316  fImpactParamHistograms->Add(fHistoNegTrackImpactParamX);
317 
318  fHistoPosTrackImpactParamZvsPt = new TH2F("fHistoPosTrackImpactParamZvsPt","",100,0,10,480,-80,80);
319  fHistoPosTrackImpactParamZvsPt->SetYTitle("Z (cm)");
320  fHistoPosTrackImpactParamZvsPt->SetXTitle("Pt (GeV)");
321  fImpactParamHistograms->Add(fHistoPosTrackImpactParamZvsPt);
322 
323  fHistoPosTrackImpactParamYvsPt = new TH2F("fHistoPosTrackImpactParamYvsPt","",100,0,10,720,-120,120);
324  fHistoPosTrackImpactParamYvsPt->SetYTitle("Y (cm)");
325  fHistoPosTrackImpactParamYvsPt->SetXTitle("Pt (GeV)");
326  fImpactParamHistograms->Add(fHistoPosTrackImpactParamYvsPt);
327 
328  fHistoPosTrackImpactParamXvsPt = new TH2F("fHistoPosTrackImpactParamXvsPt","",100,0,10,30,-3,5);
329  fHistoPosTrackImpactParamXvsPt->SetYTitle("X (cm)");
330  fHistoPosTrackImpactParamXvsPt->SetXTitle("Pt (GeV)");
331  fImpactParamHistograms->Add(fHistoPosTrackImpactParamXvsPt);
332 
333  fHistoNegTrackImpactParamZvsPt = new TH2F("fHistoNegTrackImpactParamZvsPt","",100,0,10,480,-80,80);
334  fHistoNegTrackImpactParamZvsPt->SetYTitle("Z (cm)");
335  fHistoNegTrackImpactParamZvsPt->SetXTitle("Pt (GeV)");
336  fImpactParamHistograms->Add(fHistoNegTrackImpactParamZvsPt);
337 
338  fHistoNegTrackImpactParamYvsPt = new TH2F("fHistoNegTrackImpactParamYvsPt","",100,0,10,720,-120,120);
339  fHistoNegTrackImpactParamYvsPt->SetYTitle("Y (cm)");
340  fHistoNegTrackImpactParamYvsPt->SetXTitle("Pt (GeV)");
341  fImpactParamHistograms->Add(fHistoNegTrackImpactParamYvsPt);
342 
343  fHistoNegTrackImpactParamXvsPt = new TH2F("fHistoNegTrackImpactParamXvsPt","",100,0,10,30,-3,3);
344  fHistoNegTrackImpactParamXvsPt->SetYTitle("X (cm)");
345  fHistoNegTrackImpactParamXvsPt->SetXTitle("Pt (GeV)");
346  fImpactParamHistograms->Add(fHistoNegTrackImpactParamXvsPt);
347 
348  fHistoImpactParamZvsR = new TH2F("fHistoImpactParamZvsR","Before cuts",300,-150,150,200,0,200);
349  fHistoImpactParamZvsR->SetXTitle("Z (cm)");
350  fHistoImpactParamZvsR->SetYTitle("R (cm)");
351  fImpactParamHistograms->Add(fHistoImpactParamZvsR);
352 
353  fHistoImpactParamZvsR2 = new TH2F("fHistoImpactParamZvsR2","After cuts",300,-150,150,200,0,200);
354  fHistoImpactParamZvsR2->SetXTitle("Z (cm)");
355  fHistoImpactParamZvsR2->SetYTitle("R (cm)");
356  fImpactParamHistograms->Add(fHistoImpactParamZvsR2);
357 
358  fHistoPt = new TH1F("fHistoPt","Before all cuts",100,0,10);
359  fHistoPt->SetXTitle("Pt (GeV)");
360  fImpactParamHistograms->Add(fHistoPt);
361 
362  fHistoPt2 = new TH1F("fHistoPt2","After all cuts",100,0,10);
363  fHistoPt2->SetXTitle("Pt (GeV)");
364  fImpactParamHistograms->Add(fHistoPt2);
365 
366  fHistoDCAzPhoton = new TH1F("fHistoDCAzPhoton","Before cuts",20,-2,2);
367  fHistoDCAzPhoton->SetXTitle("DCAz photon (cm)");
368  fImpactParamHistograms->Add(fHistoDCAzPhoton);
369 
370  fHistoDCAzPhoton2 = new TH1F("fHistoDCAzPhoton2","After cuts",20,-2,2);
371  fHistoDCAzPhoton2->SetXTitle("DCAz photon (cm)");
372  fImpactParamHistograms->Add(fHistoDCAzPhoton2);
373 
374  fHistoR = new TH1F("fHistoR","",200,0,200);
375  fHistoR->SetXTitle("Conversion radius (cm)");
376  fImpactParamHistograms->Add(fHistoR);
377 
378  fHistoRrecalc = new TH1F("fHistoRrecalc","",200,0,200);
379  fHistoRrecalc->SetXTitle("conversion radius (cm)");
380  fImpactParamHistograms->Add(fHistoRrecalc);
381 
382  fHistoRviaAlpha = new TH1F("fHistoRviaAlpha","",200,0,200);
383  fHistoRviaAlpha->SetXTitle("Conversion radius (cm)");
384  fImpactParamHistograms->Add(fHistoRviaAlpha);
385 
386  fHistoRviaAlphaRecalc = new TH1F("fHistoRviaAlphaRecalc","",200,0,200);
387  fHistoRviaAlphaRecalc->SetXTitle("conversion radius (cm)");
388  fImpactParamHistograms->Add(fHistoRviaAlphaRecalc);
389 
390  fHistoRdiff = new TH1F("fHistoRdiff","",200,0,200);
391  fHistoRdiff->SetXTitle("R_conv - R_cluster conflict (cm)");
392  fImpactParamHistograms->Add(fHistoRdiff);
393 
394  fHistoImpactParameterStudy = new TH1F("fHistoImpactParameterStudy","",7,-0.5,6.5);
395  fHistoImpactParameterStudy->GetXaxis()->SetBinLabel(1,"# V0s");
396  fHistoImpactParameterStudy->GetXaxis()->SetBinLabel(2,"two TPC-only tracks");
397  fHistoImpactParameterStudy->GetXaxis()->SetBinLabel(3,"Z cut not passed");
398  fHistoImpactParameterStudy->GetXaxis()->SetBinLabel(4,"Y cut not passed");
399  fHistoImpactParameterStudy->GetXaxis()->SetBinLabel(5,"R>80cm");
400  fHistoImpactParameterStudy->GetXaxis()->SetBinLabel(6,"causality cut not p.");
401  fHistoImpactParameterStudy->GetXaxis()->SetBinLabel(7,"# removed V0s"); // of all V0s
402  fImpactParamHistograms->Add(fHistoImpactParameterStudy);
403 
404  fImpactParamTree = new TTree("fImpactParamTree","");
405  fImpactParamHistograms->Add(fImpactParamTree);
406  }
407 
408  if (fProduceV0findingEffi){
409  TH1::AddDirectory(kFALSE);
410  if(fHistograms != NULL){
411  delete fHistograms;
412  fHistograms = NULL;
413  }
414  if(fHistograms==NULL){
415  fHistograms = new TList();
416  fHistograms->SetOwner(kTRUE);
417  fHistograms->SetName(Form("V0FindingEfficiencyInput_%s_%s",fEventCuts->GetCutNumber().Data(),fConversionCuts->GetCutNumber().Data()));
418  }
419 
420  fHistoMCGammaPtvsR = new TH2F("MCconvGamma_Pt_R","MC converted gamma Pt vs R (|eta| < 0.9)",250,0.0,25,400,0,200);
421  fHistoMCGammaPtvsR->SetXTitle("p_{MC,T} (GeV/c)");
422  fHistoMCGammaPtvsR->SetYTitle("R_{MC,conv} (cm)");
423  fHistograms->Add(fHistoMCGammaPtvsR);
424 
425  fHistoMCGammaPtvsEta = new TH2F("MCconvGamma_Pt_Eta","MC converted gamma Pt vs Eta ",250,0.0,25,280,-1.4,1.4);
426  fHistoMCGammaPtvsEta->SetXTitle("p_{MC,T} (GeV/c)");
427  fHistoMCGammaPtvsEta->SetYTitle("#eta_{MC}");
428  fHistograms->Add(fHistoMCGammaPtvsEta);
429 
430  fHistoMCGammaPtvsPhi = new TH2F("MCconvGamma_Pt_Phi","MC converted gamma Pt vs Phi (|eta| < 0.9) ",250,0.0,25,400,0,2*TMath::Pi());
431  fHistoMCGammaPtvsPhi->SetXTitle("p_{MC,T} (GeV/c)");
432  fHistoMCGammaPtvsPhi->SetYTitle("#varphi_{MC} (rad)");
433  fHistograms->Add(fHistoMCGammaPtvsPhi);
434 
435  fHistoMCGammaRvsPhi = new TH2F("MCconvGamma_R_Phi","MC converted gamma R vs Phi (|eta| < 0.9) ",400,0,200,400,0,2*TMath::Pi());
436  fHistoMCGammaRvsPhi->SetXTitle("R_{MC,conv} (cm)");
437  fHistoMCGammaRvsPhi->SetYTitle("#varphi_{MC} (rad)");
438  fHistograms->Add(fHistoMCGammaRvsPhi);
439 
440  fHistoMCGammaRvsEta = new TH2F("MCconvGamma_R_Eta","MC converted gamma R vs Eta ",400,0,200,280,-1.4,1.4);
441  fHistoMCGammaRvsEta->SetXTitle("R_{MC,conv} (cm)");
442  fHistoMCGammaRvsEta->SetYTitle("#eta_{MC}");
443  fHistograms->Add(fHistoMCGammaRvsEta);
444 
445  fHistoMCGammaPhivsEta = new TH2F("MCconvGamma_Phi_Eta","MC converted gamma Phi vs Eta ",400,0,2*TMath::Pi(),280,-1.4,1.4);
446  fHistoMCGammaPhivsEta->SetXTitle("#phi_{MC} (rad)");
447  fHistoMCGammaPhivsEta->SetYTitle("#eta_{MC}");
448  fHistograms->Add(fHistoMCGammaPhivsEta);
449 
450  fHistoRecMCGammaPtvsR = new TH2F("RecMCconvGamma_Pt_R","rec MC converted gamma Pt vs R (|eta| < 0.9)",250,0.0,25,400,0,200);
451  fHistoRecMCGammaPtvsR->SetXTitle("p_{MC,T} (GeV/c)");
452  fHistoRecMCGammaPtvsR->SetYTitle("R_{MC,conv} (cm)");
453  fHistograms->Add(fHistoRecMCGammaPtvsR);
454 
455  fHistoRecMCGammaPtvsEta = new TH2F("RecMCconvGamma_Pt_Eta","rec MC converted gamma Pt vs Eta ",250,0.0,25,280,-1.4,1.4);
456  fHistoRecMCGammaPtvsEta->SetXTitle("p_{MC,T} (GeV/c)");
457  fHistoRecMCGammaPtvsEta->SetYTitle("#eta_{MC}");
458  fHistograms->Add(fHistoRecMCGammaPtvsEta);
459 
460  fHistoRecMCGammaPtvsPhi = new TH2F("RecMCconvGamma_Pt_Phi","rec MC converted gamma Pt vs Phi (|eta| < 0.9) ",250,0.0,25,400,0,2*TMath::Pi());
461  fHistoRecMCGammaPtvsPhi->SetXTitle("p_{MC,T} (GeV/c)");
462  fHistoRecMCGammaPtvsPhi->SetYTitle("#varphi_{MC} (rad)");
463  fHistograms->Add(fHistoRecMCGammaPtvsPhi);
464 
465  fHistoRecMCGammaRvsPhi = new TH2F("RecMCconvGamma_R_Phi","rec MC converted gamma R vs Phi (|eta| < 0.9) ",400,0,200,400,0,2*TMath::Pi());
466  fHistoRecMCGammaRvsPhi->SetXTitle("R_{MC,conv} (cm)");
467  fHistoRecMCGammaRvsPhi->SetYTitle("#varphi_{MC} (rad)");
468  fHistograms->Add(fHistoRecMCGammaRvsPhi);
469 
470  fHistoRecMCGammaRvsEta = new TH2F("RecMCconvGamma_R_Eta","rec MC converted gamma R vs Eta ",400,0,200,280,-1.4,1.4);
471  fHistoRecMCGammaRvsEta->SetXTitle("R_{MC,conv} (cm)");
472  fHistoRecMCGammaRvsEta->SetYTitle("#eta_{MC}");
473  fHistograms->Add(fHistoRecMCGammaRvsEta);
474 
475  fHistoRecMCGammaPhivsEta = new TH2F("RecMCconvGamma_Phi_Eta","rec MC converted gamma Phi vs Eta ",400,0,2*TMath::Pi(),280,-1.4,1.4);
476  fHistoRecMCGammaPhivsEta->SetXTitle("#phi_{MC} (rad)");
477  fHistoRecMCGammaPhivsEta->SetYTitle("#eta_{MC}");
478  fHistograms->Add(fHistoRecMCGammaPhivsEta);
479 
480  fHistoRecMCGammaMultiPtvsEta = new TH2F("RecMCconvGammaMulti_Pt_Eta","rec MC converted gamma (at least double counted) Pt vs Eta ",250,0.0,25,280,-1.4,1.4);
481  fHistoRecMCGammaMultiPtvsEta->SetXTitle("p_{MC,T} (GeV/c)");
482  fHistoRecMCGammaMultiPtvsEta->SetYTitle("#eta_{MC}");
483  fHistograms->Add(fHistoRecMCGammaMultiPtvsEta);
484 
485  fHistoRecMCGammaMultiPt = new TH1F("RecMCconvGammaMulti_Pt","rec MC converted gamma (at least double counted) Pt (|eta| < 0.9)",250,0.0,25);
486  fHistoRecMCGammaMultiPt->SetXTitle("p_{MC,T} (GeV/c)");
487  fHistograms->Add(fHistoRecMCGammaMultiPt);
488 
489  fHistoRecMCGammaMultiR = new TH1F("RecMCconvGammaMulti_R","rec MC converted gamma (at least double counted) R (|eta| < 0.9)",400,0,200);
490  fHistoRecMCGammaMultiR->SetXTitle("R_{MC,conv} (cm)");
491  fHistograms->Add(fHistoRecMCGammaMultiR);
492 
493  fHistoRecMCGammaMultiPhi = new TH1F("RecMCconvGammaMulti_Phi","rec MC converted gamma (at least double counted) Phi (|eta| < 0.9)",400,0,2*TMath::Pi());
494  fHistoRecMCGammaMultiPhi->SetXTitle("#phi_{MC} (rad)");
495  fHistograms->Add(fHistoRecMCGammaMultiPhi);
496 
497  fVectorFoundGammas.clear();
498  }
499 
500 }
501 //________________________________________________________________________
503 
504  // set file name to empty & reset check flag
505  fCurrentFileName = "";
506  fMCFileChecked = kFALSE;
507 
508  // obtain file name from analysis manager
509  AliAnalysisManager *man=AliAnalysisManager::GetAnalysisManager();
510  if(man) {
511  AliInputEventHandler* inputHandler = (AliInputEventHandler*) (man->GetInputEventHandler());
512  if (inputHandler){
513  TTree* tree = (TTree*) inputHandler->GetTree();
514  TFile* file = (TFile*) tree->GetCurrentFile();
515  fCurrentFileName = file->GetName();
516  }
517  }
518 
519  // try to extract period name from file name if not given
520  if (fPeriodName.CompareTo("") == 0){
521  TObjArray *arr = fCurrentFileName.Tokenize("/");
522  fPtHardBin = -1;
523  for (Int_t i = 0; i < arr->GetEntriesFast();i++ ){
524  TObjString* testObjString = (TObjString*)arr->At(i);
525  if (testObjString->GetString().BeginsWith("LHC")){
526  fPeriodName = testObjString->GetString();
527  i = arr->GetEntriesFast();
528  }
529  }
530  if (fPeriodName.CompareTo("")==0){
531  TObjArray *arr2 = fCurrentFileName.Tokenize("__");
532  for (Int_t i = 0; i < arr->GetEntriesFast();i++ ){
533  TObjString* testObjString = (TObjString*)arr2->At(i);
534  if (testObjString->GetString().BeginsWith("LHC")){
535  fPeriodName = testObjString->GetString();
536  i = arr2->GetEntriesFast();
537  }
538  }
539  }
540  if (fPeriodName.CompareTo("") != 0 && fEventCuts->GetPeriodEnum() == AliConvEventCuts::kNoPeriod ){
541  fEventCuts->SetPeriodEnum (fPeriodName);
542  }
543  } else {
544  if(fEventCuts->GetPeriodEnum() == AliConvEventCuts::kNoPeriod ){
545  fEventCuts->SetPeriodEnum (fPeriodName);
546  }
547  }
548  // set pt hard bin from file name
549  if ( fEventCuts->GetPeriodEnum() == AliConvEventCuts::kLHC15a3a || fEventCuts->GetPeriodEnum() == AliConvEventCuts::kLHC15a3a_plus ||
550  fEventCuts->GetPeriodEnum() == AliConvEventCuts::kLHC15a3b ||
551  fEventCuts->GetPeriodEnum() == AliConvEventCuts::kLHC15g1a || fEventCuts->GetPeriodEnum() == AliConvEventCuts::kLHC15g1b ||
552  fEventCuts->GetPeriodEnum() == AliConvEventCuts::kLHC12P2JJ ||
553  fEventCuts->GetPeriodEnum() == AliConvEventCuts::kLHC16c3a || fEventCuts->GetPeriodEnum() == AliConvEventCuts::kLHC16c3b ||
554  fEventCuts->GetPeriodEnum() == AliConvEventCuts::kLHC16c3c ||
555  fEventCuts->GetPeriodEnum() == AliConvEventCuts::kLHC16P1JJ || fEventCuts->GetPeriodEnum() == AliConvEventCuts::kLHC16P1JJLowB ||
556  fEventCuts->GetPeriodEnum() == AliConvEventCuts::kLHC18b8
557  ){
558  TObjArray *arr = fCurrentFileName.Tokenize("/");
559  fPtHardBin = -1;
560  for (Int_t i = 0; i < arr->GetEntriesFast();i++ ){
561  TObjString* testObjString = (TObjString*)arr->At(i);
562  if (testObjString->GetString().BeginsWith("LHC")){
563  TObjString* testObjString2 = (TObjString*)arr->At(i+1);
564  fPtHardBin = testObjString2->GetString().Atoi();
565  i = arr->GetEntriesFast();
566  }
567  }
568  }
569 
570  if(!fEventCuts->GetDoEtaShift()) return kTRUE; // No Eta Shift requested, continue
571  if(fEventCuts->GetEtaShift() == 0.0){ // Eta Shift requested but not set, get shift automatically
572  fEventCuts->GetCorrectEtaShiftFromPeriod();
573  fEventCuts->DoEtaShift(kFALSE); // Eta Shift Set, make sure that it is called only once
574  return kTRUE;
575  } else {
576  printf(" Gamma Conversion Reader %s_%s :: Eta Shift Manually Set to %f \n\n",
577  (fEventCuts->GetCutNumber()).Data(),(fConversionCuts->GetCutNumber()).Data(),fEventCuts->GetEtaShift());
578  fEventCuts->DoEtaShift(kFALSE); // Eta Shift Set, make sure that it is called only once
579  }
580 
581  return kTRUE;
582 }
583 //________________________________________________________________________
585 
586  AliESDEvent * esdEvent = dynamic_cast<AliESDEvent*>(fInputEvent);
587  if(esdEvent) {
588  if (!TGeoGlobalMagField::Instance()->GetField()) esdEvent->InitMagneticField();
589  }
590 
591  // Check if correctly initialized
592  if(!fConversionGammas)Init();
593 
594  // User Exec
595  fEventIsSelected=ProcessEvent(fInputEvent,fMCEvent);
596 
597  // AOD Output
598  FillAODOutput();
599 
600 }
601 
602 //________________________________________________________________________
603 Bool_t AliV0ReaderV1::ProcessEvent(AliVEvent *inputEvent,AliMCEvent *mcEvent)
604 {
605 
606  //Reset the TClonesArray
607  fConversionGammas->Delete();
608 
609  //Clear TBits object with accepted v0s from previous event
610  if (kAddv0sInESDFilter){fPCMv0BitField->Clear();}
611 
612  fInputEvent = inputEvent;
613  fMCEvent = mcEvent;
614 
615  if(!fInputEvent){
616  AliError("No Input event");
617  return kFALSE;
618  }
619  if(!fEventCuts){AliError("No EventCuts");return kFALSE;}
620  if(!fConversionCuts){AliError("No ConversionCuts");return kFALSE;}
621 
622 
623  // Count Primary Tracks Event
624  CountTracks();
625 
626  //Count Tracks with TPCout flag
627  CountTPCoutTracks();
628 
629  // Event Cuts
630  if(!fEventCuts->EventIsSelected(fInputEvent,fMCEvent)){
631  if (fEventCuts->GetEventQuality() == 2 && !fMCFileChecked ){
632  cout << "ERROR with MC reading for: "<< fCurrentFileName.Data() << endl;
633  fMCFileChecked = kTRUE;
634  }
635  return kFALSE;
636  }
637  // Set Magnetic Field
638  AliKFParticle::SetField(fInputEvent->GetMagneticField());
639 
640  if(fInputEvent->IsA()==AliAODEvent::Class() && fProduceV0findingEffi){
641  fProduceV0findingEffi = kFALSE;
642  AliWarning("V0finding effi cannot be run on AODs ");
643  }
644 
645  if(fProduceV0findingEffi){
646  CreatePureMCHistosForV0FinderEffiESD();
647  fVectorFoundGammas.clear();
648  }
649 
650  if(fInputEvent->IsA()==AliESDEvent::Class()){
651  ProcessESDV0s();
652  }
653  if(fInputEvent->IsA()==AliAODEvent::Class()){
654  GetAODConversionGammas();
655  }
656 
657  return kTRUE;
658 }
661 {
662  // Fill AOD Output with reconstructed Photons
663 
664  if(fInputEvent->IsA()==AliESDEvent::Class()){
666  if(fCreateAOD) {
667  PostData(1, fPCMv0BitField);
668  AliAODHandler * aodhandler = dynamic_cast<AliAODHandler*>(AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler());
669  if (aodhandler && aodhandler->GetFillAOD()) {
670  AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler()->SetFillExtension(kTRUE);
671  //PostData(0, fConversionGammas);
672 
673  }
674  }
675  }
676 }
677 
679 const AliExternalTrackParam *AliV0ReaderV1::GetExternalTrackParam(AliESDv0 *fCurrentV0,Int_t &tracklabel,Int_t charge){
680 
681  // Get External Track Parameter with given charge
682 
683  if(!(charge==1||charge==-1)){AliError("Charge not defined");return 0x0;}
684 
685  if(fConversionCuts->GetV0FinderSameSign()==1){
686  if(fCurrentV0){
687  if((fConversionCuts->GetTrack(fInputEvent,fCurrentV0->GetPindex()))->Charge()!=(fConversionCuts->GetTrack(fInputEvent,fCurrentV0->GetNindex()))->Charge())return 0x0;
688  if((fConversionCuts->GetTrack(fInputEvent,fCurrentV0->GetPindex()))->Charge()==(fConversionCuts->GetTrack(fInputEvent,fCurrentV0->GetNindex()))->Charge()){
689  if(charge==1){
690  tracklabel=fCurrentV0->GetPindex();
691  return fCurrentV0->GetParamP();
692  }else{
693  tracklabel=fCurrentV0->GetNindex();
694  return fCurrentV0->GetParamN();
695  }
696  }
697  }
698  }else if(fConversionCuts->GetV0FinderSameSign()==2){
699  if(fCurrentV0){
700  if(charge==1){
701  tracklabel=fCurrentV0->GetPindex();
702  return fCurrentV0->GetParamP();
703  }else{
704  tracklabel=fCurrentV0->GetNindex();
705  return fCurrentV0->GetParamN();
706  }
707  }
708  }else{
709  // Check for sign flip
710  if(fCurrentV0){
711  if(!fCurrentV0->GetParamN()||!fCurrentV0->GetParamP())return 0x0;
712  if(!fConversionCuts->GetTrack(fInputEvent,fCurrentV0->GetNindex())||!fConversionCuts->GetTrack(fInputEvent,fCurrentV0->GetPindex()))return 0x0;
713  if((fConversionCuts->GetTrack(fInputEvent,fCurrentV0->GetPindex()))->Charge()==charge){
714  tracklabel=fCurrentV0->GetPindex();
715  return fCurrentV0->GetParamP();}
716  if((fConversionCuts->GetTrack(fInputEvent,fCurrentV0->GetNindex()))->Charge()==charge){
717  tracklabel=fCurrentV0->GetNindex();
718  return fCurrentV0->GetParamN();}
719  }
720  }
721  return 0x0;
722 }
723 
726 {
727 
728  // Process ESD V0s for conversion photon reconstruction
729  AliESDEvent *fESDEvent=dynamic_cast<AliESDEvent*>(fInputEvent);
730 
731  AliKFConversionPhoton *fCurrentMotherKFCandidate=NULL;
732 
733  if(fESDEvent){
734  for(Int_t currentV0Index=0;currentV0Index<fESDEvent->GetNumberOfV0s();currentV0Index++){
735  AliESDv0 *fCurrentV0=(AliESDv0*)(fESDEvent->GetV0(currentV0Index));
736  if(!fCurrentV0){
737  printf("Requested V0 does not exist");
738  continue;
739  }
740 
741  fCurrentMotherKFCandidate=ReconstructV0(fCurrentV0,currentV0Index);
742 
743  if(fCurrentMotherKFCandidate){
744  // Add Gamma to the TClonesArray
745 
746  if(kUseAODConversionPhoton){
747  new((*fConversionGammas)[fConversionGammas->GetEntriesFast()]) AliAODConversionPhoton(fCurrentMotherKFCandidate);
748  AliAODConversionPhoton * currentConversionPhoton = (AliAODConversionPhoton*)(fConversionGammas->At(fConversionGammas->GetEntriesFast()-1));
749  currentConversionPhoton->SetMass(fCurrentMotherKFCandidate->M());
750  if (fUseMassToZero) currentConversionPhoton->SetMassToZero();
751  currentConversionPhoton->SetInvMassPair(fCurrentInvMassPair);
752  if(kAddv0sInESDFilter){fPCMv0BitField->SetBitNumber(currentV0Index, kTRUE);}
753  } else {
754  new((*fConversionGammas)[fConversionGammas->GetEntriesFast()]) AliKFConversionPhoton(*fCurrentMotherKFCandidate);
755  }
756 
757  delete fCurrentMotherKFCandidate;
758  fCurrentMotherKFCandidate=NULL;
759  }
760  }
761  if(kAddv0sInESDFilter){fPCMv0BitField->Compact();}
762  }
763  return kTRUE;
764 }
765 
767 AliKFConversionPhoton *AliV0ReaderV1::ReconstructV0(AliESDv0 *fCurrentV0,Int_t currentV0Index)
768 {
769  // cout << currentV0Index << endl;
770  // Reconstruct conversion photon from ESD v0
771  fConversionCuts->FillPhotonCutIndex(AliConversionPhotonCuts::kPhotonIn);
772 
773  //checks if on the fly mode is set
774  if(!fConversionCuts->SelectV0Finder(fCurrentV0->GetOnFlyStatus())){
775  fConversionCuts->FillPhotonCutIndex(AliConversionPhotonCuts::kOnFly);
776  return 0x0;
777  }
778 
779  if (fMCEvent && fProduceV0findingEffi ) FillRecMCHistosForV0FinderEffiESD(fCurrentV0);
780 
781  // TrackLabels
782  Int_t currentTrackLabels[2]={-1,-1};
783 
784  // Get Daughter KF Particles
785 
786  const AliExternalTrackParam *fCurrentExternalTrackParamPositive=GetExternalTrackParamP(fCurrentV0,currentTrackLabels[0]);
787  // cout << fCurrentExternalTrackParamPositive << "\t" << currentTrackLabels[0] << endl;
788  const AliExternalTrackParam *fCurrentExternalTrackParamNegative=GetExternalTrackParamN(fCurrentV0,currentTrackLabels[1]);
789  // cout << fCurrentExternalTrackParamNegative << "\t" << currentTrackLabels[1] << endl;
790  if(!fCurrentExternalTrackParamPositive||!fCurrentExternalTrackParamNegative)return 0x0;
791 
792  // Apply some Cuts before Reconstruction
793 
794  AliVTrack * posTrack = fConversionCuts->GetTrack(fInputEvent,currentTrackLabels[0]);
795  AliVTrack * negTrack = fConversionCuts->GetTrack(fInputEvent,currentTrackLabels[1]);
796  if(!negTrack || !posTrack) {
797  fConversionCuts->FillPhotonCutIndex(AliConversionPhotonCuts::kNoTracks);
798  return 0x0;
799  }
800  // Track Cuts
801  if(!fConversionCuts->TracksAreSelected(negTrack, posTrack)){
802  fConversionCuts->FillPhotonCutIndex(AliConversionPhotonCuts::kTrackCuts);
803  return 0x0;
804  }
805 
806  fConversionCuts->FillV0EtaBeforedEdxCuts(fCurrentV0->Eta());
807  if (!fConversionCuts->dEdxCuts(posTrack)) {
808  fConversionCuts->FillPhotonCutIndex(AliConversionPhotonCuts::kdEdxCuts);
809  return 0x0;
810  }
811  // PID Cuts
812  if(!fConversionCuts->dEdxCuts(negTrack)) {
813  fConversionCuts->FillPhotonCutIndex(AliConversionPhotonCuts::kdEdxCuts);
814  return 0x0;
815  }
816  fConversionCuts->FillV0EtaAfterdEdxCuts(fCurrentV0->Eta());
817  // Reconstruct Photon
818  AliKFConversionPhoton *fCurrentMotherKF=NULL;
819  // fUseConstructGamma = kFALSE;
820  // cout << "construct gamma " << endl;
821  AliKFParticle fCurrentNegativeKFParticle(*(fCurrentExternalTrackParamNegative),11);
822  // cout << fCurrentExternalTrackParamNegative << "\t" << endl;
823  AliKFParticle fCurrentPositiveKFParticle(*(fCurrentExternalTrackParamPositive),-11);
824  // cout << fCurrentExternalTrackParamPositive << "\t" << endl;
825  // cout << currentTrackLabels[0] << "\t" << currentTrackLabels[1] << endl;
826  // cout << "construct gamma " <<fUseConstructGamma << endl;
827 
828  // Reconstruct Gamma
829  if(fUseConstructGamma){
830  fCurrentMotherKF = new AliKFConversionPhoton();
831  fCurrentMotherKF->ConstructGamma(fCurrentNegativeKFParticle,fCurrentPositiveKFParticle);
832  }else{
833  fCurrentMotherKF = new AliKFConversionPhoton(fCurrentNegativeKFParticle,fCurrentPositiveKFParticle);
834  fCurrentMotherKF->SetMassConstraint(0,0.0001);
835  }
836 
837  // Set Track Labels
838 
839  fCurrentMotherKF->SetTrackLabels(currentTrackLabels[0],currentTrackLabels[1]);
840 
841  // Set V0 index
842 
843  fCurrentMotherKF->SetV0Index(currentV0Index);
844 
845  //Set MC Label
846  if(fMCEvent){
847 
848  Int_t labelp=TMath::Abs(fConversionCuts->GetTrack(fInputEvent,fCurrentMotherKF->GetTrackLabelPositive())->GetLabel());
849  Int_t labeln=TMath::Abs(fConversionCuts->GetTrack(fInputEvent,fCurrentMotherKF->GetTrackLabelNegative())->GetLabel());
850 
851 // cout << "rec: " << currentTrackLabels[0] << "\t" << currentTrackLabels[1] << endl;
852 // cout << "recProp: " << fCurrentMotherKF->GetTrackLabelPositive() << "\t" << fCurrentMotherKF->GetTrackLabelNegative() << endl;
853 // cout << "MC: " << labeln << "\t" << labelp << endl;
854 
855  TParticle *fNegativeMCParticle = 0x0;
856  if(labeln>-1) fNegativeMCParticle = fMCEvent->Particle(labeln);
857  TParticle *fPositiveMCParticle = 0x0;
858  if(labelp>-1) fPositiveMCParticle = fMCEvent->Particle(labelp);
859 
860  if(fPositiveMCParticle&&fNegativeMCParticle){
861  fCurrentMotherKF->SetMCLabelPositive(labelp);
862  fCurrentMotherKF->SetMCLabelNegative(labeln);
863  }
864  }
865 
866  // Update Vertex (moved for same eta compared to old)
867  // cout << currentV0Index <<" \t before: \t" << fCurrentMotherKF->GetPx() << "\t" << fCurrentMotherKF->GetPy() << "\t" << fCurrentMotherKF->GetPz() << endl;
868  if(fUseImprovedVertex == kTRUE){
869  AliKFVertex primaryVertexImproved(*fInputEvent->GetPrimaryVertex());
870  // cout << "Prim Vtx: " << primaryVertexImproved.GetX() << "\t" << primaryVertexImproved.GetY() << "\t" << primaryVertexImproved.GetZ() << endl;
871  primaryVertexImproved+=*fCurrentMotherKF;
872  fCurrentMotherKF->SetProductionVertex(primaryVertexImproved);
873  }
874  // SetPsiPair
875  Double_t convpos[3]={0,0,0};
876  if (fImprovedPsiPair == 0){
877  Double_t PsiPair=GetPsiPair(fCurrentV0,fCurrentExternalTrackParamPositive,fCurrentExternalTrackParamNegative, convpos);
878  fCurrentMotherKF->SetPsiPair(PsiPair);
879  }
880 
881  // Recalculate ConversionPoint
882  Double_t dca[2]={0,0};
883  if(fUseOwnXYZCalculation){
884  // Double_t convpos[3]={0,0,0};
885  if(!GetConversionPoint(fCurrentExternalTrackParamPositive,fCurrentExternalTrackParamNegative,convpos,dca)){
886  fConversionCuts->FillPhotonCutIndex(AliConversionPhotonCuts::kConvPointFail);
887  delete fCurrentMotherKF;
888  fCurrentMotherKF=NULL;
889  return 0x0;
890  }
891 
892  fCurrentMotherKF->SetConversionPoint(convpos);
893  }
894 
895 
896  // SetPsiPair
897  if (fImprovedPsiPair >= 1){
898  // the propagation can be more precise after the precise conversion point calculation
899  Double_t PsiPair=GetPsiPair(fCurrentV0,fCurrentExternalTrackParamPositive,fCurrentExternalTrackParamNegative,convpos);
900  fCurrentMotherKF->SetPsiPair(PsiPair);
901  //cout<<" GetPsiPair::"<<fCurrentMotherKF->GetPsiPair() <<endl;
902  }
903 
904  if(fCurrentMotherKF->GetNDF() > 0.)
905  fCurrentMotherKF->SetChi2perNDF(fCurrentMotherKF->GetChi2()/fCurrentMotherKF->GetNDF()); //->Photon is created before all chi2 relevant changes are performed, set it "by hand"
906 
907 
908  // Set Dilepton Mass (moved down for same eta compared to old)
909  fCurrentMotherKF->SetMass(fCurrentMotherKF->M());
910 
911  // Calculating invariant mass
912  Double_t mass=-99.0, mass_width=-99.0, Pt=-99.0, Pt_width=-99.0;
913  AliKFParticle fCurrentMotherKFForMass(fCurrentNegativeKFParticle,fCurrentPositiveKFParticle);
914  fCurrentMotherKFForMass.GetMass(mass,mass_width);
915  fCurrentMotherKFForMass.GetPt(Pt,Pt_width);
916  fCurrentInvMassPair=mass;
917 
918  // apply possible Kappa cut
919  if (!fConversionCuts->KappaCuts(fCurrentMotherKF,fInputEvent)){
920  fConversionCuts->FillPhotonCutIndex(AliConversionPhotonCuts::kdEdxCuts);
921  delete fCurrentMotherKF;
922  fCurrentMotherKF=NULL;
923  return 0x0;
924  }
925 
926  // Apply Photon Cuts
927  if(!fConversionCuts->PhotonCuts(fCurrentMotherKF,fInputEvent)){
928  fConversionCuts->FillPhotonCutIndex(AliConversionPhotonCuts::kPhotonCuts);
929  delete fCurrentMotherKF;
930  fCurrentMotherKF=NULL;
931  return 0x0;
932  }
933 
934  // cout << currentV0Index <<" \t after: \t" <<fCurrentMotherKF->GetPx() << "\t" << fCurrentMotherKF->GetPy() << "\t" << fCurrentMotherKF->GetPz() << endl;
935 
936  if(fProduceImpactParamHistograms) FillImpactParamHistograms(posTrack, negTrack, fCurrentV0, fCurrentMotherKF);
937 
938  fConversionCuts->FillPhotonCutIndex(AliConversionPhotonCuts::kPhotonOut);
939  return fCurrentMotherKF;
940 }
941 
943 Double_t AliV0ReaderV1::GetPsiPair(const AliESDv0* v0, const AliExternalTrackParam *positiveparam,const AliExternalTrackParam *negativeparam,const Double_t convpos[3]) const {
944  //
945  // Angle between daughter momentum plane and plane
946  //
947 
948  AliExternalTrackParam nt(*negativeparam);
949  AliExternalTrackParam pt(*positiveparam);
950 
951  Float_t magField = fInputEvent->GetMagneticField();
952 
953  Double_t xyz[3] = {0.,0.,0.};
954  if (fImprovedPsiPair==0 ) {
955  v0->GetXYZ(xyz[0],xyz[1],xyz[2]);
956  } else if (fImprovedPsiPair>=1 ) {
957  xyz[0]= convpos[0];
958  xyz[1]= convpos[1];
959  xyz[2]= convpos[2];
960  }
961 
962 
963  // Double_t pPlus[3] = {pt.Px(),pt.Py(),pt.Pz()};
964  // Double_t pMinus[3] = {nt.Px(),nt.Py(),nt.Pz()};
965 
966  // Double_t u[3] = {pPlus[0]+pMinus[0],pPlus[1]+pMinus[1],pPlus[2]+pMinus[2]};
967  // Double_t normu = sqrt( (u[0]*u[0]) + (u[1]*u[1]) + (u[2]*u[2]) );
968 
969  // u[0] = u[0] / normu;
970  // u[1] = u[1] / normu;
971  // u[2] = u[2] / normu;
972 
973  // Double_t normpPlus = sqrt( (pPlus[0]*pPlus[0]) + (pPlus[1]*pPlus[1]) + (pPlus[2]*pPlus[2]) );
974  // Double_t normpMinus = sqrt( (pMinus[0]*pMinus[0]) + (pMinus[1]*pMinus[1]) + (pMinus[2]*pMinus[2]) );
975 
976  // pPlus[0] = pPlus[0] / normpPlus;
977  // pPlus[1] = pPlus[1] / normpPlus;
978  // pPlus[2] = pPlus[2] / normpPlus;
979 
980  // pMinus[0] = pMinus[0] / normpMinus;
981  // pMinus[1] = pMinus[1] / normpMinus;
982  // pMinus[2] = pMinus[2] / normpMinus;
983 
984  // Double_t v[3] = {0,0,0}; // pPlus X pMinus
985  // v[0] = (pPlus[1]*pMinus[2]) - (pPlus[2]*pMinus[1]);
986  // v[1] = (pPlus[2]*pMinus[0]) - (pPlus[0]*pMinus[2]);
987  // v[2] = (pPlus[0]*pMinus[1]) - (pPlus[1]*pMinus[0]);
988 
989  // Double_t w[3] = {0,0,0}; // u X v
990  // w[0] = (u[1]*v[2]) - (u[2]*v[1]);
991  // w[1] = (u[2]*v[0]) - (u[0]*v[2]);
992  // w[2] = (u[0]*v[1]) - (u[1]*v[0]);
993 
994  // Double_t z[3] = {0,0,1};
995  // Double_t wc[3] = {0,0,0}; // u X v
996  // wc[0] = (u[1]*z[2]) - (u[2]*z[1]);
997  // wc[1] = (u[2]*z[0]) - (u[0]*z[2]);
998  // wc[2] = (u[0]*z[1]) - (u[1]*z[0]);
999 
1000  // Double_t PhiV = TMath::ACos((w[0]*wc[0]) + (w[1]*wc[1]) + (w[2]*wc[2]));
1001  //return TMath::Abs(PhiV);
1002 
1003 
1004  // TVector3 pPlus(pt.Px(),pt.Py(),pt.Pz());
1005  // TVector3 pMinus(nt.Px(),nt.Py(),nt.Pz());
1006 
1007  // TVector3 u = pMinus + pPlus;
1008  // u = u*(1/u.Mag());
1009 
1010  // TVector3 pHPlus = pPlus*(1/pPlus.Mag());
1011  // TVector3 pHMinus = pMinus*(1/pMinus.Mag());
1012 
1013  // TVector3 v = pHPlus.Cross(pHMinus);
1014  // TVector3 w = u.Cross(v);
1015  // TVector3 z(0,0,1);
1016  // TVector3 wc = u.Cross(z);
1017 
1018  // Double_t PhiV = w * wc;
1019 
1020  Double_t mn[3] = {0,0,0};
1021  Double_t mp[3] = {0,0,0};
1022 
1023  v0->GetNPxPyPz(mn[0],mn[1],mn[2]);//reconstructed cartesian momentum components of negative daughter;
1024  v0->GetPPxPyPz(mp[0],mp[1],mp[2]);//reconstructed cartesian momentum components of positive daughter;
1025 
1026  Double_t deltat = 1.;
1027  deltat = TMath::ATan(mp[2]/(TMath::Sqrt(mp[0]*mp[0] + mp[1]*mp[1])+1.e-13)) - TMath::ATan(mn[2]/(TMath::Sqrt(mn[0]*mn[0] + mn[1]*mn[1])+1.e-13));//difference of angles of the two daughter tracks with z-axis
1028  Double_t radiussum = TMath::Sqrt(xyz[0]*xyz[0] + xyz[1]*xyz[1]) + 50;//radius to which tracks shall be propagated
1029 
1030  Double_t momPosProp[3] = {0,0,0};
1031  Double_t momNegProp[3] = {0,0,0};
1032 
1033  Double_t psiPair = 4.;
1034  // cout<< "Momentum before propagation at radius ::"<< TMath::Sqrt(xyz[0]*xyz[0] + xyz[1]*xyz[1])<< endl;
1035  // cout<< mn[0]<< " " <<mn[1]<<" "<< mn[2]<<" "<<TMath::Sqrt(mn[0]*mn[0] + mn[1]*mn[1] +mn[2]*mn[2])<< endl;
1036  // cout<< mp[0]<< " " <<mp[1]<<" "<< mp[2]<<" "<<TMath::Sqrt(mp[0]*mp[0] + mp[1]*mp[1] +mp[2]*mp[2])<< endl;
1037  Double_t pEle,pPos;
1038 
1039  if (fImprovedPsiPair==1 || fImprovedPsiPair==0 ){
1040  if(nt.PropagateTo(radiussum,magField) == 0) return psiPair; //propagate tracks to the outside -> Better Purity and Efficiency
1041  if(pt.PropagateTo(radiussum,magField) == 0) return psiPair; //propagate tracks to the outside -> Better Purity and Efficiency
1042 
1043  pt.GetPxPyPz(momPosProp);//Get momentum vectors of tracks after propagation
1044  nt.GetPxPyPz(momNegProp);
1045  } else if (fImprovedPsiPair>=2) {
1046  momPosProp[0] = pt.GetParameterAtRadius(radiussum,magField,3);
1047  momPosProp[1] = pt.GetParameterAtRadius(radiussum,magField,4);
1048  momPosProp[2] = pt.GetParameterAtRadius(radiussum,magField,5);
1049 
1050  momNegProp[0] = nt.GetParameterAtRadius(radiussum,magField,3);
1051  momNegProp[1] = nt.GetParameterAtRadius(radiussum,magField,4);
1052  momNegProp[2] = nt.GetParameterAtRadius(radiussum,magField,5);
1053  pEle = TMath::Sqrt(momNegProp[0]*momNegProp[0]+momNegProp[1]*momNegProp[1]+momNegProp[2]*momNegProp[2]);//absolute momentum value of negative daughter
1054 
1055  pPos = TMath::Sqrt(momPosProp[0]*momPosProp[0]+momPosProp[1]*momPosProp[1]+momPosProp[2]*momPosProp[2]);//absolute momentum value of positive daughter
1056 
1057  if ( (pEle==0 || pPos==0) && fImprovedPsiPair==3) {
1058  radiussum = TMath::Sqrt(xyz[0]*xyz[0] + xyz[1]*xyz[1]) + 30;
1059  momPosProp[0] = pt.GetParameterAtRadius(radiussum,magField,3);
1060  momPosProp[1] = pt.GetParameterAtRadius(radiussum,magField,4);
1061  momPosProp[2] = pt.GetParameterAtRadius(radiussum,magField,5);
1062 
1063  momNegProp[0] = nt.GetParameterAtRadius(radiussum,magField,3);
1064  momNegProp[1] = nt.GetParameterAtRadius(radiussum,magField,4);
1065  momNegProp[2] = nt.GetParameterAtRadius(radiussum,magField,5);
1066 
1067  }
1068  }
1069 
1070  pEle = TMath::Sqrt(momNegProp[0]*momNegProp[0]+momNegProp[1]*momNegProp[1]+momNegProp[2]*momNegProp[2]);//absolute momentum value of negative daughter
1071  pPos = TMath::Sqrt(momPosProp[0]*momPosProp[0]+momPosProp[1]*momPosProp[1]+momPosProp[2]*momPosProp[2]);//absolute momentum value of positive daughter
1072 
1073  Double_t scalarproduct =
1074  momPosProp[0]*momNegProp[0]+momPosProp[1]*momNegProp[1]+momPosProp[2]*momNegProp[2];//scalar product of propagated positive and negative daughters' momenta
1075 
1076  if (pEle==0 || pPos==0) return psiPair;
1077 
1078 
1079  Double_t chipair = TMath::ACos(scalarproduct/(pEle*pPos));//Angle between propagated daughter tracks
1080 
1081  // psiPair = TMath::Abs(TMath::ASin(deltat/chipair));
1082  psiPair = TMath::ASin(deltat/chipair);
1083  return psiPair;
1084 }
1085 
1087 Bool_t AliV0ReaderV1::GetHelixCenter(const AliExternalTrackParam *track,Double_t center[2]){
1088 
1089  // Get Center of the helix track parametrization
1090 
1091  Int_t charge=track->Charge();
1092  Double_t b=fInputEvent->GetMagneticField();
1093 
1094  Double_t helix[6];
1095  track->GetHelixParameters(helix,b);
1096 
1097  Double_t xpos = helix[5];
1098  Double_t ypos = helix[0];
1099  Double_t radius = TMath::Abs(1./helix[4]);
1100  Double_t phi = helix[2];
1101 
1102  if(phi < 0){
1103  phi = phi + 2*TMath::Pi();
1104  }
1105 
1106  phi -= TMath::Pi()/2.;
1107  Double_t xpoint = radius * TMath::Cos(phi);
1108  Double_t ypoint = radius * TMath::Sin(phi);
1109 
1110  if(b<0){
1111  if(charge > 0){
1112  xpoint = - xpoint;
1113  ypoint = - ypoint;
1114  }
1115  }
1116  if(b>0){
1117  if(charge < 0){
1118  xpoint = - xpoint;
1119  ypoint = - ypoint;
1120  }
1121  }
1122  center[0] = xpos + xpoint;
1123  center[1] = ypos + ypoint;
1124 
1125  return 1;
1126 }
1128 Bool_t AliV0ReaderV1::GetConversionPoint(const AliExternalTrackParam *pparam,const AliExternalTrackParam *nparam,Double_t convpos[3],Double_t dca[2]){
1129 
1130  // Recalculate Conversion Point
1131 
1132  if(!pparam||!nparam)return kFALSE;
1133 
1134  Double_t helixcenterpos[2];
1135  GetHelixCenter(pparam,helixcenterpos);
1136 
1137  Double_t helixcenterneg[2];
1138  GetHelixCenter(nparam,helixcenterneg);
1139 
1140  Double_t helixpos[6];
1141  pparam->GetHelixParameters(helixpos,fInputEvent->GetMagneticField());
1142  Double_t posradius = TMath::Abs(1./helixpos[4]);
1143 
1144  Double_t helixneg[6];
1145  nparam->GetHelixParameters(helixneg,fInputEvent->GetMagneticField());
1146  Double_t negradius = TMath::Abs(1./helixneg[4]);
1147 
1148  // Calculate xy-position
1149 
1150  Double_t xpos = helixcenterpos[0];
1151  Double_t ypos = helixcenterpos[1];
1152  Double_t xneg = helixcenterneg[0];
1153  Double_t yneg = helixcenterneg[1];
1154 
1155  convpos[0] = (xpos*negradius + xneg*posradius)/(negradius+posradius);
1156  convpos[1] = (ypos*negradius+ yneg*posradius)/(negradius+posradius);
1157 
1158 
1159  // Calculate Track XY vertex position
1160 
1161  Double_t deltaXPos = convpos[0] - xpos;
1162  Double_t deltaYPos = convpos[1] - ypos;
1163 
1164  Double_t deltaXNeg = convpos[0] - xneg;
1165  Double_t deltaYNeg = convpos[1] - yneg;
1166 
1167  Double_t alphaPos = TMath::Pi() + TMath::ATan2(-deltaYPos,-deltaXPos);
1168  Double_t alphaNeg = TMath::Pi() + TMath::ATan2(-deltaYNeg,-deltaXNeg);
1169 
1170  Double_t vertexXNeg = xneg + TMath::Abs(negradius)*TMath::Cos(alphaNeg);
1171  Double_t vertexYNeg = yneg + TMath::Abs(negradius)*TMath::Sin(alphaNeg);
1172 
1173  Double_t vertexXPos = xpos + TMath::Abs(posradius)*TMath::Cos(alphaPos);
1174  Double_t vertexYPos = ypos + TMath::Abs(posradius)*TMath::Sin(alphaPos);
1175 
1176  AliExternalTrackParam p(*pparam);
1177  AliExternalTrackParam n(*nparam);
1178 
1179  TVector2 vertexPos(vertexXPos,vertexYPos);
1180  TVector2 vertexNeg(vertexXNeg,vertexYNeg);
1181 
1182  // Convert to local coordinate system
1183  TVector2 vertexPosRot=vertexPos.Rotate(-p.GetAlpha());
1184  TVector2 vertexNegRot=vertexNeg.Rotate(-n.GetAlpha());
1185 
1186  // Propagate Track Params to Vertex
1187 
1188  if(!p.PropagateTo(vertexPosRot.X(),fInputEvent->GetMagneticField()))return kFALSE;
1189  if(!n.PropagateTo(vertexNegRot.X(),fInputEvent->GetMagneticField()))return kFALSE;
1190 
1191  // Check whether propagation was sucessful
1192 
1193  if(TMath::Abs(vertexPos.Mod()-TMath::Sqrt(p.GetX()*p.GetX()+p.GetY()*p.GetY()))>0.01)return kFALSE;
1194  if(TMath::Abs(vertexNeg.Mod()-TMath::Sqrt(n.GetX()*n.GetX()+n.GetY()*n.GetY()))>0.01)return kFALSE;
1195 
1196  // Calculate z position
1197 
1198  convpos[2] = (p.GetZ()*negradius+n.GetZ()*posradius)/(negradius+posradius);
1199 
1200  // Calculate DCA
1201  TVector2 vdca=vertexPos-vertexNeg;
1202  dca[0]=vdca.Mod();
1203  dca[1]=TMath::Abs(n.GetZ()-p.GetZ());
1204 
1205  return kTRUE;
1206 }
1207 
1208 //________________________________________________________________________
1210 
1211  // Get reconstructed conversion photons from satellite AOD file
1212 
1213  AliAODEvent *fAODEvent=dynamic_cast<AliAODEvent*>(fInputEvent);
1214 
1215  if(fAODEvent){
1216 
1217  if(fConversionGammas == NULL){
1218  fConversionGammas = new TClonesArray("AliAODConversionPhoton",100);
1219  }
1220  fConversionGammas->Delete();//Reset the TClonesArray
1221 
1222  //Get Gammas from satellite AOD gamma branch
1223 
1224  AliAODConversionPhoton *gamma=0x0;
1225 
1226  TClonesArray *fInputGammas=dynamic_cast<TClonesArray*>(fAODEvent->FindListObject(fDeltaAODBranchName.Data()));
1227 
1228  if(!fInputGammas){
1229  FindDeltaAODBranchName();
1230  fInputGammas=dynamic_cast<TClonesArray*>(fAODEvent->FindListObject(fDeltaAODBranchName.Data()));}
1231  if(!fInputGammas){AliError("No Gamma Satellites found");return kFALSE;}
1232  // Apply Selection Cuts to Gammas and create local working copy
1233  if(fInputGammas){
1234  for(Int_t i=0;i<fInputGammas->GetEntriesFast();i++){
1235  gamma=dynamic_cast<AliAODConversionPhoton*>(fInputGammas->At(i));
1236  if(gamma){
1237  if(fRelabelAODs)RelabelAODPhotonCandidates(gamma);
1238  if(fConversionCuts->PhotonIsSelected(gamma,fInputEvent)){
1239  new((*fConversionGammas)[fConversionGammas->GetEntriesFast()]) AliAODConversionPhoton(*gamma);}
1240  }
1241  }
1242  }
1243  }
1244 
1245  if(fConversionGammas->GetEntries()){return kTRUE;}
1246 
1247  return kFALSE;
1248 }
1249 
1250 //________________________________________________________________________
1252 
1253  // Find delta AOD branchname containing reconstructed photons
1254 
1255  TList *list=fInputEvent->GetList();
1256  for(Int_t ii=0;ii<list->GetEntries();ii++){
1257  TString name((list->At(ii))->GetName());
1258  if(name.BeginsWith(fDeltaAODBranchName)&&name.EndsWith("gamma")){
1259  fDeltaAODBranchName=name;
1260  AliInfo(Form("Set DeltaAOD BranchName to: %s",fDeltaAODBranchName.Data()));
1261  return;
1262  }
1263  }
1264 }
1265 
1266 //________________________________________________________________________
1268 
1269  if(fPreviousV0ReaderPerformsAODRelabeling == 2) return;
1270  else if(fPreviousV0ReaderPerformsAODRelabeling == 0){
1271  printf("Running AODs! Determine if V0Reader '%s' should perform relabeling\n",this->GetName());
1272  TObjArray* obj = (TObjArray*)AliAnalysisManager::GetAnalysisManager()->GetTasks();
1273  Int_t iPosition = obj->IndexOf(this);
1274  Bool_t prevV0ReaderRunningButNotRelabeling = kFALSE;
1275  for(Int_t i=iPosition-1; i>=0; i--){
1276  if( (obj->At(i))->IsA() == AliV0ReaderV1::Class()){
1277  AliV0ReaderV1* tempReader = (AliV0ReaderV1*) obj->At(i);
1278  if( tempReader->AreAODsRelabeled() && tempReader->IsReaderPerformingRelabeling() == 1){
1279  fPreviousV0ReaderPerformsAODRelabeling = 2;
1280  prevV0ReaderRunningButNotRelabeling = kFALSE;
1281  printf("V0Reader '%s' is running before this V0Reader '%s': do _NOT_ relabel AODs by current reader!\n",tempReader->GetName(),this->GetName());
1282  break;
1283  }else prevV0ReaderRunningButNotRelabeling = kTRUE;
1284  }
1285  }
1286  if(prevV0ReaderRunningButNotRelabeling) AliFatal(Form("There are V0Readers before '%s', but none of them is relabeling!",this->GetName()));
1287 
1288  if(fPreviousV0ReaderPerformsAODRelabeling == 2) return;
1289  else{
1290  printf("This V0Reader '%s' is first to be processed: do relabel AODs by current reader!\n",this->GetName());
1291  fPreviousV0ReaderPerformsAODRelabeling = 1;
1292  }
1293  }
1294 
1295  if(fPreviousV0ReaderPerformsAODRelabeling != 1) AliFatal(Form("In %s: fPreviousV0ReaderPerformsAODRelabeling = '%i' - while it should be impossible it is something different than '1'!",this->GetName(),fPreviousV0ReaderPerformsAODRelabeling));
1296 
1297  // Relabeling For AOD Event
1298  // ESDiD -> AODiD
1299  // MCLabel -> AODMCLabel
1300  Bool_t AODLabelPos = kFALSE;
1301  Bool_t AODLabelNeg = kFALSE;
1302 
1303  for(Int_t i = 0; i<fInputEvent->GetNumberOfTracks();i++){
1304  AliAODTrack *tempDaughter = static_cast<AliAODTrack*>(fInputEvent->GetTrack(i));
1305  if(!AODLabelPos){
1306  if( tempDaughter->GetID() == PhotonCandidate->GetTrackLabelPositive() ){
1307  PhotonCandidate->SetMCLabelPositive(TMath::Abs(tempDaughter->GetLabel()));
1308  PhotonCandidate->SetLabelPositive(i);
1309  AODLabelPos = kTRUE;
1310  }
1311  }
1312  if(!AODLabelNeg){
1313  if( tempDaughter->GetID() == PhotonCandidate->GetTrackLabelNegative()){
1314  PhotonCandidate->SetMCLabelNegative(TMath::Abs(tempDaughter->GetLabel()));
1315  PhotonCandidate->SetLabelNegative(i);
1316  AODLabelNeg = kTRUE;
1317  }
1318  }
1319  if(AODLabelNeg && AODLabelPos){
1320  return;
1321  }
1322  }
1323  if(!AODLabelPos || !AODLabelNeg){
1324  AliError(Form("NO AOD Daughters Found Pos: %i %i Neg: %i %i, setting all labels to -999999",AODLabelPos,PhotonCandidate->GetTrackLabelPositive(),AODLabelNeg,PhotonCandidate->GetTrackLabelNegative()));
1325  if(!AODLabelNeg){
1326  PhotonCandidate->SetMCLabelNegative(-999999);
1327  PhotonCandidate->SetLabelNegative(-999999);
1328  }
1329  if(!AODLabelPos){
1330  PhotonCandidate->SetMCLabelPositive(-999999);
1331  PhotonCandidate->SetLabelPositive(-999999);
1332  }
1333  }
1334 
1335 }
1336 
1337 //************************************************************************
1338 // This function counts the number of primary tracks in the event, the
1339 // implementation for ESD and AOD had to be different due to the different
1340 // filters which are already applied on AOD level the so-called filter
1341 // bits, we tried to replicate the conditions in both but an exact match
1342 // could not be reached. The cuts are similar to the ones used by the
1343 // jet-group
1344 //************************************************************************
1345 //________________________________________________________________________
1347 
1348  // In Case of MC count only MB tracks
1349  // if(fMCEvent && fInputEvent->IsA()==AliESDEvent::Class()){
1350  // fEventCuts->GetNotRejectedParticles(1,NULL,fMCEvent);
1351  // }
1352  // if(fMCEvent && fInputEvent->IsA()==AliAODEvent::Class()){
1353  // fEventCuts->GetNotRejectedParticles(1,NULL,fInputEvent);
1354  // }
1355 
1356  if(fInputEvent->IsA()==AliESDEvent::Class()){
1357  static AliESDtrackCuts *EsdTrackCuts = 0x0;
1358  static int prevRun = -1;
1359  // Using standard function for setting Cuts
1360  Int_t runNumber = fInputEvent->GetRunNumber();
1361 
1362  if (prevRun!=runNumber) {
1363  delete EsdTrackCuts;
1364  EsdTrackCuts = 0;
1365  prevRun = runNumber;
1366  }
1367  if (!EsdTrackCuts) {
1368  // if LHC11a or earlier or if LHC13g or if LHC12a-i -> use 2010 cuts
1369  if( (runNumber<=146860) || (runNumber>=197470 && runNumber<=197692) || (runNumber>=172440 && runNumber<=193766) ){
1370  EsdTrackCuts = AliESDtrackCuts::GetStandardITSTPCTrackCuts2010();
1371  // else if run2 data use 2015 PbPb cuts
1372  } else if (runNumber>=209122){
1373  //EsdTrackCuts = AliESDtrackCuts::GetStandardITSTPCTrackCuts2015PbPb();
1374  // hard coded track cuts for the moment, because AliESDtrackCuts::GetStandardITSTPCTrackCuts2015PbPb() gives spams warnings
1375  EsdTrackCuts = new AliESDtrackCuts();
1376  // TPC; clusterCut = 1, cutAcceptanceEdges = kTRUE, removeDistortedRegions = kFALSE
1377  EsdTrackCuts->AliESDtrackCuts::SetMinNCrossedRowsTPC(70);
1378  EsdTrackCuts->AliESDtrackCuts::SetMinRatioCrossedRowsOverFindableClustersTPC(0.8);
1379  EsdTrackCuts->SetCutGeoNcrNcl(2., 130., 1.5, 0.0, 0.0); // only dead zone and not clusters per length
1380  //EsdTrackCuts->AliESDtrackCuts::SetCutOutDistortedRegionsTPC(kTRUE);
1381  EsdTrackCuts->AliESDtrackCuts::SetMaxChi2PerClusterTPC(4);
1382  EsdTrackCuts->AliESDtrackCuts::SetAcceptKinkDaughters(kFALSE);
1383  EsdTrackCuts->AliESDtrackCuts::SetRequireTPCRefit(kTRUE);
1384  // ITS; selPrimaries = 1
1385  EsdTrackCuts->AliESDtrackCuts::SetRequireITSRefit(kTRUE);
1386  EsdTrackCuts->AliESDtrackCuts::SetClusterRequirementITS(AliESDtrackCuts::kSPD,
1387  AliESDtrackCuts::kAny);
1388  EsdTrackCuts->AliESDtrackCuts::SetMaxDCAToVertexXYPtDep("0.0105+0.0350/pt^1.1");
1389  EsdTrackCuts->AliESDtrackCuts::SetMaxChi2TPCConstrainedGlobal(36);
1390  EsdTrackCuts->AliESDtrackCuts::SetMaxDCAToVertexZ(2);
1391  EsdTrackCuts->AliESDtrackCuts::SetDCAToVertex2D(kFALSE);
1392  EsdTrackCuts->AliESDtrackCuts::SetRequireSigmaToVertex(kFALSE);
1393  EsdTrackCuts->AliESDtrackCuts::SetMaxChi2PerClusterITS(36);
1394  // else use 2011 version of track cuts
1395  }else{
1396  EsdTrackCuts = AliESDtrackCuts::GetStandardITSTPCTrackCuts2011();
1397  }
1398  EsdTrackCuts->SetMaxDCAToVertexZ(2);
1399  EsdTrackCuts->SetEtaRange(-0.8, 0.8);
1400  EsdTrackCuts->SetPtRange(0.15);
1401  }
1402  fNumberOfPrimaryTracks = 0;
1403  for(Int_t iTracks = 0; iTracks < fInputEvent->GetNumberOfTracks(); iTracks++){
1404  AliESDtrack* curTrack = (AliESDtrack*) fInputEvent->GetTrack(iTracks);
1405  if(!curTrack) continue;
1406  if(!EsdTrackCuts->AcceptTrack(curTrack)) continue;
1407  //if(fMCEvent && !(fEventCuts->IsParticleFromBGEvent(TMath::Abs(curTrack->GetLabel()),fMCEvent->Stack(),fInputEvent))) continue;
1408  fNumberOfPrimaryTracks++;
1409  }
1410  }
1411  else if(fInputEvent->IsA()==AliAODEvent::Class()){
1412  fNumberOfPrimaryTracks = 0;
1413  for(Int_t iTracks = 0; iTracks<fInputEvent->GetNumberOfTracks(); iTracks++){
1414  AliAODTrack* curTrack = (AliAODTrack*) fInputEvent->GetTrack(iTracks);
1415  if(curTrack->GetID()<0) continue; // Avoid double counting of tracks
1416  if(!curTrack->IsHybridGlobalConstrainedGlobal()) continue;
1417  if(TMath::Abs(curTrack->Eta())>0.8) continue;
1418  if(curTrack->Pt()<0.15) continue;
1419  //if(fMCEvent && !(fEventCuts->IsParticleFromBGEvent(TMath::Abs(curTrack->GetLabel()),NULL,fInputEvent))) continue;
1420  //if(TMath::Abs(curTrack->ZAtDCA())>2) continue; // Only Set For TPCOnly tracks
1421  fNumberOfPrimaryTracks++;
1422  }
1423  }
1424 
1425  return;
1426 }
1427 
1430  fNumberOfTPCoutTracks = 0;
1431 
1432  for (Int_t itrk = 0; itrk < fInputEvent->GetNumberOfTracks(); itrk++) {
1433  AliVTrack *trk = dynamic_cast<AliVTrack*>(fInputEvent->GetTrack(itrk));
1434 
1435  if (trk != NULL) {
1436  /* the initial method of counting TPC out tracks */
1437  if (!(trk->Pt() < 0.15) && (TMath::Abs(trk->Eta()) < 0.8)) {
1438  if ((trk->GetStatus() & AliVTrack::kTPCout) == AliVTrack::kTPCout) {
1439  fNumberOfTPCoutTracks++;
1440  }
1441  }
1442  }
1443  }
1444  return;
1445 }
1446 
1447 
1449 Bool_t AliV0ReaderV1::ParticleIsConvertedPhoton(AliMCEvent *mcEvent, TParticle *particle, Double_t etaMax, Double_t rMax, Double_t zMax){
1450  // MonteCarlo Photon Selection
1451  if(!mcEvent)return kFALSE;
1452 
1453  if (particle->GetPdgCode() == 22){
1454  // check whether particle is within eta range
1455  if( TMath::Abs(particle->Eta()) > etaMax ) return kFALSE;
1456  // check if particle doesn't have a photon as mother
1457  if(particle->GetMother(0) >-1 && mcEvent->Particle(particle->GetMother(0))->GetPdgCode() == 22){
1458  return kFALSE; // no photon as mothers!
1459  }
1460  // looking for conversion gammas (electron + positron from pairbuilding (= 5) )
1461  TParticle* ePos = NULL;
1462  TParticle* eNeg = NULL;
1463  if(particle->GetNDaughters() >= 2){
1464  for(Int_t daughterIndex=particle->GetFirstDaughter();daughterIndex<=particle->GetLastDaughter();daughterIndex++){
1465  if(daughterIndex<0) continue;
1466  TParticle *tmpDaughter = mcEvent->Particle(daughterIndex);
1467  if(tmpDaughter->GetUniqueID() == 5){
1468  if(tmpDaughter->GetPdgCode() == 11){
1469  eNeg = tmpDaughter;
1470  } else if(tmpDaughter->GetPdgCode() == -11){
1471  ePos = tmpDaughter;
1472  }
1473  }
1474  }
1475  }
1476  if(ePos == NULL || eNeg == NULL){ // means we do not have two daughters from pair production
1477  return kFALSE;
1478  }
1479  // check if electrons are in correct eta window
1480  if( TMath::Abs(ePos->Eta()) > etaMax ||
1481  TMath::Abs(eNeg->Eta()) > etaMax )
1482  return kFALSE;
1483 
1484  // check if photons have converted in reconstructable range
1485  if(ePos->R() > rMax){
1486  return kFALSE; // cuts on distance from collision point
1487  }
1488  if(TMath::Abs(ePos->Vz()) > zMax){
1489  return kFALSE; // outside material
1490  }
1491  if(TMath::Abs(eNeg->Vz()) > zMax){
1492  return kFALSE; // outside material
1493  }
1494 
1495 
1496  Double_t lineCutZRSlope = tan(2*atan(exp(-etaMax)));
1497  Double_t lineCutZValue = 7.;
1498  if( ePos->R() <= ((TMath::Abs(ePos->Vz()) * lineCutZRSlope) - lineCutZValue)){
1499  return kFALSE; // line cut to exclude regions where we do not reconstruct
1500  }
1501  if( eNeg->R() <= ((TMath::Abs(eNeg->Vz()) * lineCutZRSlope) - lineCutZValue)){
1502  return kFALSE; // line cut to exclude regions where we do not reconstruct
1503  }
1504  if (ePos->Pt() < 0.05 || eNeg->Pt() < 0.05){
1505  return kFALSE;
1506  }
1507 
1508  return kTRUE;
1509  }
1510  return kFALSE;
1511 }
1512 
1515  const AliVVertex* primVtxMC = fMCEvent->GetPrimaryVertex();
1516  Double_t mcProdVtxX = primVtxMC->GetX();
1517  Double_t mcProdVtxY = primVtxMC->GetY();
1518  Double_t mcProdVtxZ = primVtxMC->GetZ();
1519  //cout << mcProdVtxX <<"\t" << mcProdVtxY << "\t" << mcProdVtxZ << endl;
1520 
1521  // Loop over all primary MC particle
1522  for(Long_t i = 0; i < fMCEvent->GetNumberOfTracks(); i++) {
1523  if (fEventCuts->IsConversionPrimaryESD( fMCEvent, i, mcProdVtxX, mcProdVtxY, mcProdVtxZ)){
1524  // fill primary histogram
1525  TParticle* particle = (TParticle *)fMCEvent->Particle(i);
1526  if (!particle) continue;
1527  if (ParticleIsConvertedPhoton(fMCEvent, particle, 0.9, 180.,250. )){
1528  if(particle->GetFirstDaughter()<0) continue;
1529  TParticle *tmpDaughter = fMCEvent->Particle(particle->GetFirstDaughter());
1530  if (!tmpDaughter) continue;
1531  fHistoMCGammaPtvsR->Fill(particle->Pt(),tmpDaughter->R());
1532  fHistoMCGammaPtvsPhi->Fill(particle->Pt(),particle->Phi());
1533  fHistoMCGammaRvsPhi->Fill(tmpDaughter->R(),particle->Phi());
1534  }
1535  if (ParticleIsConvertedPhoton(fMCEvent, particle, 1.4, 180.,250. )){
1536  if(particle->GetFirstDaughter()<0) continue;
1537  TParticle *tmpDaughter = fMCEvent->Particle(particle->GetFirstDaughter());
1538  if (!tmpDaughter) continue;
1539  fHistoMCGammaPtvsEta->Fill(particle->Pt(),particle->Eta());
1540  fHistoMCGammaRvsEta->Fill(tmpDaughter->R(),particle->Eta());
1541  fHistoMCGammaPhivsEta->Fill(particle->Phi(),particle->Eta());
1542  }
1543  }
1544  }
1545 }
1546 
1549 
1550  const AliVVertex* primVtxMC = fMCEvent->GetPrimaryVertex();
1551  Double_t mcProdVtxX = primVtxMC->GetX();
1552  Double_t mcProdVtxY = primVtxMC->GetY();
1553  Double_t mcProdVtxZ = primVtxMC->GetZ();
1554 // cout << mcProdVtxX <<"\t" << mcProdVtxY << "\t" << mcProdVtxZ << endl;
1555 
1556  Int_t tracklabelPos=currentV0->GetPindex();
1557  Int_t tracklabelNeg=currentV0->GetNindex();
1558 
1559  Int_t labelp=TMath::Abs(fConversionCuts->GetTrack(fInputEvent,tracklabelPos)->GetLabel());
1560  Int_t labeln=TMath::Abs(fConversionCuts->GetTrack(fInputEvent,tracklabelNeg)->GetLabel());
1561 
1562  TParticle* negPart = 0x0;
1563  if(labeln>-1) negPart = (TParticle *)fMCEvent->Particle(labeln);
1564  TParticle* posPart = 0x0;
1565  if(labelp>-1) posPart = (TParticle *)fMCEvent->Particle(labelp);
1566 
1567  if ( negPart == NULL || posPart == NULL ) return;
1568 // if (!(negPart->GetPdgCode() == 11)) return;
1569 // if (!(posPart->GetPdgCode() == -11)) return;
1570  Long_t motherlabelNeg = negPart->GetFirstMother();
1571  Long_t motherlabelPos = posPart->GetFirstMother();
1572 
1573 // cout << "mother neg " << motherlabelNeg << " mother pos " << motherlabelPos << endl;
1574  if (motherlabelNeg>-1 && motherlabelNeg == motherlabelPos && negPart->GetFirstMother() != -1){
1575  if (fEventCuts->IsConversionPrimaryESD( fMCEvent, negPart->GetFirstMother(), mcProdVtxX, mcProdVtxY, mcProdVtxZ)){
1576 
1577  TParticle* mother = (TParticle *)fMCEvent->Particle(motherlabelNeg);
1578  if (mother->GetPdgCode() == 22 ){
1579  if (!CheckVectorForDoubleCount(fVectorFoundGammas,motherlabelNeg ) ){
1580  if (ParticleIsConvertedPhoton(fMCEvent, mother, 0.9, 180.,250. )){
1581  fHistoRecMCGammaPtvsR->Fill(mother->Pt(),negPart->R());
1582  fHistoRecMCGammaPtvsPhi->Fill(mother->Pt(),mother->Phi());
1583  fHistoRecMCGammaRvsPhi->Fill(negPart->R(),mother->Phi());
1584  }
1585  if (ParticleIsConvertedPhoton(fMCEvent, mother, 1.4, 180.,250. )){
1586  fHistoRecMCGammaPtvsEta->Fill(mother->Pt(),mother->Eta());
1587  fHistoRecMCGammaRvsEta->Fill(negPart->R(),mother->Eta());
1588  fHistoRecMCGammaPhivsEta->Fill(mother->Phi(),mother->Eta());
1589  }
1590 // cout << "new gamma found" << endl;
1591  } else {
1592  if (ParticleIsConvertedPhoton(fMCEvent, mother, 0.9, 180.,250. )){
1593  fHistoRecMCGammaMultiPt->Fill(mother->Pt());
1594  fHistoRecMCGammaMultiPhi->Fill(mother->Phi());
1595  fHistoRecMCGammaMultiR->Fill(negPart->R());
1596  }
1597  if (ParticleIsConvertedPhoton(fMCEvent, mother, 1.4, 180.,250. )){
1598  fHistoRecMCGammaMultiPtvsEta->Fill(mother->Pt(),mother->Eta());
1599  }
1600 // cout << "this one I had already: " << motherlabelNeg << endl << "-----------------------" << endl;
1601  }
1602 // cout << "event gammas: " << endl;
1603 // for(Int_t iGamma=0; iGamma<fVectorFoundGammas.size(); iGamma++){cout << fVectorFoundGammas.at(iGamma) << ", ";}
1604  }
1605  }
1606  }
1607 }
1608 
1609 //_________________________________________________________________________________
1610 void AliV0ReaderV1::FillImpactParamHistograms( AliVTrack* pTrack, AliVTrack* nTrack, AliESDv0 *fCurrentV0, AliKFConversionPhoton *fCurrentMotherKF){
1611 
1612  // values of cuts to be introduced in 2016 to reduce ESD size
1613  Float_t fZmax = 25; //cm
1614  Float_t fYmax = 25; //cm
1615  Double_t kTPCMargin = 1.0; //cm
1616  Int_t kMaxTPCV0Conflicts = 1; //# conflicting clusters tolerated
1617 
1618  //conversion point
1619  Double_t convX, convY, convZ;
1620  fCurrentV0->GetXYZ(convX,convY,convZ);
1621  Double_t convR = TMath::Sqrt(convX*convX+convY*convY);
1622  //recalculated conversion point
1623  Double_t convXrecalc=fCurrentMotherKF->GetConversionX();
1624  Double_t convYrecalc=fCurrentMotherKF->GetConversionY();
1625  Double_t convZrecalc=fCurrentMotherKF->GetConversionZ();
1626  Double_t convRrecalc = fCurrentMotherKF->GetConversionRadius();
1627 
1628  //count V0s
1629  fHistoImpactParameterStudy->Fill(0);
1630 
1631  //count V0s with two TPC-only tracks
1632  AliESDtrack* positiveTrack = (AliESDtrack*) pTrack;
1633  AliESDtrack* negativeTrack = (AliESDtrack*) nTrack;
1634  Bool_t TPConly = kFALSE;
1635  if (!positiveTrack->IsOn(AliESDtrack::kITSin) && !negativeTrack->IsOn(AliESDtrack::kITSin)){
1636  fHistoImpactParameterStudy->Fill(1);
1637  TPConly = kTRUE;
1638  }
1639 
1640  Bool_t RemovedByZcut = kFALSE;
1641  Bool_t RemovedByYcut = kFALSE;
1642 
1643  //count V0s which have...
1644  if(TPConly){
1645  //not passed z cut: pos.tr. or neg.tr.
1646  if(((TMath::Abs(positiveTrack->GetZ()))>fZmax) || ((TMath::Abs(negativeTrack->GetZ()))>fZmax)){
1647  fHistoImpactParameterStudy->Fill(2);
1648  RemovedByZcut=kTRUE;
1649  }
1650 
1651  //not passed y cut: pos.tr. or neg.tr.
1652  if(((TMath::Abs(positiveTrack->GetY()))>fYmax) || ((TMath::Abs(negativeTrack->GetY()))>fYmax)){
1653  fHistoImpactParameterStudy->Fill(3);
1654  RemovedByYcut=kTRUE;
1655  }
1656  }
1657 
1658  //causality cut
1659  Bool_t RemovedByCausality=kFALSE;
1660  AliESDEvent *fESDEvent=dynamic_cast<AliESDEvent*>(fInputEvent);
1661  const Float_t rTPC[159]={
1662  85.225, 85.975, 86.725, 87.475, 88.225, 88.975, 89.725, 90.475, 91.225, 91.975, 92.725, 93.475, 94.225, 94.975, 95.725,
1663  96.475, 97.225, 97.975, 98.725, 99.475,100.225,100.975,101.725,102.475,103.225,103.975,104.725,105.475,106.225,106.975,
1664  107.725,108.475,109.225,109.975,110.725,111.475,112.225,112.975,113.725,114.475,115.225,115.975,116.725,117.475,118.225,
1665  118.975,119.725,120.475,121.225,121.975,122.725,123.475,124.225,124.975,125.725,126.475,127.225,127.975,128.725,129.475,
1666  130.225,130.975,131.725,135.100,136.100,137.100,138.100,139.100,140.100,141.100,142.100,143.100,144.100,145.100,146.100,
1667  147.100,148.100,149.100,150.100,151.100,152.100,153.100,154.100,155.100,156.100,157.100,158.100,159.100,160.100,161.100,
1668  162.100,163.100,164.100,165.100,166.100,167.100,168.100,169.100,170.100,171.100,172.100,173.100,174.100,175.100,176.100,
1669  177.100,178.100,179.100,180.100,181.100,182.100,183.100,184.100,185.100,186.100,187.100,188.100,189.100,190.100,191.100,
1670  192.100,193.100,194.100,195.100,196.100,197.100,198.100,199.350,200.850,202.350,203.850,205.350,206.850,208.350,209.850,
1671  211.350,212.850,214.350,215.850,217.350,218.850,220.350,221.850,223.350,224.850,226.350,227.850,229.350,230.850,232.350,
1672  233.850,235.350,236.850,238.350,239.850,241.350,242.850,244.350,245.850};
1673 
1674  // fill conversion radius histograms
1675  fHistoR->Fill(convR);
1676  fHistoRrecalc->Fill(convRrecalc);
1677  Double_t alpha = TMath::ATan2(convY,convX);
1678  if (alpha<0) alpha += TMath::Pi()*2;
1679  Int_t sec = alpha/(TMath::Pi()/9);
1680  alpha = (10.+sec*20.)*TMath::DegToRad();
1681  Double_t cs = TMath::Cos(alpha);
1682  Double_t sn = TMath::Sin(alpha);
1683  Double_t xsV0 = convX*cs - convY*sn;
1684  fHistoRviaAlpha->Fill(xsV0);
1685  Double_t alpha_r = TMath::ATan2(convYrecalc,convXrecalc);
1686  if (alpha_r<0) alpha_r += TMath::Pi()*2;
1687  Int_t sec_r = alpha_r/(TMath::Pi()/9);
1688  alpha_r = (10.+sec_r*20.)*TMath::DegToRad();
1689  Double_t cs_r = TMath::Cos(alpha_r);
1690  Double_t sn_r = TMath::Sin(alpha_r);
1691  Double_t xsV0_r = convXrecalc*cs_r - convYrecalc*sn_r;
1692  fHistoRviaAlphaRecalc->Fill(xsV0_r);
1693 
1694  if (convR > 80) { // conversion within TPC <-> TPC-only tracks
1695  fHistoImpactParameterStudy->Fill(4);
1696 
1697  for (Int_t it=2;it--;) {
1698  Int_t trId = fCurrentV0->GetIndex(it);
1699  AliESDtrack* tr = fESDEvent->GetTrack(trId);
1700  const TBits& bits = tr->GetTPCClusterMap();
1701  Int_t nConflict = 0;
1702  for (Int_t ic=0;ic<159;ic++) {
1703  if (rTPC[ic]>(xsV0-kTPCMargin)) break;
1704  if (bits.TestBitNumber(ic)){
1705  nConflict++;
1706  fHistoRdiff->Fill(xsV0-rTPC[ic]);
1707  }
1708  if (nConflict>kMaxTPCV0Conflicts) {
1709  fHistoImpactParameterStudy->Fill(5);
1710  RemovedByCausality=kTRUE;
1711  break;
1712  }
1713  }
1714  }
1715  }
1716 
1717  //not passed y or z o causality cut:
1718  Bool_t RemovedByAnyCut=kFALSE;
1719  if(RemovedByZcut||RemovedByYcut||RemovedByCausality){
1720  fHistoImpactParameterStudy->Fill(6);
1721  RemovedByAnyCut=kTRUE;
1722  }
1723 
1724  //Fill tree for further analysis
1725  Float_t posZ;
1726  Float_t posY;
1727  Float_t posX;
1728  Float_t posPt;
1729  Float_t negZ;
1730  Float_t negY;
1731  Float_t negX;
1732  Float_t negPt;
1733  Float_t R;
1734  fImpactParamTree->Branch("posPt",&posPt,"posPt/F");
1735  fImpactParamTree->Branch("posY",&posY,"posY/F");
1736  fImpactParamTree->Branch("R",&R,"R/F");
1737  posZ = positiveTrack->GetZ();
1738  posY = positiveTrack->GetY();
1739  posX = positiveTrack->GetX();
1740  posPt = positiveTrack->Pt();
1741  negZ = negativeTrack->GetZ();
1742  negY = negativeTrack->GetY();
1743  negX = negativeTrack->GetX();
1744  negPt = negativeTrack->Pt();
1745  R=convRrecalc;
1746  fImpactParamTree->Fill();
1747 
1748  // fill impact parameter histograms
1749  fHistoPosTrackImpactParamZ->Fill(posZ);
1750  fHistoPosTrackImpactParamY->Fill(posY);
1751  fHistoPosTrackImpactParamX->Fill(posX);
1752  fHistoPosTrackImpactParamZvsPt->Fill(posPt, posZ);
1753  fHistoPosTrackImpactParamYvsPt->Fill(posPt, posY);
1754  fHistoPosTrackImpactParamXvsPt->Fill(posPt, posX);
1755  fHistoNegTrackImpactParamZ->Fill(negZ);
1756  fHistoNegTrackImpactParamY->Fill(negY);
1757  fHistoNegTrackImpactParamX->Fill(negX);
1758  fHistoNegTrackImpactParamZvsPt->Fill(negPt, negZ);
1759  fHistoNegTrackImpactParamYvsPt->Fill(negPt, negY);
1760  fHistoNegTrackImpactParamXvsPt->Fill(negPt, negX);
1761 
1762  // fill comparison histos before / after new cuts
1763  fHistoPt->Fill(positiveTrack->Pt());
1764  fHistoImpactParamZvsR->Fill(convZrecalc,convRrecalc);
1765  //Float_t *DCAzPhoton;
1766  Float_t DCAzPhoton;
1767  //const AliVVertex *PrimVertex = fInputEvent->GetPrimaryVertex();
1768  AliAODConversionPhoton* fCurrentMother=(AliAODConversionPhoton*)fCurrentMotherKF;
1769  //AliAODConversionPhoton *fCurrentMother=dynamic_cast<AliAODConversionPhoton*>(fCurrentMotherKF);
1770  //fCurrentMotherKF->GetDistanceOfClossetApproachToPrimVtx(PrimVertex, DCAzPhoton);
1771  DCAzPhoton = fCurrentMother->GetDCAzToPrimVtx();
1772  fHistoDCAzPhoton->Fill(DCAzPhoton);
1773  if(!RemovedByAnyCut) {
1774  fHistoPt2->Fill(positiveTrack->Pt());
1775  fHistoImpactParamZvsR2->Fill(convZrecalc,convRrecalc);
1776  fHistoDCAzPhoton2->Fill(DCAzPhoton);
1777  }
1778 
1779  return;
1780 }
1781 
1782 //_________________________________________________________________________________
1783 Bool_t AliV0ReaderV1::CheckVectorOnly(vector<Int_t> &vec, Int_t tobechecked)
1784 {
1785  if(tobechecked > -1)
1786  {
1787  vector<Int_t>::iterator it;
1788  it = find (vec.begin(), vec.end(), tobechecked);
1789  if (it != vec.end()) return true;
1790  else return false;
1791  }
1792  return false;
1793 }
1794 
1795 //_________________________________________________________________________________
1797 {
1798  if(tobechecked > -1)
1799  {
1800  vector<Int_t>::iterator it;
1801  it = find (vec.begin(), vec.end(), tobechecked);
1802  if (it != vec.end()) return true;
1803  else{
1804  vec.push_back(tobechecked);
1805  return false;
1806  }
1807  }
1808  return false;
1809 }
1810 
1811 //________________________________________________________________________
1813 {
1814 
1815 }
1816 
1818  std::iterator<std::bidirectional_iterator_tag, AliConversionPhotonBase>(),
1819  fkData(reader),
1820  fCurrentIndex(position),
1821  fDirection(dir)
1822 {
1823 }
1824 
1826  std::iterator<std::bidirectional_iterator_tag, AliConversionPhotonBase>(other),
1827  fkData(other.fkData),
1829  fDirection(other.fDirection)
1830 {
1831 }
1832 
1834  if(this != &other){
1835  fkData = other.fkData;
1836  fCurrentIndex = other.fCurrentIndex;
1837  fDirection = other.fDirection;
1838  }
1839  return *this;
1840 }
1841 
1843  if(fkData != other.fkData) return true;
1844  return fCurrentIndex != other.fCurrentIndex;
1845 }
1846 
1848  iterator tmp(*this);
1849  operator++();
1850  return tmp;
1851 }
1852 
1855  fCurrentIndex++;
1856  } else {
1857  fCurrentIndex--;
1858  }
1859  return *this;
1860 }
1861 
1863  iterator tmp(*this);
1864  operator--();
1865  return tmp;
1866 }
1867 
1870  fCurrentIndex++;
1871  } else {
1872  fCurrentIndex--;
1873  }
1874  return *this;
1875 }
1876 
1878  return (*fkData)[fCurrentIndex];
1879 }
Int_t IsReaderPerformingRelabeling()
Int_t charge
Direction_t fDirection
Iterator in forward direction.
Definition: AliV0ReaderV1.h:69
void FillAODOutput()
double Double_t
Definition: External.C:58
Bool_t ParticleIsConvertedPhoton(AliMCEvent *mcEvent, TParticle *particle, Double_t etaMax, Double_t rMax, Double_t zMax)
anchored LHC13[d-e] pass 2 - JJ
Definition: External.C:236
void FindDeltaAODBranchName()
void SetConversionPoint(Double_t convpoint[3])
Bool_t CheckVectorForDoubleCount(vector< Int_t > &vec, Int_t tobechecked)
void SetMCLabelNegative(Int_t label)
Bool_t AreAODsRelabeled()
const AliV0ReaderV1 * fkData
V0 reader used to iterate over.
Definition: AliV0ReaderV1.h:67
void CreatePureMCHistosForV0FinderEffiESD()
Double_t mass
AliConversionPhotonBase * operator*()
Float_t GetDCAzToPrimVtx() const
virtual void Init()
Double_t GetPsiPair(const AliESDv0 *v0, const AliExternalTrackParam *positiveparam, const AliExternalTrackParam *negativeparam, const Double_t convpos[3]) const
void FillRecMCHistosForV0FinderEffiESD(AliESDv0 *currentV0)
Bool_t CheckVectorOnly(vector< Int_t > &vec, Int_t tobechecked)
anchored LHC16f pass 1 low B-field - Pythia8 JJ
anchored LHC13g pass 1 - JJ
anchored LHC13g pass 1 - JJ
int fCurrentIndex
Index of the current element.
Definition: AliV0ReaderV1.h:68
anchored LHC17p/q pass 1 - jet-jet MC w/GEANT3,
int Int_t
Definition: External.C:63
virtual ~AliV0ReaderV1()
Bool_t ProcessEvent(AliVEvent *inputEvent, AliMCEvent *mcEvent=NULL)
float Float_t
Definition: External.C:68
void SetMCLabelPositive(Int_t label)
void ConstructGamma(const AliKFParticle &fCurrentNegativeKFParticle, const AliKFParticle &fCurrentPositiveKFParticle)
void SetInvMassPair(Float_t mass)
iterator(const AliV0ReaderV1 *reader, Direction_t dir, int position)
Bool_t GetHelixCenter(const AliExternalTrackParam *track, Double_t center[2])
void RelabelAODPhotonCandidates(AliAODConversionPhoton *PhotonCandidate)
void SetChi2perNDF(Float_t chi2)
void CountTPCoutTracks()
anchored LHC13[d-e] pass 2 - GJ
Bool_t Data(TH1F *h, Double_t *rangefit, Bool_t writefit, Double_t &sgn, Double_t &errsgn, Double_t &bkg, Double_t &errbkg, Double_t &sgnf, Double_t &errsgnf, Double_t &sigmafit, Int_t &status)
anchored LHC12[a-h] pass 2 - JJ
virtual Bool_t Notify()
iterator & operator=(const iterator &other)
anchored LHC13[d-e] pass 2 - JJ
void SetTrackLabels(Int_t label1, Int_t label2)
anchored LHC11a pass 4 - JJ
virtual void Terminate(Option_t *)
const AliExternalTrackParam * GetExternalTrackParam(AliESDv0 *fCurrentV0, Int_t &tracklabel, Int_t charge)
anchored LHC13g pass 1 - JJ
AliKFConversionPhoton * ReconstructV0(AliESDv0 *fCurrentV0, Int_t currentV0Index)
void SetLabelPositive(Int_t label)
Track labels.
void SetLabelNegative(Int_t label)
bool operator!=(iterator &other) const
TFile * file
TList with histograms for a given trigger.
void SetPsiPair(Float_t PsiPair)
void FillImpactParamHistograms(AliVTrack *ptrack, AliVTrack *ntrack, AliESDv0 *fCurrentV0, AliKFConversionPhoton *fCurrentMotherKF)
AliConversionPhotonBase * operator[](int index) const
Array index operator.
const char Option_t
Definition: External.C:48
anchored LHC11a pass 4 - JJ
Bool_t ProcessESDV0s()
bool Bool_t
Definition: External.C:53
anchored LHC16x pass 1 nom B-field - Pythia8 JJ
Bool_t GetConversionPoint(const AliExternalTrackParam *pparam, const AliExternalTrackParam *nparam, Double_t convpos[3], Double_t dca[2])
Bool_t GetAODConversionGammas()
void UserCreateOutputObjects()
virtual void UserExec(Option_t *option)
Double_t GetConversionRadius() const
TDirectoryFile * dir