1)Terminate() method implemented in the frame. Simple examples on what to do with...
[u/mrichter/AliRoot.git] / PWG4 / PartCorrDep / AliAnaPi0EbE.cxx
CommitLineData
1c5acb87 1/**************************************************************************\r
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *\r
3 * *\r
4 * Author: The ALICE Off-line Project. *\r
5 * Contributors are mentioned in the code where appropriate. *\r
6 * *\r
7 * Permission to use, copy, modify and distribute this software and its *\r
8 * documentation strictly for non-commercial purposes iGetEntriesFast(s hereby granted *\r
9 * without fee, provided that the above copyright notice appears in all *\r
10 * copies and that both the copyright notice and this permission notice *\r
11 * appear in the supporting documentation. The authors make no claims *\r
12 * about the suitability of this software for any purpose. It is *\r
13 * provided "as is" without express or implied warranty. *\r
14 **************************************************************************/\r
15/* $Id: AliAnaPi0EbE.cxx 28688 2008-09-11 15:04:07Z gconesab $ */\r
16\r
17//_________________________________________________________________________\r
18// Class for the analysis of high pT pi0 event by event\r
19// Pi0 identified by one of the following:\r
20// -Invariant mass of 2 cluster in calorimeter\r
21// -Shower shape analysis in calorimeter\r
22// -Invariant mass of one cluster in calorimeter and one photon reconstructed in TPC (in near future)\r
23//\r
24// -- Author: Gustavo Conesa (LNF-INFN) & Raphaelle Ichou (SUBATECH)\r
25//////////////////////////////////////////////////////////////////////////////\r
26 \r
27 \r
28// --- ROOT system --- \r
29#include <TList.h>\r
30//#include "Riostream.h"\r
31\r
32// --- Analysis system --- \r
33#include "AliAnaPi0EbE.h" \r
34#include "AliLog.h"\r
35#include "AliCaloTrackReader.h"\r
36#include "AliIsolationCut.h"\r
37#include "AliNeutralMesonSelection.h"\r
38#include "AliCaloPID.h"\r
6639984f 39#include "AliMCAnalysisUtils.h"\r
1c5acb87 40#include "AliAODPWG4ParticleCorrelation.h"\r
6639984f 41 #include "AliStack.h"\r
1c5acb87 42ClassImp(AliAnaPi0EbE)\r
43 \r
44//____________________________________________________________________________\r
45AliAnaPi0EbE::AliAnaPi0EbE() : \r
46AliAnaPartCorrBaseClass(), fAnaType(kIMCalo),fCalorimeter(""),\r
47fMinDist(0.),fMinDist2(0.),fMinDist3(0.), \r
48fInputAODGammaConv(0x0),fInputAODGammaConvName(""),\r
49fhPtPi0(0),fhPhiPi0(0),fhEtaPi0(0), \r
50fhPtMCNoPi0(0),fhPhiMCNoPi0(0),fhEtaMCNoPi0(0), \r
51fhPtMCPi0(0),fhPhiMCPi0(0),fhEtaMCPi0(0)\r
52{\r
53 //default ctor\r
54 \r
55 //Initialize parameters\r
56 InitParameters();\r
57 \r
58}\r
59\r
60//____________________________________________________________________________\r
61AliAnaPi0EbE::AliAnaPi0EbE(const AliAnaPi0EbE & p) : \r
62AliAnaPartCorrBaseClass(p), fAnaType(p.fAnaType), fCalorimeter(p.fCalorimeter),\r
63fMinDist(p.fMinDist),fMinDist2(p.fMinDist2), fMinDist3(p.fMinDist3),\r
64fInputAODGammaConv(new TClonesArray(*p.fInputAODGammaConv)),\r
65fInputAODGammaConvName(p.fInputAODGammaConvName),\r
66fhPtPi0(p.fhPtPi0),fhPhiPi0(p.fhPhiPi0),fhEtaPi0(p.fhEtaPi0), \r
67fhPtMCNoPi0(p.fhPtMCNoPi0),fhPhiMCNoPi0(p.fhPhiMCNoPi0),fhEtaMCNoPi0(p.fhEtaMCNoPi0), \r
68fhPtMCPi0(p.fhPtMCPi0),fhPhiMCPi0(p.fhPhiMCPi0),fhEtaMCPi0(p.fhEtaMCPi0) \r
69{\r
70 // cpy ctor\r
71}\r
72\r
73//_________________________________________________________________________\r
74AliAnaPi0EbE & AliAnaPi0EbE::operator = (const AliAnaPi0EbE & p)\r
75{\r
76 // assignment operator\r
77 \r
78 if(&p == this) return *this;\r
79 \r
80 fAnaType = p.fAnaType ;\r
81 fCalorimeter = p.fCalorimeter ;\r
82 fMinDist = p.fMinDist;\r
83 fMinDist2 = p.fMinDist2;\r
84 fMinDist3 = p.fMinDist3;\r
85 \r
86 fInputAODGammaConv = new TClonesArray(*p.fInputAODGammaConv) ;\r
87 fInputAODGammaConvName = p.fInputAODGammaConvName ;\r
88 \r
89 fhPtPi0 = p.fhPtPi0; fhPhiPi0 = p.fhPhiPi0; fhEtaPi0 = p.fhEtaPi0; \r
90 fhPtMCNoPi0 = p.fhPtMCNoPi0; fhPhiMCNoPi0 = p.fhPhiMCNoPi0; fhEtaMCNoPi0 = p.fhEtaMCPi0; \r
91 fhPtMCPi0 = p.fhPtMCPi0; fhPhiMCPi0 = p.fhPhiMCPi0; fhEtaMCPi0 = p.fhEtaMCPi0; \r
92 \r
93 return *this;\r
94 \r
95}\r
96\r
97//____________________________________________________________________________\r
98AliAnaPi0EbE::~AliAnaPi0EbE() \r
99{\r
100 //dtor\r
101 if(fInputAODGammaConv){\r
102 fInputAODGammaConv->Clear() ; \r
103 delete fInputAODGammaConv ;\r
104 }\r
105}\r
106\r
107//________________________________________________________________________\r
108TList * AliAnaPi0EbE::GetCreateOutputObjects()\r
109{ \r
110 // Create histograms to be saved in output file and \r
111 // store them in outputContainer\r
112 TList * outputContainer = new TList() ; \r
113 outputContainer->SetName("Pi0EbEHistos") ; \r
114 \r
115 Int_t nptbins = GetHistoNPtBins();\r
116 Int_t nphibins = GetHistoNPhiBins();\r
117 Int_t netabins = GetHistoNEtaBins();\r
118 Float_t ptmax = GetHistoPtMax();\r
119 Float_t phimax = GetHistoPhiMax();\r
120 Float_t etamax = GetHistoEtaMax();\r
121 Float_t ptmin = GetHistoPtMin();\r
122 Float_t phimin = GetHistoPhiMin();\r
123 Float_t etamin = GetHistoEtaMin(); \r
124 \r
125 fhPtPi0 = new TH1F("hPtPi0","Number of identified #pi^{0} decay",nptbins,ptmin,ptmax); \r
126 fhPtPi0->SetYTitle("N");\r
127 fhPtPi0->SetXTitle("p_{T #pi^{0}}(GeV/c)");\r
128 outputContainer->Add(fhPtPi0) ; \r
129 \r
130 fhPhiPi0 = new TH2F\r
131 ("hPhiPi0","#phi_{#pi^{0}}",nptbins,ptmin,ptmax,nphibins,phimin,phimax); \r
132 fhPhiPi0->SetYTitle("#phi");\r
133 fhPhiPi0->SetXTitle("p_{T #pi^{0}} (GeV/c)");\r
134 outputContainer->Add(fhPhiPi0) ; \r
135 \r
136 fhEtaPi0 = new TH2F\r
137 ("hEtaPi0","#phi_{#pi^{0}}",nptbins,ptmin,ptmax,netabins,etamin,etamax); \r
138 fhEtaPi0->SetYTitle("#eta");\r
139 fhEtaPi0->SetXTitle("p_{T #pi^{0}} (GeV/c)");\r
140 outputContainer->Add(fhEtaPi0) ;\r
141 \r
142 if(IsDataMC()){\r
143 \r
144 fhPtMCPi0 = new TH1F("hPtMCPi0","Identified pi0 from pi0",nptbins,ptmin,ptmax); \r
145 fhPtMCPi0->SetYTitle("N");\r
146 fhPtMCPi0->SetXTitle("p_{T #pi^{0}}(GeV/c)");\r
147 outputContainer->Add(fhPtMCPi0) ; \r
148 \r
149 fhPhiMCPi0 = new TH2F\r
150 ("hPhiMCPi0","Identified pi0 from pi0",nptbins,ptmin,ptmax,nphibins,phimin,phimax); \r
151 fhPhiMCPi0->SetYTitle("#phi");\r
152 fhPhiMCPi0->SetXTitle("p_{T #pi^{0}} (GeV/c)");\r
153 outputContainer->Add(fhPhiMCPi0) ; \r
154 \r
155 fhEtaMCPi0 = new TH2F\r
156 ("hEtaMCPi0","Identified pi0 from pi0",nptbins,ptmin,ptmax,netabins,etamin,etamax); \r
157 fhEtaMCPi0->SetYTitle("#eta");\r
158 fhEtaMCPi0->SetXTitle("p_{T #pi^{0}} (GeV/c)");\r
159 outputContainer->Add(fhEtaMCPi0) ;\r
160 \r
161 fhPtMCNoPi0 = new TH1F("hPtMCNoPi0","Identified pi0 not from pi0",nptbins,ptmin,ptmax); \r
162 fhPtMCNoPi0->SetYTitle("N");\r
163 fhPtMCNoPi0->SetXTitle("p_{T #pi^{0}}(GeV/c)");\r
164 outputContainer->Add(fhPtMCNoPi0) ; \r
165 \r
166 fhPhiMCNoPi0 = new TH2F\r
167 ("hPhiMCNoPi0","Identified pi0 not from pi0",nptbins,ptmin,ptmax,nphibins,phimin,phimax); \r
168 fhPhiMCNoPi0->SetYTitle("#phi");\r
169 fhPhiMCNoPi0->SetXTitle("p_{T #pi^{0}} (GeV/c)");\r
170 outputContainer->Add(fhPhiMCNoPi0) ; \r
171 \r
172 fhEtaMCNoPi0 = new TH2F\r
173 ("hEtaMCNoPi0","Identified pi0 not from pi0",nptbins,ptmin,ptmax,netabins,etamin,etamax); \r
174 fhEtaMCNoPi0->SetYTitle("#eta");\r
175 fhEtaMCNoPi0->SetXTitle("p_{T #pi^{0}} (GeV/c)");\r
176 outputContainer->Add(fhEtaMCNoPi0) ;\r
177 \r
178 }//Histos with MC\r
179 \r
180 \r
181 //Keep neutral meson selection histograms if requiered\r
182 //Setting done in AliNeutralMesonSelection\r
183 \r
184 if(fAnaType!=kSSCalo && GetNeutralMesonSelection()){\r
185 \r
186 TList * nmsHistos = GetNeutralMesonSelection()->GetCreateOutputObjects() ;\r
187 if(GetNeutralMesonSelection()->AreNeutralMesonSelectionHistosKept())\r
188 for(Int_t i = 0; i < nmsHistos->GetEntries(); i++) outputContainer->Add(nmsHistos->At(i)) ;\r
189 }\r
190 \r
191 //Save parameters used for analysis\r
192 TString parList ; //this will be list of parameters used for this analysis.\r
193 char onePar[255] ;\r
194 \r
195 sprintf(onePar,"--- AliAnaPi0EbE ---\n") ;\r
196 parList+=onePar ; \r
197 sprintf(onePar,"fAnaType=%d (Pi0 selection type) \n",fAnaType) ;\r
198 parList+=onePar ;\r
199 \r
200 if(fAnaType == kSSCalo){\r
201 sprintf(onePar,"Calorimeter: %s\n",fCalorimeter.Data()) ;\r
202 parList+=onePar ;\r
203 sprintf(onePar,"fMinDist =%2.2f (Minimal distance to bad channel to accept cluster) \n",fMinDist) ;\r
204 parList+=onePar ;\r
205 sprintf(onePar,"fMinDist2=%2.2f (Cuts on Minimal distance to study acceptance evaluation) \n",fMinDist2) ;\r
206 parList+=onePar ;\r
207 sprintf(onePar,"fMinDist3=%2.2f (One more cut on distance used for acceptance-efficiency study) \n",fMinDist3) ;\r
208 parList+=onePar ;\r
209 }\r
210 \r
211 //Get parameters set in base class.\r
212 parList += GetBaseParametersList() ;\r
213 \r
214 //Get parameters set in PID class.\r
215 if(fAnaType == kSSCalo) parList += GetCaloPID()->GetPIDParametersList() ;\r
216 \r
217 TObjString *oString= new TObjString(parList) ;\r
218 outputContainer->Add(oString);\r
219 \r
220 return outputContainer ;\r
221 \r
222}\r
223\r
224//__________________________________________________________________\r
225void AliAnaPi0EbE::MakeAnalysisFillAOD() \r
226{\r
227 //Do analysis and fill aods\r
228 \r
229 switch(fAnaType) \r
230 {\r
231 case kIMCalo:\r
232 MakeInvMassInCalorimeter();\r
233 break;\r
234 \r
235 case kSSCalo:\r
236 MakeShowerShapeIdentification();\r
237 break;\r
238 \r
239 case kIMCaloTracks:\r
240 MakeInvMassInCalorimeterAndCTS();\r
241 break;\r
242 \r
243 }\r
244}\r
245\r
246//__________________________________________________________________\r
247void AliAnaPi0EbE::MakeInvMassInCalorimeter() \r
248{\r
249 //Do analysis and fill aods\r
250 //Search for the photon decay in calorimeters\r
251 //Read photon list from AOD, produced in class AliAnaPhoton\r
252 //Check if 2 photons have the mass of the pi0.\r
253 \r
254 TLorentzVector mom1;\r
255 TLorentzVector mom2;\r
256 TLorentzVector mom ;\r
257 Int_t tag1 =-1;\r
258 Int_t tag2 =-1;\r
259 Int_t tag = -1;\r
260 \r
261 if(!GetInputAODBranch())\r
262 AliFatal(Form("MakeInvMassInCalo: No input calo photons in AOD with name branch < %s > \n",GetInputAODName().Data()));\r
263 \r
264 for(Int_t iphoton = 0; iphoton < GetInputAODBranch()->GetEntriesFast(); iphoton++){\r
265 AliAODPWG4Particle * photon1 = (AliAODPWG4Particle*) (GetInputAODBranch()->At(iphoton));\r
266 mom1 = *(photon1->Momentum());\r
267 \r
6639984f 268 \r
1c5acb87 269 for(Int_t jphoton = iphoton+1; jphoton < GetInputAODBranch()->GetEntriesFast()-1; jphoton++){\r
270 \r
271 AliAODPWG4ParticleCorrelation * photon2 = (AliAODPWG4ParticleCorrelation*) (GetInputAODBranch()->At(jphoton));\r
272 mom2 = *(photon2->Momentum());\r
6639984f 273 \r
1c5acb87 274 //Select good pair (good phi, pt cuts, aperture and invariant mass)\r
275 if(GetNeutralMesonSelection()->SelectPair(mom1, mom2))\r
276 {\r
6639984f 277 if(GetDebug()>1) \r
278 printf("Selected gamma pair: pt %f, phi %f, eta%f \n",(mom1+mom2).Pt(), (mom1+mom2).Phi()*180./3.1416, (mom1+mom2).Eta());\r
1c5acb87 279 \r
6639984f 280 //Play with the MC stack if available\r
1c5acb87 281 if(IsDataMC()){\r
282 //Check origin of the candidates\r
6639984f 283 tag1 = GetMCAnalysisUtils()->CheckOrigin(photon1->GetLabel(), GetMCStack());\r
284 tag2 = GetMCAnalysisUtils()->CheckOrigin(photon2->GetLabel(), GetMCStack());\r
285 \r
1c5acb87 286 if(GetDebug() > 0) printf("Origin of: photon1 %d; photon2 %d \n",tag1, tag2);\r
6639984f 287 if(tag1 == AliMCAnalysisUtils::kMCPi0Decay && tag2 == AliMCAnalysisUtils::kMCPi0Decay){\r
288 \r
289 //Check if pi0 mother is the same\r
290 Int_t label1 = photon1->GetLabel();\r
291 TParticle * mother1 = GetMCStack()->Particle(label1);//photon in kine tree\r
292 label1 = mother1->GetFirstMother();\r
293 //mother1 = GetMCStack()->Particle(label1);//pi0\r
294 \r
295 Int_t label2 = photon2->GetLabel();\r
296 TParticle * mother2 = GetMCStack()->Particle(label2);//photon in kine tree\r
297 label2 = mother2->GetFirstMother();\r
298 //mother2 = GetMCStack()->Particle(label2);//pi0\r
299 \r
300 //printf("mother1 %d, mother2 %d\n",label1,label2);\r
301 if(label1 == label2)\r
302 tag = AliMCAnalysisUtils::kMCPi0;\r
303 }\r
1c5acb87 304 }//Work with stack also \r
305 \r
306 //Create AOD for analysis\r
307 mom = mom1+mom2;\r
308 AliAODPWG4Particle pi0 = AliAODPWG4Particle(mom);\r
309 //pi0.SetLabel(calo->GetLabel(0));\r
310 pi0.SetPdg(AliCaloPID::kPi0);\r
311 pi0.SetDetector(photon1->GetDetector());\r
312 pi0.SetTag(tag); \r
313 //Set the indeces of the original caloclusters \r
314 pi0.SetCaloLabel(photon1->GetCaloLabel(0), photon2->GetCaloLabel(0));\r
315 AddAODParticle(pi0);\r
316 }//pi0\r
317 }//2n photon loop\r
318 \r
319 }//1st photon loop\r
320 \r
321 if(GetDebug() > 1) printf("End fill AODs \n"); \r
322 \r
323}\r
324\r
325//__________________________________________________________________\r
326void AliAnaPi0EbE::MakeInvMassInCalorimeterAndCTS() \r
327{\r
328 //Do analysis and fill aods\r
329 //Search for the photon decay in calorimeters\r
330 //Read photon list from AOD, produced in class AliAnaPhoton and AliGammaConversion\r
331 //Check if 2 photons have the mass of the pi0.\r
332 \r
333 TLorentzVector mom1;\r
334 TLorentzVector mom2;\r
335 TLorentzVector mom ;\r
336 Int_t tag1 =-1;\r
337 Int_t tag2 =-1;\r
338 Int_t tag = -1;\r
339 \r
340 if(!GetInputAODBranch())\r
341 AliFatal(Form("MakeInvMassInCalo: No input calo photons in AOD branch with name < %s > \n",GetInputAODName().Data()));\r
342 \r
343 for(Int_t iphoton = 0; iphoton < GetInputAODBranch()->GetEntriesFast(); iphoton++){\r
344 AliAODPWG4Particle * photon1 = (AliAODPWG4Particle*) (GetInputAODBranch()->At(iphoton));\r
345 mom1 = *(photon1->Momentum());\r
346 \r
347 //Play with the MC stack if available\r
348 fInputAODGammaConv = (TClonesArray *) GetReader()->GetAOD()->FindListObject(fInputAODGammaConvName);\r
349 if(!fInputAODGammaConv) \r
350 AliFatal(Form("MakeInvMassInCaloAndCTS: No input gamma conversions in AOD branch with name < %s >",fInputAODGammaConvName.Data())); \r
351 \r
352 for(Int_t jphoton = iphoton+1; jphoton < fInputAODGammaConv->GetEntriesFast()-1; jphoton++){\r
353 AliAODPWG4ParticleCorrelation * photon2 = (AliAODPWG4ParticleCorrelation*) (fInputAODGammaConv->At(jphoton));\r
354 mom2 = *(photon2->Momentum());\r
355 //Select good pair (good phi, pt cuts, aperture and invariant mass)\r
356 if(GetNeutralMesonSelection()->SelectPair(mom1, mom2)){\r
357 if(GetDebug() > 1) printf("Selected gamma pair: pt %f, phi %f, eta%f\n",(mom1+mom2).Pt(), (mom1+mom2).Phi()*180./3.1416, (mom1+mom2).Eta());\r
358 \r
359 if(IsDataMC()){\r
360 //Check origin of the candidates\r
6639984f 361 tag1 = GetMCAnalysisUtils()->CheckOrigin(photon1->GetLabel(), GetMCStack());\r
362 tag2 = GetMCAnalysisUtils()->CheckOrigin(photon2->GetLabel(), GetMCStack());\r
1c5acb87 363 if(GetDebug() > 0) printf("Origin of: photon1 %d; photon2 %d \n",tag1, tag2);\r
6639984f 364 if(tag1 == AliMCAnalysisUtils::kMCPi0Decay && tag2 == AliMCAnalysisUtils::kMCPi0Decay){\r
365 //Check if pi0 mother is the same\r
366 Int_t label1 = photon1->GetLabel();\r
367 TParticle * mother1 = GetMCStack()->Particle(label1);//photon in kine tree\r
368 label1 = mother1->GetFirstMother();\r
369 //mother1 = GetMCStack()->Particle(label1);//pi0\r
370 \r
371 Int_t label2 = photon2->GetLabel();\r
372 TParticle * mother2 = GetMCStack()->Particle(label2);//photon in kine tree\r
373 label2 = mother2->GetFirstMother();\r
374 //mother2 = GetMCStack()->Particle(label2);//pi0\r
375 \r
376 //printf("mother1 %d, mother2 %d\n",label1,label2);\r
377 if(label1 == label2)\r
378 tag = AliMCAnalysisUtils::kMCPi0;\r
379 }\r
1c5acb87 380 }//Work with stack also \r
381 \r
382 //Create AOD for analysis\r
383 mom = mom1+mom2;\r
384 AliAODPWG4ParticleCorrelation pi0 = AliAODPWG4ParticleCorrelation(mom);\r
385 //pi0.SetLabel(calo->GetLabel(0));\r
386 pi0.SetPdg(AliCaloPID::kPi0);\r
387 pi0.SetDetector(photon1->GetDetector());\r
388 pi0.SetTag(tag);\r
389 //Set the indeces of the original tracks or caloclusters \r
390 pi0.SetCaloLabel(photon1->GetCaloLabel(0), -1);\r
391 pi0.SetTrackLabel(photon2->GetTrackLabel(0), photon2->GetTrackLabel(1));\r
392 AddAODParticle(pi0);\r
393 }//pi0\r
394 }//2n photon loop\r
395 \r
396 }//1st photon loop\r
397 \r
398 if(GetDebug() > 1) printf("End fill AODs \n"); \r
399 \r
400}\r
401\r
402\r
403//__________________________________________________________________\r
404void AliAnaPi0EbE::MakeShowerShapeIdentification() \r
405{\r
406 //Search for pi0 in fCalorimeter with shower shape analysis \r
407 \r
408 TClonesArray * pl = new TClonesArray; \r
409 \r
410 //Get vertex for photon momentum calculation\r
411 Double_t vertex[]={0,0,0} ; //vertex ;\r
412 if(!GetReader()->GetDataType()== AliCaloTrackReader::kMC) GetReader()->GetVertex(vertex);\r
413 \r
414 //Select the Calorimeter of the photon\r
415 if(fCalorimeter == "PHOS")\r
416 pl = GetAODPHOS();\r
417 else if (fCalorimeter == "EMCAL")\r
418 pl = GetAODEMCAL();\r
419 //Fill AODCaloClusters and AODParticle with PHOS aods\r
420 TLorentzVector mom ;\r
421 for(Int_t icalo = 0; icalo < pl->GetEntriesFast(); icalo++){\r
422 AliAODCaloCluster * calo = (AliAODCaloCluster*) (pl->At(icalo)); \r
423 \r
424 //Cluster selection, not charged, with photon id and in fidutial cut\r
425 //Get Momentum vector, \r
426 calo->GetMomentum(mom,vertex);//Assume that come from vertex in straight line\r
427 //If too small or big pt, skip it\r
428 if(mom.Pt() < GetMinPt() || mom.Pt() > GetMaxPt() ) continue ; \r
429 //Check acceptance selection\r
430 if(IsFidutialCutOn()){\r
431 Bool_t in = GetFidutialCut()->IsInFidutialCut(mom,fCalorimeter) ;\r
432 if(! in ) continue ;\r
433 }\r
434 \r
435 //Create AOD for analysis\r
436 AliAODPWG4Particle aodph = AliAODPWG4Particle(mom);\r
437 aodph.SetLabel(calo->GetLabel(0));\r
438 aodph.SetDetector(fCalorimeter);\r
439 if(GetDebug() > 1) printf("AliAnaPhoton::FillAOD: Min pt cut and fidutial cut passed: pt %3.2f, phi %2.2f, eta %1.2f\n",aodph.Pt(),aodph.Phi(),aodph.Eta()); \r
440 \r
441 //Check Distance to Bad channel, set bit.\r
442 Double_t distBad=calo->GetDistToBadChannel() ; //Distance to bad channel\r
443 if(distBad < 0.) distBad=9999. ; //workout strange convension dist = -1. ;\r
444 if(distBad < fMinDist) //In bad channel (PHOS cristal size 2.2x2.2 cm)\r
445 continue ;\r
446 \r
447 if(GetDebug() > 1) printf("AliAnaPhoton::FillAOD: Bad channel cut passed %4.2f\n",distBad);\r
448 \r
449 if(distBad > fMinDist3) aodph.SetDistToBad(2) ;\r
450 else if(distBad > fMinDist2) aodph.SetDistToBad(1) ; \r
451 else aodph.SetDistToBad(0) ;\r
452 \r
453 //Check PID\r
454 //PID selection or bit setting\r
455 if(GetReader()->GetDataType() == AliCaloTrackReader::kMC){\r
456 //Get most probable PID, check PID weights (in MC this option is mandatory)\r
457 aodph.SetPdg(GetCaloPID()->GetPdg(fCalorimeter,calo->PID(),mom.E()));//PID with weights\r
458 if(GetDebug() > 1) printf("AliAnaPhoton::FillAOD: PDG of identified particle %d\n",aodph.GetPdg());\r
459 //If primary is not photon, skip it.\r
460 if(aodph.GetPdg() != AliCaloPID::kPhoton) continue ;\r
461 } \r
462 else if(IsCaloPIDOn()){\r
463 //Skip matched clusters with tracks\r
464 if(calo->GetNTracksMatched() > 0) continue ;\r
465 \r
466 //Get most probable PID, 2 options check PID weights \r
467 //or redo PID, recommended option for EMCal. \r
468 if(!IsCaloPIDRecalculationOn())\r
469 aodph.SetPdg(GetCaloPID()->GetPdg(fCalorimeter,calo->PID(),mom.E()));//PID with weights\r
470 else\r
471 aodph.SetPdg(GetCaloPID()->GetPdg(fCalorimeter,mom,calo));//PID recalculated\r
472 \r
473 if(GetDebug() > 1) printf("AliAnaPhoton::FillAOD: PDG of identified particle %d\n",aodph.GetPdg());\r
474 \r
475 //If cluster does not pass pid, not photon, skip it.\r
476 if(aodph.GetPdg() != AliCaloPID::kPhoton) continue ; \r
477 \r
478 }\r
479 else{\r
480 //Set PID bits for later selection (AliAnaPi0 for example)\r
481 //GetPDG already called in SetPIDBits.\r
482 GetCaloPID()->SetPIDBits(fCalorimeter,calo,&aodph);\r
483 if(GetDebug() > 1) printf("AliAnaPhoton::FillAOD: PID Bits set \n"); \r
484 }\r
485 \r
486 if(GetDebug() > 1) printf("AliAnaPhoton::FillAOD: Photon selection cuts passed: pT %3.2f, pdg %d\n",aodph.Pt(), aodph.GetPdg());\r
487 \r
488 //Play with the MC stack if available\r
489 //Check origin of the candidates\r
490 if(IsDataMC()){\r
6639984f 491 aodph.SetTag(GetMCAnalysisUtils()->CheckOrigin(calo->GetLabel(0),GetMCStack()));\r
1c5acb87 492 if(GetDebug() > 0) printf("AliAnaPi0EbE::FillAOD: Origin of candidate %d\n",aodph.GetTag());\r
493 }//Work with stack also \r
494 \r
495 //Add AOD with photon object to aod branch\r
496 AddAODParticle(aodph);\r
497 \r
498 }//loop\r
499 \r
500 if(GetDebug() > 1) printf("AliAnaPhoton::FillAOD: End fill AODs \n"); \r
501 \r
502}\r
503//__________________________________________________________________\r
504void AliAnaPi0EbE::MakeAnalysisFillHistograms() \r
505{\r
506 //Do analysis and fill histograms\r
507 \r
508 if(!GetOutputAODBranch())\r
509 AliFatal(Form("MakeInvMassInCalo: No output pi0 in AOD branch with name < %s > \n",GetOutputAODName().Data()));\r
510 \r
511 //Loop on stored AOD pi0\r
512 Int_t naod = GetOutputAODBranch()->GetEntriesFast();\r
513 if(GetDebug() > 0) printf("pi0 aod branch entries %d\n", naod);\r
514 \r
515 for(Int_t iaod = 0; iaod < naod ; iaod++){\r
516 \r
517 AliAODPWG4Particle* pi0 = (AliAODPWG4Particle*) (GetOutputAODBranch()->At(iaod));\r
518 Int_t pdg = pi0->GetPdg();\r
519 \r
520 if(pdg != AliCaloPID::kPi0) continue; \r
521 \r
522 //Fill pi0 histograms \r
523 Float_t pt = pi0->Pt();\r
524 Float_t phi = pi0->Phi();\r
525 Float_t eta = pi0->Eta();\r
526 \r
527 \r
528 fhPtPi0 ->Fill(pt);\r
529 fhPhiPi0 ->Fill(pt,phi);\r
530 fhEtaPi0 ->Fill(pt,eta);\r
531 \r
532 if(IsDataMC()){\r
6639984f 533 if(pi0->GetTag()== AliMCAnalysisUtils::kMCPi0){\r
1c5acb87 534 fhPtMCPi0 ->Fill(pt);\r
535 fhPhiMCPi0 ->Fill(pt,phi);\r
536 fhEtaMCPi0 ->Fill(pt,eta);\r
537 }\r
538 else{\r
539 fhPtMCNoPi0 ->Fill(pt);\r
540 fhPhiMCNoPi0 ->Fill(pt,phi);\r
541 fhEtaMCNoPi0 ->Fill(pt,eta);\r
542 }\r
543 }//Histograms with MC\r
544 \r
545 }// aod loop\r
546 \r
547}\r
548\r
549\r
550//____________________________________________________________________________\r
551void AliAnaPi0EbE::InitParameters()\r
552{\r
553 \r
554 //Initialize the parameters of the analysis.\r
555 SetOutputAODClassName("AliAODPWG4Particle");\r
556 SetOutputAODName("pi0s");\r
557 fInputAODGammaConvName = "gammaconv" ;\r
558 fAnaType = kIMCalo ;\r
559 fCalorimeter = "PHOS" ;\r
560 fMinDist = 2.;\r
561 fMinDist2 = 4.;\r
562 fMinDist3 = 5.;\r
563 \r
564}\r
565\r
566//__________________________________________________________________\r
567void AliAnaPi0EbE::Print(const Option_t * opt) const\r
568{\r
569 \r
570 //Print some relevant parameters set for the analysis\r
571 if(! opt)\r
572 return;\r
573 \r
574 printf("**** Print %s %s ****\n", GetName(), GetTitle() ) ;\r
575 AliAnaPartCorrBaseClass::Print("");\r
576 printf("Analysis Type = %d \n", fAnaType) ;\r
577 if(fAnaType == kSSCalo){ \r
578 printf("Calorimeter = %s\n", fCalorimeter.Data()) ;\r
579 printf("Min Distance to Bad Channel = %2.1f\n",fMinDist);\r
580 printf("Min Distance to Bad Channel 2 = %2.1f\n",fMinDist2);\r
581 printf("Min Distance to Bad Channel 3 = %2.1f\n",fMinDist3); \r
582 } \r
583 printf(" \n") ;\r
584 \r
585} \r