AliPhysics  e6c8d43 (e6c8d43)
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::kLHC16rP1JJ || fEventCuts->GetPeriodEnum() == AliConvEventCuts::kLHC16sP1JJ ||
556  fEventCuts->GetPeriodEnum() == AliConvEventCuts::kLHC16P1JJ || fEventCuts->GetPeriodEnum() == AliConvEventCuts::kLHC16P1JJLowB ||
557  fEventCuts->GetPeriodEnum() == AliConvEventCuts::kLHC18b8
558  ){
559  TObjArray *arr = fCurrentFileName.Tokenize("/");
560  fPtHardBin = -1;
561  for (Int_t i = 0; i < arr->GetEntriesFast();i++ ){
562  TObjString* testObjString = (TObjString*)arr->At(i);
563  if (testObjString->GetString().BeginsWith("LHC")){
564  TObjString* testObjString2 = (TObjString*)arr->At(i+1);
565  fPtHardBin = testObjString2->GetString().Atoi();
566  i = arr->GetEntriesFast();
567  }
568  }
569  }
570 
571  if(!fEventCuts->GetDoEtaShift()) return kTRUE; // No Eta Shift requested, continue
572  if(fEventCuts->GetEtaShift() == 0.0){ // Eta Shift requested but not set, get shift automatically
573  fEventCuts->GetCorrectEtaShiftFromPeriod();
574  fEventCuts->DoEtaShift(kFALSE); // Eta Shift Set, make sure that it is called only once
575  return kTRUE;
576  } else {
577  printf(" Gamma Conversion Reader %s_%s :: Eta Shift Manually Set to %f \n\n",
578  (fEventCuts->GetCutNumber()).Data(),(fConversionCuts->GetCutNumber()).Data(),fEventCuts->GetEtaShift());
579  fEventCuts->DoEtaShift(kFALSE); // Eta Shift Set, make sure that it is called only once
580  }
581 
582  return kTRUE;
583 }
584 //________________________________________________________________________
586 
587  AliESDEvent * esdEvent = dynamic_cast<AliESDEvent*>(fInputEvent);
588  if(esdEvent) {
589  if (!TGeoGlobalMagField::Instance()->GetField()) esdEvent->InitMagneticField();
590  }
591 
592  // Check if correctly initialized
593  if(!fConversionGammas)Init();
594 
595  // User Exec
596  fEventIsSelected=ProcessEvent(fInputEvent,fMCEvent);
597 
598  // AOD Output
599  FillAODOutput();
600 
601 }
602 
603 //________________________________________________________________________
604 Bool_t AliV0ReaderV1::ProcessEvent(AliVEvent *inputEvent,AliMCEvent *mcEvent)
605 {
606 
607  //Reset the TClonesArray
608  fConversionGammas->Delete();
609 
610  //Clear TBits object with accepted v0s from previous event
611  if (kAddv0sInESDFilter){fPCMv0BitField->Clear();}
612 
613  fInputEvent = inputEvent;
614  fMCEvent = mcEvent;
615 
616  if(!fInputEvent){
617  AliError("No Input event");
618  return kFALSE;
619  }
620  if(!fEventCuts){AliError("No EventCuts");return kFALSE;}
621  if(!fConversionCuts){AliError("No ConversionCuts");return kFALSE;}
622 
623 
624  // Count Primary Tracks Event
625  CountTracks();
626 
627  //Count Tracks with TPCout flag
628  CountTPCoutTracks();
629 
630  // Event Cuts
631  if(!fEventCuts->EventIsSelected(fInputEvent,fMCEvent)){
632  if (fEventCuts->GetEventQuality() == 2 && !fMCFileChecked ){
633  cout << "ERROR with MC reading for: "<< fCurrentFileName.Data() << endl;
634  fMCFileChecked = kTRUE;
635  }
636  return kFALSE;
637  }
638  // Set Magnetic Field
639  AliKFParticle::SetField(fInputEvent->GetMagneticField());
640 
641  if(fInputEvent->IsA()==AliAODEvent::Class() && fProduceV0findingEffi){
642  fProduceV0findingEffi = kFALSE;
643  AliWarning("V0finding effi cannot be run on AODs ");
644  }
645 
646  if(fProduceV0findingEffi){
647  CreatePureMCHistosForV0FinderEffiESD();
648  fVectorFoundGammas.clear();
649  }
650 
651  if(fInputEvent->IsA()==AliESDEvent::Class()){
652  ProcessESDV0s();
653  }
654  if(fInputEvent->IsA()==AliAODEvent::Class()){
655  GetAODConversionGammas();
656  }
657 
658  return kTRUE;
659 }
662 {
663  // Fill AOD Output with reconstructed Photons
664 
665  if(fInputEvent->IsA()==AliESDEvent::Class()){
667  if(fCreateAOD) {
668  PostData(1, fPCMv0BitField);
669  AliAODHandler * aodhandler = dynamic_cast<AliAODHandler*>(AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler());
670  if (aodhandler && aodhandler->GetFillAOD()) {
671  AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler()->SetFillExtension(kTRUE);
672  //PostData(0, fConversionGammas);
673 
674  }
675  }
676  }
677 }
678 
680 const AliExternalTrackParam *AliV0ReaderV1::GetExternalTrackParam(AliESDv0 *fCurrentV0,Int_t &tracklabel,Int_t charge){
681 
682  // Get External Track Parameter with given charge
683 
684  if(!(charge==1||charge==-1)){AliError("Charge not defined");return 0x0;}
685 
686  if(fConversionCuts->GetV0FinderSameSign()==1){
687  if(fCurrentV0){
688  if((fConversionCuts->GetTrack(fInputEvent,fCurrentV0->GetPindex()))->Charge()!=(fConversionCuts->GetTrack(fInputEvent,fCurrentV0->GetNindex()))->Charge())return 0x0;
689  if((fConversionCuts->GetTrack(fInputEvent,fCurrentV0->GetPindex()))->Charge()==(fConversionCuts->GetTrack(fInputEvent,fCurrentV0->GetNindex()))->Charge()){
690  if(charge==1){
691  tracklabel=fCurrentV0->GetPindex();
692  return fCurrentV0->GetParamP();
693  }else{
694  tracklabel=fCurrentV0->GetNindex();
695  return fCurrentV0->GetParamN();
696  }
697  }
698  }
699  }else if(fConversionCuts->GetV0FinderSameSign()==2){
700  if(fCurrentV0){
701  if(charge==1){
702  tracklabel=fCurrentV0->GetPindex();
703  return fCurrentV0->GetParamP();
704  }else{
705  tracklabel=fCurrentV0->GetNindex();
706  return fCurrentV0->GetParamN();
707  }
708  }
709  }else{
710  // Check for sign flip
711  if(fCurrentV0){
712  if(!fCurrentV0->GetParamN()||!fCurrentV0->GetParamP())return 0x0;
713  if(!fConversionCuts->GetTrack(fInputEvent,fCurrentV0->GetNindex())||!fConversionCuts->GetTrack(fInputEvent,fCurrentV0->GetPindex()))return 0x0;
714  if((fConversionCuts->GetTrack(fInputEvent,fCurrentV0->GetPindex()))->Charge()==charge){
715  tracklabel=fCurrentV0->GetPindex();
716  return fCurrentV0->GetParamP();}
717  if((fConversionCuts->GetTrack(fInputEvent,fCurrentV0->GetNindex()))->Charge()==charge){
718  tracklabel=fCurrentV0->GetNindex();
719  return fCurrentV0->GetParamN();}
720  }
721  }
722  return 0x0;
723 }
724 
727 {
728 
729  // Process ESD V0s for conversion photon reconstruction
730  AliESDEvent *fESDEvent=dynamic_cast<AliESDEvent*>(fInputEvent);
731 
732  AliKFConversionPhoton *fCurrentMotherKFCandidate=NULL;
733 
734  if(fESDEvent){
735  for(Int_t currentV0Index=0;currentV0Index<fESDEvent->GetNumberOfV0s();currentV0Index++){
736  AliESDv0 *fCurrentV0=(AliESDv0*)(fESDEvent->GetV0(currentV0Index));
737  if(!fCurrentV0){
738  printf("Requested V0 does not exist");
739  continue;
740  }
741 
742  fCurrentMotherKFCandidate=ReconstructV0(fCurrentV0,currentV0Index);
743 
744  if(fCurrentMotherKFCandidate){
745  // Add Gamma to the TClonesArray
746 
747  if(kUseAODConversionPhoton){
748  new((*fConversionGammas)[fConversionGammas->GetEntriesFast()]) AliAODConversionPhoton(fCurrentMotherKFCandidate);
749  AliAODConversionPhoton * currentConversionPhoton = (AliAODConversionPhoton*)(fConversionGammas->At(fConversionGammas->GetEntriesFast()-1));
750  currentConversionPhoton->SetMass(fCurrentMotherKFCandidate->M());
751  if (fUseMassToZero) currentConversionPhoton->SetMassToZero();
752  currentConversionPhoton->SetInvMassPair(fCurrentInvMassPair);
753  if(kAddv0sInESDFilter){fPCMv0BitField->SetBitNumber(currentV0Index, kTRUE);}
754  } else {
755  new((*fConversionGammas)[fConversionGammas->GetEntriesFast()]) AliKFConversionPhoton(*fCurrentMotherKFCandidate);
756  }
757 
758  delete fCurrentMotherKFCandidate;
759  fCurrentMotherKFCandidate=NULL;
760  }
761  }
762  if(kAddv0sInESDFilter){fPCMv0BitField->Compact();}
763  }
764  return kTRUE;
765 }
766 
768 AliKFConversionPhoton *AliV0ReaderV1::ReconstructV0(AliESDv0 *fCurrentV0,Int_t currentV0Index)
769 {
770  // cout << currentV0Index << endl;
771  // Reconstruct conversion photon from ESD v0
772  fConversionCuts->FillPhotonCutIndex(AliConversionPhotonCuts::kPhotonIn);
773 
774  //checks if on the fly mode is set
775  if(!fConversionCuts->SelectV0Finder(fCurrentV0->GetOnFlyStatus())){
776  fConversionCuts->FillPhotonCutIndex(AliConversionPhotonCuts::kOnFly);
777  return 0x0;
778  }
779 
780  if (fMCEvent && fProduceV0findingEffi ) FillRecMCHistosForV0FinderEffiESD(fCurrentV0);
781 
782  // TrackLabels
783  Int_t currentTrackLabels[2]={-1,-1};
784 
785  // Get Daughter KF Particles
786 
787  const AliExternalTrackParam *fCurrentExternalTrackParamPositive=GetExternalTrackParamP(fCurrentV0,currentTrackLabels[0]);
788  // cout << fCurrentExternalTrackParamPositive << "\t" << currentTrackLabels[0] << endl;
789  const AliExternalTrackParam *fCurrentExternalTrackParamNegative=GetExternalTrackParamN(fCurrentV0,currentTrackLabels[1]);
790  // cout << fCurrentExternalTrackParamNegative << "\t" << currentTrackLabels[1] << endl;
791  if(!fCurrentExternalTrackParamPositive||!fCurrentExternalTrackParamNegative)return 0x0;
792 
793  // Apply some Cuts before Reconstruction
794 
795  AliVTrack * posTrack = fConversionCuts->GetTrack(fInputEvent,currentTrackLabels[0]);
796  AliVTrack * negTrack = fConversionCuts->GetTrack(fInputEvent,currentTrackLabels[1]);
797  if(!negTrack || !posTrack) {
798  fConversionCuts->FillPhotonCutIndex(AliConversionPhotonCuts::kNoTracks);
799  return 0x0;
800  }
801  // Track Cuts
802  if(!fConversionCuts->TracksAreSelected(negTrack, posTrack)){
803  fConversionCuts->FillPhotonCutIndex(AliConversionPhotonCuts::kTrackCuts);
804  return 0x0;
805  }
806 
807  fConversionCuts->FillV0EtaBeforedEdxCuts(fCurrentV0->Eta());
808  if (!fConversionCuts->dEdxCuts(posTrack)) {
809  fConversionCuts->FillPhotonCutIndex(AliConversionPhotonCuts::kdEdxCuts);
810  return 0x0;
811  }
812  // PID Cuts
813  if(!fConversionCuts->dEdxCuts(negTrack)) {
814  fConversionCuts->FillPhotonCutIndex(AliConversionPhotonCuts::kdEdxCuts);
815  return 0x0;
816  }
817  fConversionCuts->FillV0EtaAfterdEdxCuts(fCurrentV0->Eta());
818  // Reconstruct Photon
819  AliKFConversionPhoton *fCurrentMotherKF=NULL;
820  // fUseConstructGamma = kFALSE;
821  // cout << "construct gamma " << endl;
822  AliKFParticle fCurrentNegativeKFParticle(*(fCurrentExternalTrackParamNegative),11);
823  // cout << fCurrentExternalTrackParamNegative << "\t" << endl;
824  AliKFParticle fCurrentPositiveKFParticle(*(fCurrentExternalTrackParamPositive),-11);
825  // cout << fCurrentExternalTrackParamPositive << "\t" << endl;
826  // cout << currentTrackLabels[0] << "\t" << currentTrackLabels[1] << endl;
827  // cout << "construct gamma " <<fUseConstructGamma << endl;
828 
829  // Reconstruct Gamma
830  if(fUseConstructGamma){
831  fCurrentMotherKF = new AliKFConversionPhoton();
832  fCurrentMotherKF->ConstructGamma(fCurrentNegativeKFParticle,fCurrentPositiveKFParticle);
833  }else{
834  fCurrentMotherKF = new AliKFConversionPhoton(fCurrentNegativeKFParticle,fCurrentPositiveKFParticle);
835  fCurrentMotherKF->SetMassConstraint(0,0.0001);
836  }
837 
838  // Set Track Labels
839 
840  fCurrentMotherKF->SetTrackLabels(currentTrackLabels[0],currentTrackLabels[1]);
841 
842  // Set V0 index
843 
844  fCurrentMotherKF->SetV0Index(currentV0Index);
845 
846  //Set MC Label
847  if(fMCEvent){
848 
849  Int_t labelp=TMath::Abs(fConversionCuts->GetTrack(fInputEvent,fCurrentMotherKF->GetTrackLabelPositive())->GetLabel());
850  Int_t labeln=TMath::Abs(fConversionCuts->GetTrack(fInputEvent,fCurrentMotherKF->GetTrackLabelNegative())->GetLabel());
851 
852 // cout << "rec: " << currentTrackLabels[0] << "\t" << currentTrackLabels[1] << endl;
853 // cout << "recProp: " << fCurrentMotherKF->GetTrackLabelPositive() << "\t" << fCurrentMotherKF->GetTrackLabelNegative() << endl;
854 // cout << "MC: " << labeln << "\t" << labelp << endl;
855 
856  TParticle *fNegativeMCParticle = 0x0;
857  if(labeln>-1) fNegativeMCParticle = fMCEvent->Particle(labeln);
858  TParticle *fPositiveMCParticle = 0x0;
859  if(labelp>-1) fPositiveMCParticle = fMCEvent->Particle(labelp);
860 
861  if(fPositiveMCParticle&&fNegativeMCParticle){
862  fCurrentMotherKF->SetMCLabelPositive(labelp);
863  fCurrentMotherKF->SetMCLabelNegative(labeln);
864  }
865  }
866 
867  // Update Vertex (moved for same eta compared to old)
868  // cout << currentV0Index <<" \t before: \t" << fCurrentMotherKF->GetPx() << "\t" << fCurrentMotherKF->GetPy() << "\t" << fCurrentMotherKF->GetPz() << endl;
869  if(fUseImprovedVertex == kTRUE){
870  AliKFVertex primaryVertexImproved(*fInputEvent->GetPrimaryVertex());
871  // cout << "Prim Vtx: " << primaryVertexImproved.GetX() << "\t" << primaryVertexImproved.GetY() << "\t" << primaryVertexImproved.GetZ() << endl;
872  primaryVertexImproved+=*fCurrentMotherKF;
873  fCurrentMotherKF->SetProductionVertex(primaryVertexImproved);
874  }
875  // SetPsiPair
876  Double_t convpos[3]={0,0,0};
877  if (fImprovedPsiPair == 0){
878  Double_t PsiPair=GetPsiPair(fCurrentV0,fCurrentExternalTrackParamPositive,fCurrentExternalTrackParamNegative, convpos);
879  fCurrentMotherKF->SetPsiPair(PsiPair);
880  }
881 
882  // Recalculate ConversionPoint
883  Double_t dca[2]={0,0};
884  if(fUseOwnXYZCalculation){
885  // Double_t convpos[3]={0,0,0};
886  if(!GetConversionPoint(fCurrentExternalTrackParamPositive,fCurrentExternalTrackParamNegative,convpos,dca)){
887  fConversionCuts->FillPhotonCutIndex(AliConversionPhotonCuts::kConvPointFail);
888  delete fCurrentMotherKF;
889  fCurrentMotherKF=NULL;
890  return 0x0;
891  }
892 
893  fCurrentMotherKF->SetConversionPoint(convpos);
894  }
895 
896 
897  // SetPsiPair
898  if (fImprovedPsiPair >= 1){
899  // the propagation can be more precise after the precise conversion point calculation
900  Double_t PsiPair=GetPsiPair(fCurrentV0,fCurrentExternalTrackParamPositive,fCurrentExternalTrackParamNegative,convpos);
901  fCurrentMotherKF->SetPsiPair(PsiPair);
902  //cout<<" GetPsiPair::"<<fCurrentMotherKF->GetPsiPair() <<endl;
903  }
904 
905  if(fCurrentMotherKF->GetNDF() > 0.)
906  fCurrentMotherKF->SetChi2perNDF(fCurrentMotherKF->GetChi2()/fCurrentMotherKF->GetNDF()); //->Photon is created before all chi2 relevant changes are performed, set it "by hand"
907 
908 
909  // Set Dilepton Mass (moved down for same eta compared to old)
910  fCurrentMotherKF->SetMass(fCurrentMotherKF->M());
911 
912  // Calculating invariant mass
913  Double_t mass=-99.0, mass_width=-99.0, Pt=-99.0, Pt_width=-99.0;
914  AliKFParticle fCurrentMotherKFForMass(fCurrentNegativeKFParticle,fCurrentPositiveKFParticle);
915  fCurrentMotherKFForMass.GetMass(mass,mass_width);
916  fCurrentMotherKFForMass.GetPt(Pt,Pt_width);
917  fCurrentInvMassPair=mass;
918 
919  // apply possible Kappa cut
920  if (!fConversionCuts->KappaCuts(fCurrentMotherKF,fInputEvent)){
921  fConversionCuts->FillPhotonCutIndex(AliConversionPhotonCuts::kdEdxCuts);
922  delete fCurrentMotherKF;
923  fCurrentMotherKF=NULL;
924  return 0x0;
925  }
926 
927  // Apply Photon Cuts
928  if(!fConversionCuts->PhotonCuts(fCurrentMotherKF,fInputEvent)){
929  fConversionCuts->FillPhotonCutIndex(AliConversionPhotonCuts::kPhotonCuts);
930  delete fCurrentMotherKF;
931  fCurrentMotherKF=NULL;
932  return 0x0;
933  }
934 
935  // cout << currentV0Index <<" \t after: \t" <<fCurrentMotherKF->GetPx() << "\t" << fCurrentMotherKF->GetPy() << "\t" << fCurrentMotherKF->GetPz() << endl;
936 
937  if(fProduceImpactParamHistograms) FillImpactParamHistograms(posTrack, negTrack, fCurrentV0, fCurrentMotherKF);
938 
939  fConversionCuts->FillPhotonCutIndex(AliConversionPhotonCuts::kPhotonOut);
940  return fCurrentMotherKF;
941 }
942 
944 Double_t AliV0ReaderV1::GetPsiPair(const AliESDv0* v0, const AliExternalTrackParam *positiveparam,const AliExternalTrackParam *negativeparam,const Double_t convpos[3]) const {
945  //
946  // Angle between daughter momentum plane and plane
947  //
948 
949  AliExternalTrackParam nt(*negativeparam);
950  AliExternalTrackParam pt(*positiveparam);
951 
952  Float_t magField = fInputEvent->GetMagneticField();
953 
954  Double_t xyz[3] = {0.,0.,0.};
955  if (fImprovedPsiPair==0 ) {
956  v0->GetXYZ(xyz[0],xyz[1],xyz[2]);
957  } else if (fImprovedPsiPair>=1 ) {
958  xyz[0]= convpos[0];
959  xyz[1]= convpos[1];
960  xyz[2]= convpos[2];
961  }
962 
963 
964  // Double_t pPlus[3] = {pt.Px(),pt.Py(),pt.Pz()};
965  // Double_t pMinus[3] = {nt.Px(),nt.Py(),nt.Pz()};
966 
967  // Double_t u[3] = {pPlus[0]+pMinus[0],pPlus[1]+pMinus[1],pPlus[2]+pMinus[2]};
968  // Double_t normu = sqrt( (u[0]*u[0]) + (u[1]*u[1]) + (u[2]*u[2]) );
969 
970  // u[0] = u[0] / normu;
971  // u[1] = u[1] / normu;
972  // u[2] = u[2] / normu;
973 
974  // Double_t normpPlus = sqrt( (pPlus[0]*pPlus[0]) + (pPlus[1]*pPlus[1]) + (pPlus[2]*pPlus[2]) );
975  // Double_t normpMinus = sqrt( (pMinus[0]*pMinus[0]) + (pMinus[1]*pMinus[1]) + (pMinus[2]*pMinus[2]) );
976 
977  // pPlus[0] = pPlus[0] / normpPlus;
978  // pPlus[1] = pPlus[1] / normpPlus;
979  // pPlus[2] = pPlus[2] / normpPlus;
980 
981  // pMinus[0] = pMinus[0] / normpMinus;
982  // pMinus[1] = pMinus[1] / normpMinus;
983  // pMinus[2] = pMinus[2] / normpMinus;
984 
985  // Double_t v[3] = {0,0,0}; // pPlus X pMinus
986  // v[0] = (pPlus[1]*pMinus[2]) - (pPlus[2]*pMinus[1]);
987  // v[1] = (pPlus[2]*pMinus[0]) - (pPlus[0]*pMinus[2]);
988  // v[2] = (pPlus[0]*pMinus[1]) - (pPlus[1]*pMinus[0]);
989 
990  // Double_t w[3] = {0,0,0}; // u X v
991  // w[0] = (u[1]*v[2]) - (u[2]*v[1]);
992  // w[1] = (u[2]*v[0]) - (u[0]*v[2]);
993  // w[2] = (u[0]*v[1]) - (u[1]*v[0]);
994 
995  // Double_t z[3] = {0,0,1};
996  // Double_t wc[3] = {0,0,0}; // u X v
997  // wc[0] = (u[1]*z[2]) - (u[2]*z[1]);
998  // wc[1] = (u[2]*z[0]) - (u[0]*z[2]);
999  // wc[2] = (u[0]*z[1]) - (u[1]*z[0]);
1000 
1001  // Double_t PhiV = TMath::ACos((w[0]*wc[0]) + (w[1]*wc[1]) + (w[2]*wc[2]));
1002  //return TMath::Abs(PhiV);
1003 
1004 
1005  // TVector3 pPlus(pt.Px(),pt.Py(),pt.Pz());
1006  // TVector3 pMinus(nt.Px(),nt.Py(),nt.Pz());
1007 
1008  // TVector3 u = pMinus + pPlus;
1009  // u = u*(1/u.Mag());
1010 
1011  // TVector3 pHPlus = pPlus*(1/pPlus.Mag());
1012  // TVector3 pHMinus = pMinus*(1/pMinus.Mag());
1013 
1014  // TVector3 v = pHPlus.Cross(pHMinus);
1015  // TVector3 w = u.Cross(v);
1016  // TVector3 z(0,0,1);
1017  // TVector3 wc = u.Cross(z);
1018 
1019  // Double_t PhiV = w * wc;
1020 
1021  Double_t mn[3] = {0,0,0};
1022  Double_t mp[3] = {0,0,0};
1023 
1024  v0->GetNPxPyPz(mn[0],mn[1],mn[2]);//reconstructed cartesian momentum components of negative daughter;
1025  v0->GetPPxPyPz(mp[0],mp[1],mp[2]);//reconstructed cartesian momentum components of positive daughter;
1026 
1027  Double_t deltat = 1.;
1028  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
1029  Double_t radiussum = TMath::Sqrt(xyz[0]*xyz[0] + xyz[1]*xyz[1]) + 50;//radius to which tracks shall be propagated
1030 
1031  Double_t momPosProp[3] = {0,0,0};
1032  Double_t momNegProp[3] = {0,0,0};
1033 
1034  Double_t psiPair = 4.;
1035  // cout<< "Momentum before propagation at radius ::"<< TMath::Sqrt(xyz[0]*xyz[0] + xyz[1]*xyz[1])<< endl;
1036  // cout<< mn[0]<< " " <<mn[1]<<" "<< mn[2]<<" "<<TMath::Sqrt(mn[0]*mn[0] + mn[1]*mn[1] +mn[2]*mn[2])<< endl;
1037  // cout<< mp[0]<< " " <<mp[1]<<" "<< mp[2]<<" "<<TMath::Sqrt(mp[0]*mp[0] + mp[1]*mp[1] +mp[2]*mp[2])<< endl;
1038  Double_t pEle,pPos;
1039 
1040  if (fImprovedPsiPair==1 || fImprovedPsiPair==0 ){
1041  if(nt.PropagateTo(radiussum,magField) == 0) return psiPair; //propagate tracks to the outside -> Better Purity and Efficiency
1042  if(pt.PropagateTo(radiussum,magField) == 0) return psiPair; //propagate tracks to the outside -> Better Purity and Efficiency
1043 
1044  pt.GetPxPyPz(momPosProp);//Get momentum vectors of tracks after propagation
1045  nt.GetPxPyPz(momNegProp);
1046  } else if (fImprovedPsiPair>=2) {
1047  momPosProp[0] = pt.GetParameterAtRadius(radiussum,magField,3);
1048  momPosProp[1] = pt.GetParameterAtRadius(radiussum,magField,4);
1049  momPosProp[2] = pt.GetParameterAtRadius(radiussum,magField,5);
1050 
1051  momNegProp[0] = nt.GetParameterAtRadius(radiussum,magField,3);
1052  momNegProp[1] = nt.GetParameterAtRadius(radiussum,magField,4);
1053  momNegProp[2] = nt.GetParameterAtRadius(radiussum,magField,5);
1054  pEle = TMath::Sqrt(momNegProp[0]*momNegProp[0]+momNegProp[1]*momNegProp[1]+momNegProp[2]*momNegProp[2]);//absolute momentum value of negative daughter
1055 
1056  pPos = TMath::Sqrt(momPosProp[0]*momPosProp[0]+momPosProp[1]*momPosProp[1]+momPosProp[2]*momPosProp[2]);//absolute momentum value of positive daughter
1057 
1058  if ( (pEle==0 || pPos==0) && fImprovedPsiPair==3) {
1059  radiussum = TMath::Sqrt(xyz[0]*xyz[0] + xyz[1]*xyz[1]) + 30;
1060  momPosProp[0] = pt.GetParameterAtRadius(radiussum,magField,3);
1061  momPosProp[1] = pt.GetParameterAtRadius(radiussum,magField,4);
1062  momPosProp[2] = pt.GetParameterAtRadius(radiussum,magField,5);
1063 
1064  momNegProp[0] = nt.GetParameterAtRadius(radiussum,magField,3);
1065  momNegProp[1] = nt.GetParameterAtRadius(radiussum,magField,4);
1066  momNegProp[2] = nt.GetParameterAtRadius(radiussum,magField,5);
1067 
1068  }
1069  }
1070 
1071  pEle = TMath::Sqrt(momNegProp[0]*momNegProp[0]+momNegProp[1]*momNegProp[1]+momNegProp[2]*momNegProp[2]);//absolute momentum value of negative daughter
1072  pPos = TMath::Sqrt(momPosProp[0]*momPosProp[0]+momPosProp[1]*momPosProp[1]+momPosProp[2]*momPosProp[2]);//absolute momentum value of positive daughter
1073 
1074  Double_t scalarproduct =
1075  momPosProp[0]*momNegProp[0]+momPosProp[1]*momNegProp[1]+momPosProp[2]*momNegProp[2];//scalar product of propagated positive and negative daughters' momenta
1076 
1077  if (pEle==0 || pPos==0) return psiPair;
1078 
1079 
1080  Double_t chipair = TMath::ACos(scalarproduct/(pEle*pPos));//Angle between propagated daughter tracks
1081 
1082  // psiPair = TMath::Abs(TMath::ASin(deltat/chipair));
1083  psiPair = TMath::ASin(deltat/chipair);
1084  return psiPair;
1085 }
1086 
1088 Bool_t AliV0ReaderV1::GetHelixCenter(const AliExternalTrackParam *track,Double_t center[2]){
1089 
1090  // Get Center of the helix track parametrization
1091 
1092  Int_t charge=track->Charge();
1093  Double_t b=fInputEvent->GetMagneticField();
1094 
1095  Double_t helix[6];
1096  track->GetHelixParameters(helix,b);
1097 
1098  Double_t xpos = helix[5];
1099  Double_t ypos = helix[0];
1100  Double_t radius = TMath::Abs(1./helix[4]);
1101  Double_t phi = helix[2];
1102 
1103  if(phi < 0){
1104  phi = phi + 2*TMath::Pi();
1105  }
1106 
1107  phi -= TMath::Pi()/2.;
1108  Double_t xpoint = radius * TMath::Cos(phi);
1109  Double_t ypoint = radius * TMath::Sin(phi);
1110 
1111  if(b<0){
1112  if(charge > 0){
1113  xpoint = - xpoint;
1114  ypoint = - ypoint;
1115  }
1116  }
1117  if(b>0){
1118  if(charge < 0){
1119  xpoint = - xpoint;
1120  ypoint = - ypoint;
1121  }
1122  }
1123  center[0] = xpos + xpoint;
1124  center[1] = ypos + ypoint;
1125 
1126  return 1;
1127 }
1129 Bool_t AliV0ReaderV1::GetConversionPoint(const AliExternalTrackParam *pparam,const AliExternalTrackParam *nparam,Double_t convpos[3],Double_t dca[2]){
1130 
1131  // Recalculate Conversion Point
1132 
1133  if(!pparam||!nparam)return kFALSE;
1134 
1135  Double_t helixcenterpos[2];
1136  GetHelixCenter(pparam,helixcenterpos);
1137 
1138  Double_t helixcenterneg[2];
1139  GetHelixCenter(nparam,helixcenterneg);
1140 
1141  Double_t helixpos[6];
1142  pparam->GetHelixParameters(helixpos,fInputEvent->GetMagneticField());
1143  Double_t posradius = TMath::Abs(1./helixpos[4]);
1144 
1145  Double_t helixneg[6];
1146  nparam->GetHelixParameters(helixneg,fInputEvent->GetMagneticField());
1147  Double_t negradius = TMath::Abs(1./helixneg[4]);
1148 
1149  // Calculate xy-position
1150 
1151  Double_t xpos = helixcenterpos[0];
1152  Double_t ypos = helixcenterpos[1];
1153  Double_t xneg = helixcenterneg[0];
1154  Double_t yneg = helixcenterneg[1];
1155 
1156  convpos[0] = (xpos*negradius + xneg*posradius)/(negradius+posradius);
1157  convpos[1] = (ypos*negradius+ yneg*posradius)/(negradius+posradius);
1158 
1159 
1160  // Calculate Track XY vertex position
1161 
1162  Double_t deltaXPos = convpos[0] - xpos;
1163  Double_t deltaYPos = convpos[1] - ypos;
1164 
1165  Double_t deltaXNeg = convpos[0] - xneg;
1166  Double_t deltaYNeg = convpos[1] - yneg;
1167 
1168  Double_t alphaPos = TMath::Pi() + TMath::ATan2(-deltaYPos,-deltaXPos);
1169  Double_t alphaNeg = TMath::Pi() + TMath::ATan2(-deltaYNeg,-deltaXNeg);
1170 
1171  Double_t vertexXNeg = xneg + TMath::Abs(negradius)*TMath::Cos(alphaNeg);
1172  Double_t vertexYNeg = yneg + TMath::Abs(negradius)*TMath::Sin(alphaNeg);
1173 
1174  Double_t vertexXPos = xpos + TMath::Abs(posradius)*TMath::Cos(alphaPos);
1175  Double_t vertexYPos = ypos + TMath::Abs(posradius)*TMath::Sin(alphaPos);
1176 
1177  AliExternalTrackParam p(*pparam);
1178  AliExternalTrackParam n(*nparam);
1179 
1180  TVector2 vertexPos(vertexXPos,vertexYPos);
1181  TVector2 vertexNeg(vertexXNeg,vertexYNeg);
1182 
1183  // Convert to local coordinate system
1184  TVector2 vertexPosRot=vertexPos.Rotate(-p.GetAlpha());
1185  TVector2 vertexNegRot=vertexNeg.Rotate(-n.GetAlpha());
1186 
1187  // Propagate Track Params to Vertex
1188 
1189  if(!p.PropagateTo(vertexPosRot.X(),fInputEvent->GetMagneticField()))return kFALSE;
1190  if(!n.PropagateTo(vertexNegRot.X(),fInputEvent->GetMagneticField()))return kFALSE;
1191 
1192  // Check whether propagation was sucessful
1193 
1194  if(TMath::Abs(vertexPos.Mod()-TMath::Sqrt(p.GetX()*p.GetX()+p.GetY()*p.GetY()))>0.01)return kFALSE;
1195  if(TMath::Abs(vertexNeg.Mod()-TMath::Sqrt(n.GetX()*n.GetX()+n.GetY()*n.GetY()))>0.01)return kFALSE;
1196 
1197  // Calculate z position
1198 
1199  convpos[2] = (p.GetZ()*negradius+n.GetZ()*posradius)/(negradius+posradius);
1200 
1201  // Calculate DCA
1202  TVector2 vdca=vertexPos-vertexNeg;
1203  dca[0]=vdca.Mod();
1204  dca[1]=TMath::Abs(n.GetZ()-p.GetZ());
1205 
1206  return kTRUE;
1207 }
1208 
1209 //________________________________________________________________________
1211 
1212  // Get reconstructed conversion photons from satellite AOD file
1213 
1214  AliAODEvent *fAODEvent=dynamic_cast<AliAODEvent*>(fInputEvent);
1215 
1216  if(fAODEvent){
1217 
1218  if(fConversionGammas == NULL){
1219  fConversionGammas = new TClonesArray("AliAODConversionPhoton",100);
1220  }
1221  fConversionGammas->Delete();//Reset the TClonesArray
1222 
1223  //Get Gammas from satellite AOD gamma branch
1224 
1225  AliAODConversionPhoton *gamma=0x0;
1226 
1227  TClonesArray *fInputGammas=dynamic_cast<TClonesArray*>(fAODEvent->FindListObject(fDeltaAODBranchName.Data()));
1228 
1229  if(!fInputGammas){
1230  FindDeltaAODBranchName();
1231  fInputGammas=dynamic_cast<TClonesArray*>(fAODEvent->FindListObject(fDeltaAODBranchName.Data()));}
1232  if(!fInputGammas){AliError("No Gamma Satellites found");return kFALSE;}
1233  // Apply Selection Cuts to Gammas and create local working copy
1234  if(fInputGammas){
1235  for(Int_t i=0;i<fInputGammas->GetEntriesFast();i++){
1236  gamma=dynamic_cast<AliAODConversionPhoton*>(fInputGammas->At(i));
1237  if(gamma){
1238  if(fRelabelAODs)RelabelAODPhotonCandidates(gamma);
1239  if(fConversionCuts->PhotonIsSelected(gamma,fInputEvent)){
1240  new((*fConversionGammas)[fConversionGammas->GetEntriesFast()]) AliAODConversionPhoton(*gamma);}
1241  }
1242  }
1243  }
1244  }
1245 
1246  if(fConversionGammas->GetEntries()){return kTRUE;}
1247 
1248  return kFALSE;
1249 }
1250 
1251 //________________________________________________________________________
1253 
1254  // Find delta AOD branchname containing reconstructed photons
1255 
1256  TList *list=fInputEvent->GetList();
1257  for(Int_t ii=0;ii<list->GetEntries();ii++){
1258  TString name((list->At(ii))->GetName());
1259  if(name.BeginsWith(fDeltaAODBranchName)&&name.EndsWith("gamma")){
1260  fDeltaAODBranchName=name;
1261  AliInfo(Form("Set DeltaAOD BranchName to: %s",fDeltaAODBranchName.Data()));
1262  return;
1263  }
1264  }
1265 }
1266 
1267 //________________________________________________________________________
1269 
1270  if(fPreviousV0ReaderPerformsAODRelabeling == 2) return;
1271  else if(fPreviousV0ReaderPerformsAODRelabeling == 0){
1272  printf("Running AODs! Determine if V0Reader '%s' should perform relabeling\n",this->GetName());
1273  TObjArray* obj = (TObjArray*)AliAnalysisManager::GetAnalysisManager()->GetTasks();
1274  Int_t iPosition = obj->IndexOf(this);
1275  Bool_t prevV0ReaderRunningButNotRelabeling = kFALSE;
1276  for(Int_t i=iPosition-1; i>=0; i--){
1277  if( (obj->At(i))->IsA() == AliV0ReaderV1::Class()){
1278  AliV0ReaderV1* tempReader = (AliV0ReaderV1*) obj->At(i);
1279  if( tempReader->AreAODsRelabeled() && tempReader->IsReaderPerformingRelabeling() == 1){
1280  fPreviousV0ReaderPerformsAODRelabeling = 2;
1281  prevV0ReaderRunningButNotRelabeling = kFALSE;
1282  printf("V0Reader '%s' is running before this V0Reader '%s': do _NOT_ relabel AODs by current reader!\n",tempReader->GetName(),this->GetName());
1283  break;
1284  }else prevV0ReaderRunningButNotRelabeling = kTRUE;
1285  }
1286  }
1287  if(prevV0ReaderRunningButNotRelabeling) AliFatal(Form("There are V0Readers before '%s', but none of them is relabeling!",this->GetName()));
1288 
1289  if(fPreviousV0ReaderPerformsAODRelabeling == 2) return;
1290  else{
1291  printf("This V0Reader '%s' is first to be processed: do relabel AODs by current reader!\n",this->GetName());
1292  fPreviousV0ReaderPerformsAODRelabeling = 1;
1293  }
1294  }
1295 
1296  if(fPreviousV0ReaderPerformsAODRelabeling != 1) AliFatal(Form("In %s: fPreviousV0ReaderPerformsAODRelabeling = '%i' - while it should be impossible it is something different than '1'!",this->GetName(),fPreviousV0ReaderPerformsAODRelabeling));
1297 
1298  // Relabeling For AOD Event
1299  // ESDiD -> AODiD
1300  // MCLabel -> AODMCLabel
1301  Bool_t AODLabelPos = kFALSE;
1302  Bool_t AODLabelNeg = kFALSE;
1303 
1304  for(Int_t i = 0; i<fInputEvent->GetNumberOfTracks();i++){
1305  AliAODTrack *tempDaughter = static_cast<AliAODTrack*>(fInputEvent->GetTrack(i));
1306  if(!AODLabelPos){
1307  if( tempDaughter->GetID() == PhotonCandidate->GetTrackLabelPositive() ){
1308  PhotonCandidate->SetMCLabelPositive(TMath::Abs(tempDaughter->GetLabel()));
1309  PhotonCandidate->SetLabelPositive(i);
1310  AODLabelPos = kTRUE;
1311  }
1312  }
1313  if(!AODLabelNeg){
1314  if( tempDaughter->GetID() == PhotonCandidate->GetTrackLabelNegative()){
1315  PhotonCandidate->SetMCLabelNegative(TMath::Abs(tempDaughter->GetLabel()));
1316  PhotonCandidate->SetLabelNegative(i);
1317  AODLabelNeg = kTRUE;
1318  }
1319  }
1320  if(AODLabelNeg && AODLabelPos){
1321  return;
1322  }
1323  }
1324  if(!AODLabelPos || !AODLabelNeg){
1325  AliError(Form("NO AOD Daughters Found Pos: %i %i Neg: %i %i, setting all labels to -999999",AODLabelPos,PhotonCandidate->GetTrackLabelPositive(),AODLabelNeg,PhotonCandidate->GetTrackLabelNegative()));
1326  if(!AODLabelNeg){
1327  PhotonCandidate->SetMCLabelNegative(-999999);
1328  PhotonCandidate->SetLabelNegative(-999999);
1329  }
1330  if(!AODLabelPos){
1331  PhotonCandidate->SetMCLabelPositive(-999999);
1332  PhotonCandidate->SetLabelPositive(-999999);
1333  }
1334  }
1335 
1336 }
1337 
1338 //************************************************************************
1339 // This function counts the number of primary tracks in the event, the
1340 // implementation for ESD and AOD had to be different due to the different
1341 // filters which are already applied on AOD level the so-called filter
1342 // bits, we tried to replicate the conditions in both but an exact match
1343 // could not be reached. The cuts are similar to the ones used by the
1344 // jet-group
1345 //************************************************************************
1346 //________________________________________________________________________
1348 
1349  // In Case of MC count only MB tracks
1350  // if(fMCEvent && fInputEvent->IsA()==AliESDEvent::Class()){
1351  // fEventCuts->GetNotRejectedParticles(1,NULL,fMCEvent);
1352  // }
1353  // if(fMCEvent && fInputEvent->IsA()==AliAODEvent::Class()){
1354  // fEventCuts->GetNotRejectedParticles(1,NULL,fInputEvent);
1355  // }
1356 
1357  if(fInputEvent->IsA()==AliESDEvent::Class()){
1358  static AliESDtrackCuts *EsdTrackCuts = 0x0;
1359  static int prevRun = -1;
1360  // Using standard function for setting Cuts
1361  Int_t runNumber = fInputEvent->GetRunNumber();
1362 
1363  if (prevRun!=runNumber) {
1364  delete EsdTrackCuts;
1365  EsdTrackCuts = 0;
1366  prevRun = runNumber;
1367  }
1368  if (!EsdTrackCuts) {
1369  // if LHC11a or earlier or if LHC13g or if LHC12a-i -> use 2010 cuts
1370  if( (runNumber<=146860) || (runNumber>=197470 && runNumber<=197692) || (runNumber>=172440 && runNumber<=193766) ){
1371  EsdTrackCuts = AliESDtrackCuts::GetStandardITSTPCTrackCuts2010();
1372  // else if run2 data use 2015 PbPb cuts
1373  } else if (runNumber>=209122){
1374  //EsdTrackCuts = AliESDtrackCuts::GetStandardITSTPCTrackCuts2015PbPb();
1375  // hard coded track cuts for the moment, because AliESDtrackCuts::GetStandardITSTPCTrackCuts2015PbPb() gives spams warnings
1376  EsdTrackCuts = new AliESDtrackCuts();
1377  // TPC; clusterCut = 1, cutAcceptanceEdges = kTRUE, removeDistortedRegions = kFALSE
1378  EsdTrackCuts->AliESDtrackCuts::SetMinNCrossedRowsTPC(70);
1379  EsdTrackCuts->AliESDtrackCuts::SetMinRatioCrossedRowsOverFindableClustersTPC(0.8);
1380  EsdTrackCuts->SetCutGeoNcrNcl(2., 130., 1.5, 0.0, 0.0); // only dead zone and not clusters per length
1381  //EsdTrackCuts->AliESDtrackCuts::SetCutOutDistortedRegionsTPC(kTRUE);
1382  EsdTrackCuts->AliESDtrackCuts::SetMaxChi2PerClusterTPC(4);
1383  EsdTrackCuts->AliESDtrackCuts::SetAcceptKinkDaughters(kFALSE);
1384  EsdTrackCuts->AliESDtrackCuts::SetRequireTPCRefit(kTRUE);
1385  // ITS; selPrimaries = 1
1386  EsdTrackCuts->AliESDtrackCuts::SetRequireITSRefit(kTRUE);
1387  EsdTrackCuts->AliESDtrackCuts::SetClusterRequirementITS(AliESDtrackCuts::kSPD,
1388  AliESDtrackCuts::kAny);
1389  EsdTrackCuts->AliESDtrackCuts::SetMaxDCAToVertexXYPtDep("0.0105+0.0350/pt^1.1");
1390  EsdTrackCuts->AliESDtrackCuts::SetMaxChi2TPCConstrainedGlobal(36);
1391  EsdTrackCuts->AliESDtrackCuts::SetMaxDCAToVertexZ(2);
1392  EsdTrackCuts->AliESDtrackCuts::SetDCAToVertex2D(kFALSE);
1393  EsdTrackCuts->AliESDtrackCuts::SetRequireSigmaToVertex(kFALSE);
1394  EsdTrackCuts->AliESDtrackCuts::SetMaxChi2PerClusterITS(36);
1395  // else use 2011 version of track cuts
1396  }else{
1397  EsdTrackCuts = AliESDtrackCuts::GetStandardITSTPCTrackCuts2011();
1398  }
1399  EsdTrackCuts->SetMaxDCAToVertexZ(2);
1400  EsdTrackCuts->SetEtaRange(-0.8, 0.8);
1401  EsdTrackCuts->SetPtRange(0.15);
1402  }
1403  fNumberOfPrimaryTracks = 0;
1404  for(Int_t iTracks = 0; iTracks < fInputEvent->GetNumberOfTracks(); iTracks++){
1405  AliESDtrack* curTrack = (AliESDtrack*) fInputEvent->GetTrack(iTracks);
1406  if(!curTrack) continue;
1407  if(!EsdTrackCuts->AcceptTrack(curTrack)) continue;
1408  //if(fMCEvent && !(fEventCuts->IsParticleFromBGEvent(TMath::Abs(curTrack->GetLabel()),fMCEvent->Stack(),fInputEvent))) continue;
1409  fNumberOfPrimaryTracks++;
1410  }
1411  }
1412  else if(fInputEvent->IsA()==AliAODEvent::Class()){
1413  fNumberOfPrimaryTracks = 0;
1414  for(Int_t iTracks = 0; iTracks<fInputEvent->GetNumberOfTracks(); iTracks++){
1415  AliAODTrack* curTrack = (AliAODTrack*) fInputEvent->GetTrack(iTracks);
1416  if(curTrack->GetID()<0) continue; // Avoid double counting of tracks
1417  if(!curTrack->IsHybridGlobalConstrainedGlobal()) continue;
1418  if(TMath::Abs(curTrack->Eta())>0.8) continue;
1419  if(curTrack->Pt()<0.15) continue;
1420  //if(fMCEvent && !(fEventCuts->IsParticleFromBGEvent(TMath::Abs(curTrack->GetLabel()),NULL,fInputEvent))) continue;
1421  //if(TMath::Abs(curTrack->ZAtDCA())>2) continue; // Only Set For TPCOnly tracks
1422  fNumberOfPrimaryTracks++;
1423  }
1424  }
1425 
1426  return;
1427 }
1428 
1431  fNumberOfTPCoutTracks = 0;
1432 
1433  for (Int_t itrk = 0; itrk < fInputEvent->GetNumberOfTracks(); itrk++) {
1434  AliVTrack *trk = dynamic_cast<AliVTrack*>(fInputEvent->GetTrack(itrk));
1435 
1436  if (trk != NULL) {
1437  /* the initial method of counting TPC out tracks */
1438  if (!(trk->Pt() < 0.15) && (TMath::Abs(trk->Eta()) < 0.8)) {
1439  if ((trk->GetStatus() & AliVTrack::kTPCout) == AliVTrack::kTPCout) {
1440  fNumberOfTPCoutTracks++;
1441  }
1442  }
1443  }
1444  }
1445  return;
1446 }
1447 
1448 
1450 Bool_t AliV0ReaderV1::ParticleIsConvertedPhoton(AliMCEvent *mcEvent, TParticle *particle, Double_t etaMax, Double_t rMax, Double_t zMax){
1451  // MonteCarlo Photon Selection
1452  if(!mcEvent)return kFALSE;
1453 
1454  if (particle->GetPdgCode() == 22){
1455  // check whether particle is within eta range
1456  if( TMath::Abs(particle->Eta()) > etaMax ) return kFALSE;
1457  // check if particle doesn't have a photon as mother
1458  if(particle->GetMother(0) >-1 && mcEvent->Particle(particle->GetMother(0))->GetPdgCode() == 22){
1459  return kFALSE; // no photon as mothers!
1460  }
1461  // looking for conversion gammas (electron + positron from pairbuilding (= 5) )
1462  TParticle* ePos = NULL;
1463  TParticle* eNeg = NULL;
1464  if(particle->GetNDaughters() >= 2){
1465  for(Int_t daughterIndex=particle->GetFirstDaughter();daughterIndex<=particle->GetLastDaughter();daughterIndex++){
1466  if(daughterIndex<0) continue;
1467  TParticle *tmpDaughter = mcEvent->Particle(daughterIndex);
1468  if(tmpDaughter->GetUniqueID() == 5){
1469  if(tmpDaughter->GetPdgCode() == 11){
1470  eNeg = tmpDaughter;
1471  } else if(tmpDaughter->GetPdgCode() == -11){
1472  ePos = tmpDaughter;
1473  }
1474  }
1475  }
1476  }
1477  if(ePos == NULL || eNeg == NULL){ // means we do not have two daughters from pair production
1478  return kFALSE;
1479  }
1480  // check if electrons are in correct eta window
1481  if( TMath::Abs(ePos->Eta()) > etaMax ||
1482  TMath::Abs(eNeg->Eta()) > etaMax )
1483  return kFALSE;
1484 
1485  // check if photons have converted in reconstructable range
1486  if(ePos->R() > rMax){
1487  return kFALSE; // cuts on distance from collision point
1488  }
1489  if(TMath::Abs(ePos->Vz()) > zMax){
1490  return kFALSE; // outside material
1491  }
1492  if(TMath::Abs(eNeg->Vz()) > zMax){
1493  return kFALSE; // outside material
1494  }
1495 
1496 
1497  Double_t lineCutZRSlope = tan(2*atan(exp(-etaMax)));
1498  Double_t lineCutZValue = 7.;
1499  if( ePos->R() <= ((TMath::Abs(ePos->Vz()) * lineCutZRSlope) - lineCutZValue)){
1500  return kFALSE; // line cut to exclude regions where we do not reconstruct
1501  }
1502  if( eNeg->R() <= ((TMath::Abs(eNeg->Vz()) * lineCutZRSlope) - lineCutZValue)){
1503  return kFALSE; // line cut to exclude regions where we do not reconstruct
1504  }
1505 
1506  return kTRUE;
1507  }
1508  return kFALSE;
1509 }
1510 
1513  const AliVVertex* primVtxMC = fMCEvent->GetPrimaryVertex();
1514  Double_t mcProdVtxX = primVtxMC->GetX();
1515  Double_t mcProdVtxY = primVtxMC->GetY();
1516  Double_t mcProdVtxZ = primVtxMC->GetZ();
1517  //cout << mcProdVtxX <<"\t" << mcProdVtxY << "\t" << mcProdVtxZ << endl;
1518 
1519  // Loop over all primary MC particle
1520  for(Long_t i = 0; i < fMCEvent->GetNumberOfTracks(); i++) {
1521  if (fEventCuts->IsConversionPrimaryESD( fMCEvent, i, mcProdVtxX, mcProdVtxY, mcProdVtxZ)){
1522  // fill primary histogram
1523  TParticle* particle = (TParticle *)fMCEvent->Particle(i);
1524  if (!particle) continue;
1525  if (ParticleIsConvertedPhoton(fMCEvent, particle, 0.9, 180.,250. )){
1526  if(particle->GetFirstDaughter()<0) continue;
1527  TParticle *tmpDaughter = fMCEvent->Particle(particle->GetFirstDaughter());
1528  if (!tmpDaughter) continue;
1529  fHistoMCGammaPtvsR->Fill(particle->Pt(),tmpDaughter->R());
1530  fHistoMCGammaPtvsPhi->Fill(particle->Pt(),particle->Phi());
1531  fHistoMCGammaRvsPhi->Fill(tmpDaughter->R(),particle->Phi());
1532  }
1533  if (ParticleIsConvertedPhoton(fMCEvent, particle, 1.4, 180.,250. )){
1534  if(particle->GetFirstDaughter()<0) continue;
1535  TParticle *tmpDaughter = fMCEvent->Particle(particle->GetFirstDaughter());
1536  if (!tmpDaughter) continue;
1537  fHistoMCGammaPtvsEta->Fill(particle->Pt(),particle->Eta());
1538  fHistoMCGammaRvsEta->Fill(tmpDaughter->R(),particle->Eta());
1539  fHistoMCGammaPhivsEta->Fill(particle->Phi(),particle->Eta());
1540  }
1541  }
1542  }
1543 }
1544 
1547 
1548  const AliVVertex* primVtxMC = fMCEvent->GetPrimaryVertex();
1549  Double_t mcProdVtxX = primVtxMC->GetX();
1550  Double_t mcProdVtxY = primVtxMC->GetY();
1551  Double_t mcProdVtxZ = primVtxMC->GetZ();
1552 // cout << mcProdVtxX <<"\t" << mcProdVtxY << "\t" << mcProdVtxZ << endl;
1553 
1554  Int_t tracklabelPos=currentV0->GetPindex();
1555  Int_t tracklabelNeg=currentV0->GetNindex();
1556 
1557  Int_t labelp=TMath::Abs(fConversionCuts->GetTrack(fInputEvent,tracklabelPos)->GetLabel());
1558  Int_t labeln=TMath::Abs(fConversionCuts->GetTrack(fInputEvent,tracklabelNeg)->GetLabel());
1559 
1560  TParticle* negPart = 0x0;
1561  if(labeln>-1) negPart = (TParticle *)fMCEvent->Particle(labeln);
1562  TParticle* posPart = 0x0;
1563  if(labelp>-1) posPart = (TParticle *)fMCEvent->Particle(labelp);
1564 
1565  if ( negPart == NULL || posPart == NULL ) return;
1566 // if (!(negPart->GetPdgCode() == 11)) return;
1567 // if (!(posPart->GetPdgCode() == -11)) return;
1568  Long_t motherlabelNeg = negPart->GetFirstMother();
1569  Long_t motherlabelPos = posPart->GetFirstMother();
1570 
1571 // cout << "mother neg " << motherlabelNeg << " mother pos " << motherlabelPos << endl;
1572  if (motherlabelNeg>-1 && motherlabelNeg == motherlabelPos && negPart->GetFirstMother() != -1){
1573  if (fEventCuts->IsConversionPrimaryESD( fMCEvent, negPart->GetFirstMother(), mcProdVtxX, mcProdVtxY, mcProdVtxZ)){
1574 
1575  TParticle* mother = (TParticle *)fMCEvent->Particle(motherlabelNeg);
1576  if (mother->GetPdgCode() == 22 ){
1577  if (!CheckVectorForDoubleCount(fVectorFoundGammas,motherlabelNeg ) ){
1578  if (ParticleIsConvertedPhoton(fMCEvent, mother, 0.9, 180.,250. )){
1579  fHistoRecMCGammaPtvsR->Fill(mother->Pt(),negPart->R());
1580  fHistoRecMCGammaPtvsPhi->Fill(mother->Pt(),mother->Phi());
1581  fHistoRecMCGammaRvsPhi->Fill(negPart->R(),mother->Phi());
1582  }
1583  if (ParticleIsConvertedPhoton(fMCEvent, mother, 1.4, 180.,250. )){
1584  fHistoRecMCGammaPtvsEta->Fill(mother->Pt(),mother->Eta());
1585  fHistoRecMCGammaRvsEta->Fill(negPart->R(),mother->Eta());
1586  fHistoRecMCGammaPhivsEta->Fill(mother->Phi(),mother->Eta());
1587  }
1588 // cout << "new gamma found" << endl;
1589  } else {
1590  if (ParticleIsConvertedPhoton(fMCEvent, mother, 0.9, 180.,250. )){
1591  fHistoRecMCGammaMultiPt->Fill(mother->Pt());
1592  fHistoRecMCGammaMultiPhi->Fill(mother->Phi());
1593  fHistoRecMCGammaMultiR->Fill(negPart->R());
1594  }
1595  if (ParticleIsConvertedPhoton(fMCEvent, mother, 1.4, 180.,250. )){
1596  fHistoRecMCGammaMultiPtvsEta->Fill(mother->Pt(),mother->Eta());
1597  }
1598 // cout << "this one I had already: " << motherlabelNeg << endl << "-----------------------" << endl;
1599  }
1600 // cout << "event gammas: " << endl;
1601 // for(Int_t iGamma=0; iGamma<fVectorFoundGammas.size(); iGamma++){cout << fVectorFoundGammas.at(iGamma) << ", ";}
1602  }
1603  }
1604  }
1605 }
1606 
1607 //_________________________________________________________________________________
1608 void AliV0ReaderV1::FillImpactParamHistograms( AliVTrack* pTrack, AliVTrack* nTrack, AliESDv0 *fCurrentV0, AliKFConversionPhoton *fCurrentMotherKF){
1609 
1610  // values of cuts to be introduced in 2016 to reduce ESD size
1611  Float_t fZmax = 25; //cm
1612  Float_t fYmax = 25; //cm
1613  Double_t kTPCMargin = 1.0; //cm
1614  Int_t kMaxTPCV0Conflicts = 1; //# conflicting clusters tolerated
1615 
1616  //conversion point
1617  Double_t convX, convY, convZ;
1618  fCurrentV0->GetXYZ(convX,convY,convZ);
1619  Double_t convR = TMath::Sqrt(convX*convX+convY*convY);
1620  //recalculated conversion point
1621  Double_t convXrecalc=fCurrentMotherKF->GetConversionX();
1622  Double_t convYrecalc=fCurrentMotherKF->GetConversionY();
1623  Double_t convZrecalc=fCurrentMotherKF->GetConversionZ();
1624  Double_t convRrecalc = fCurrentMotherKF->GetConversionRadius();
1625 
1626  //count V0s
1627  fHistoImpactParameterStudy->Fill(0);
1628 
1629  //count V0s with two TPC-only tracks
1630  AliESDtrack* positiveTrack = (AliESDtrack*) pTrack;
1631  AliESDtrack* negativeTrack = (AliESDtrack*) nTrack;
1632  Bool_t TPConly = kFALSE;
1633  if (!positiveTrack->IsOn(AliESDtrack::kITSin) && !negativeTrack->IsOn(AliESDtrack::kITSin)){
1634  fHistoImpactParameterStudy->Fill(1);
1635  TPConly = kTRUE;
1636  }
1637 
1638  Bool_t RemovedByZcut = kFALSE;
1639  Bool_t RemovedByYcut = kFALSE;
1640 
1641  //count V0s which have...
1642  if(TPConly){
1643  //not passed z cut: pos.tr. or neg.tr.
1644  if(((TMath::Abs(positiveTrack->GetZ()))>fZmax) || ((TMath::Abs(negativeTrack->GetZ()))>fZmax)){
1645  fHistoImpactParameterStudy->Fill(2);
1646  RemovedByZcut=kTRUE;
1647  }
1648 
1649  //not passed y cut: pos.tr. or neg.tr.
1650  if(((TMath::Abs(positiveTrack->GetY()))>fYmax) || ((TMath::Abs(negativeTrack->GetY()))>fYmax)){
1651  fHistoImpactParameterStudy->Fill(3);
1652  RemovedByYcut=kTRUE;
1653  }
1654  }
1655 
1656  //causality cut
1657  Bool_t RemovedByCausality=kFALSE;
1658  AliESDEvent *fESDEvent=dynamic_cast<AliESDEvent*>(fInputEvent);
1659  const Float_t rTPC[159]={
1660  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,
1661  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,
1662  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,
1663  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,
1664  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,
1665  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,
1666  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,
1667  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,
1668  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,
1669  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,
1670  233.850,235.350,236.850,238.350,239.850,241.350,242.850,244.350,245.850};
1671 
1672  // fill conversion radius histograms
1673  fHistoR->Fill(convR);
1674  fHistoRrecalc->Fill(convRrecalc);
1675  Double_t alpha = TMath::ATan2(convY,convX);
1676  if (alpha<0) alpha += TMath::Pi()*2;
1677  Int_t sec = alpha/(TMath::Pi()/9);
1678  alpha = (10.+sec*20.)*TMath::DegToRad();
1679  Double_t cs = TMath::Cos(alpha);
1680  Double_t sn = TMath::Sin(alpha);
1681  Double_t xsV0 = convX*cs - convY*sn;
1682  fHistoRviaAlpha->Fill(xsV0);
1683  Double_t alpha_r = TMath::ATan2(convYrecalc,convXrecalc);
1684  if (alpha_r<0) alpha_r += TMath::Pi()*2;
1685  Int_t sec_r = alpha_r/(TMath::Pi()/9);
1686  alpha_r = (10.+sec_r*20.)*TMath::DegToRad();
1687  Double_t cs_r = TMath::Cos(alpha_r);
1688  Double_t sn_r = TMath::Sin(alpha_r);
1689  Double_t xsV0_r = convXrecalc*cs_r - convYrecalc*sn_r;
1690  fHistoRviaAlphaRecalc->Fill(xsV0_r);
1691 
1692  if (convR > 80) { // conversion within TPC <-> TPC-only tracks
1693  fHistoImpactParameterStudy->Fill(4);
1694 
1695  for (Int_t it=2;it--;) {
1696  Int_t trId = fCurrentV0->GetIndex(it);
1697  AliESDtrack* tr = fESDEvent->GetTrack(trId);
1698  const TBits& bits = tr->GetTPCClusterMap();
1699  Int_t nConflict = 0;
1700  for (Int_t ic=0;ic<159;ic++) {
1701  if (rTPC[ic]>(xsV0-kTPCMargin)) break;
1702  if (bits.TestBitNumber(ic)){
1703  nConflict++;
1704  fHistoRdiff->Fill(xsV0-rTPC[ic]);
1705  }
1706  if (nConflict>kMaxTPCV0Conflicts) {
1707  fHistoImpactParameterStudy->Fill(5);
1708  RemovedByCausality=kTRUE;
1709  break;
1710  }
1711  }
1712  }
1713  }
1714 
1715  //not passed y or z o causality cut:
1716  Bool_t RemovedByAnyCut=kFALSE;
1717  if(RemovedByZcut||RemovedByYcut||RemovedByCausality){
1718  fHistoImpactParameterStudy->Fill(6);
1719  RemovedByAnyCut=kTRUE;
1720  }
1721 
1722  //Fill tree for further analysis
1723  Float_t posZ;
1724  Float_t posY;
1725  Float_t posX;
1726  Float_t posPt;
1727  Float_t negZ;
1728  Float_t negY;
1729  Float_t negX;
1730  Float_t negPt;
1731  Float_t R;
1732  fImpactParamTree->Branch("posPt",&posPt,"posPt/F");
1733  fImpactParamTree->Branch("posY",&posY,"posY/F");
1734  fImpactParamTree->Branch("R",&R,"R/F");
1735  posZ = positiveTrack->GetZ();
1736  posY = positiveTrack->GetY();
1737  posX = positiveTrack->GetX();
1738  posPt = positiveTrack->Pt();
1739  negZ = negativeTrack->GetZ();
1740  negY = negativeTrack->GetY();
1741  negX = negativeTrack->GetX();
1742  negPt = negativeTrack->Pt();
1743  R=convRrecalc;
1744  fImpactParamTree->Fill();
1745 
1746  // fill impact parameter histograms
1747  fHistoPosTrackImpactParamZ->Fill(posZ);
1748  fHistoPosTrackImpactParamY->Fill(posY);
1749  fHistoPosTrackImpactParamX->Fill(posX);
1750  fHistoPosTrackImpactParamZvsPt->Fill(posPt, posZ);
1751  fHistoPosTrackImpactParamYvsPt->Fill(posPt, posY);
1752  fHistoPosTrackImpactParamXvsPt->Fill(posPt, posX);
1753  fHistoNegTrackImpactParamZ->Fill(negZ);
1754  fHistoNegTrackImpactParamY->Fill(negY);
1755  fHistoNegTrackImpactParamX->Fill(negX);
1756  fHistoNegTrackImpactParamZvsPt->Fill(negPt, negZ);
1757  fHistoNegTrackImpactParamYvsPt->Fill(negPt, negY);
1758  fHistoNegTrackImpactParamXvsPt->Fill(negPt, negX);
1759 
1760  // fill comparison histos before / after new cuts
1761  fHistoPt->Fill(positiveTrack->Pt());
1762  fHistoImpactParamZvsR->Fill(convZrecalc,convRrecalc);
1763  //Float_t *DCAzPhoton;
1764  Float_t DCAzPhoton;
1765  //const AliVVertex *PrimVertex = fInputEvent->GetPrimaryVertex();
1766  AliAODConversionPhoton* fCurrentMother=(AliAODConversionPhoton*)fCurrentMotherKF;
1767  //AliAODConversionPhoton *fCurrentMother=dynamic_cast<AliAODConversionPhoton*>(fCurrentMotherKF);
1768  //fCurrentMotherKF->GetDistanceOfClossetApproachToPrimVtx(PrimVertex, DCAzPhoton);
1769  DCAzPhoton = fCurrentMother->GetDCAzToPrimVtx();
1770  fHistoDCAzPhoton->Fill(DCAzPhoton);
1771  if(!RemovedByAnyCut) {
1772  fHistoPt2->Fill(positiveTrack->Pt());
1773  fHistoImpactParamZvsR2->Fill(convZrecalc,convRrecalc);
1774  fHistoDCAzPhoton2->Fill(DCAzPhoton);
1775  }
1776 
1777  return;
1778 }
1779 
1780 //_________________________________________________________________________________
1781 Bool_t AliV0ReaderV1::CheckVectorOnly(vector<Int_t> &vec, Int_t tobechecked)
1782 {
1783  if(tobechecked > -1)
1784  {
1785  vector<Int_t>::iterator it;
1786  it = find (vec.begin(), vec.end(), tobechecked);
1787  if (it != vec.end()) return true;
1788  else return false;
1789  }
1790  return false;
1791 }
1792 
1793 //_________________________________________________________________________________
1795 {
1796  if(tobechecked > -1)
1797  {
1798  vector<Int_t>::iterator it;
1799  it = find (vec.begin(), vec.end(), tobechecked);
1800  if (it != vec.end()) return true;
1801  else{
1802  vec.push_back(tobechecked);
1803  return false;
1804  }
1805  }
1806  return false;
1807 }
1808 
1809 //________________________________________________________________________
1811 {
1812 
1813 }
1814 
1816  std::iterator<std::bidirectional_iterator_tag, AliConversionPhotonBase>(),
1817  fkData(reader),
1818  fCurrentIndex(position),
1819  fDirection(dir)
1820 {
1821 }
1822 
1824  std::iterator<std::bidirectional_iterator_tag, AliConversionPhotonBase>(other),
1825  fkData(other.fkData),
1827  fDirection(other.fDirection)
1828 {
1829 }
1830 
1832  if(this != &other){
1833  fkData = other.fkData;
1834  fCurrentIndex = other.fCurrentIndex;
1835  fDirection = other.fDirection;
1836  }
1837  return *this;
1838 }
1839 
1841  if(fkData != other.fkData) return true;
1842  return fCurrentIndex != other.fCurrentIndex;
1843 }
1844 
1846  iterator tmp(*this);
1847  operator++();
1848  return tmp;
1849 }
1850 
1853  fCurrentIndex++;
1854  } else {
1855  fCurrentIndex--;
1856  }
1857  return *this;
1858 }
1859 
1861  iterator tmp(*this);
1862  operator--();
1863  return tmp;
1864 }
1865 
1868  fCurrentIndex++;
1869  } else {
1870  fCurrentIndex--;
1871  }
1872  return *this;
1873 }
1874 
1876  return (*fkData)[fCurrentIndex];
1877 }
Int_t IsReaderPerformingRelabeling()
Int_t charge
Direction_t fDirection
Iterator in forward direction.
Definition: AliV0ReaderV1.h:69
void FillAODOutput()
anchored LHC16s pass 1 - jet-jet MC in EPOSLHC
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
anchored LHC16r pass 1 - jet-jet MC in EPOSLHC
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