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