1 /**************************************************************************
\r
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
\r
4 * Author: The ALICE Off-line Project. *
\r
5 * Contributors are mentioned in the code where appropriate. *
\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
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
24 // -- Author: Gustavo Conesa (LNF-INFN) & Raphaelle Ichou (SUBATECH)
\r
25 //////////////////////////////////////////////////////////////////////////////
\r
28 // --- ROOT system ---
\r
30 //#include "Riostream.h"
\r
32 // --- Analysis system ---
\r
33 #include "AliAnaPi0EbE.h"
\r
35 #include "AliCaloTrackReader.h"
\r
36 #include "AliIsolationCut.h"
\r
37 #include "AliNeutralMesonSelection.h"
\r
38 #include "AliCaloPID.h"
\r
39 #include "AliAODPWG4ParticleCorrelation.h"
\r
41 ClassImp(AliAnaPi0EbE)
\r
43 //____________________________________________________________________________
\r
44 AliAnaPi0EbE::AliAnaPi0EbE() :
\r
45 AliAnaPartCorrBaseClass(), fAnaType(kIMCalo),fCalorimeter(""),
\r
46 fMinDist(0.),fMinDist2(0.),fMinDist3(0.),
\r
47 fInputAODGammaConv(0x0),fInputAODGammaConvName(""),
\r
48 fhPtPi0(0),fhPhiPi0(0),fhEtaPi0(0),
\r
49 fhPtMCNoPi0(0),fhPhiMCNoPi0(0),fhEtaMCNoPi0(0),
\r
50 fhPtMCPi0(0),fhPhiMCPi0(0),fhEtaMCPi0(0)
\r
54 //Initialize parameters
\r
59 //____________________________________________________________________________
\r
60 AliAnaPi0EbE::AliAnaPi0EbE(const AliAnaPi0EbE & p) :
\r
61 AliAnaPartCorrBaseClass(p), fAnaType(p.fAnaType), fCalorimeter(p.fCalorimeter),
\r
62 fMinDist(p.fMinDist),fMinDist2(p.fMinDist2), fMinDist3(p.fMinDist3),
\r
63 fInputAODGammaConv(new TClonesArray(*p.fInputAODGammaConv)),
\r
64 fInputAODGammaConvName(p.fInputAODGammaConvName),
\r
65 fhPtPi0(p.fhPtPi0),fhPhiPi0(p.fhPhiPi0),fhEtaPi0(p.fhEtaPi0),
\r
66 fhPtMCNoPi0(p.fhPtMCNoPi0),fhPhiMCNoPi0(p.fhPhiMCNoPi0),fhEtaMCNoPi0(p.fhEtaMCNoPi0),
\r
67 fhPtMCPi0(p.fhPtMCPi0),fhPhiMCPi0(p.fhPhiMCPi0),fhEtaMCPi0(p.fhEtaMCPi0)
\r
72 //_________________________________________________________________________
\r
73 AliAnaPi0EbE & AliAnaPi0EbE::operator = (const AliAnaPi0EbE & p)
\r
75 // assignment operator
\r
77 if(&p == this) return *this;
\r
79 fAnaType = p.fAnaType ;
\r
80 fCalorimeter = p.fCalorimeter ;
\r
81 fMinDist = p.fMinDist;
\r
82 fMinDist2 = p.fMinDist2;
\r
83 fMinDist3 = p.fMinDist3;
\r
85 fInputAODGammaConv = new TClonesArray(*p.fInputAODGammaConv) ;
\r
86 fInputAODGammaConvName = p.fInputAODGammaConvName ;
\r
88 fhPtPi0 = p.fhPtPi0; fhPhiPi0 = p.fhPhiPi0; fhEtaPi0 = p.fhEtaPi0;
\r
89 fhPtMCNoPi0 = p.fhPtMCNoPi0; fhPhiMCNoPi0 = p.fhPhiMCNoPi0; fhEtaMCNoPi0 = p.fhEtaMCPi0;
\r
90 fhPtMCPi0 = p.fhPtMCPi0; fhPhiMCPi0 = p.fhPhiMCPi0; fhEtaMCPi0 = p.fhEtaMCPi0;
\r
96 //____________________________________________________________________________
\r
97 AliAnaPi0EbE::~AliAnaPi0EbE()
\r
100 if(fInputAODGammaConv){
\r
101 fInputAODGammaConv->Clear() ;
\r
102 delete fInputAODGammaConv ;
\r
106 //________________________________________________________________________
\r
107 TList * AliAnaPi0EbE::GetCreateOutputObjects()
\r
109 // Create histograms to be saved in output file and
\r
110 // store them in outputContainer
\r
111 TList * outputContainer = new TList() ;
\r
112 outputContainer->SetName("Pi0EbEHistos") ;
\r
114 Int_t nptbins = GetHistoNPtBins();
\r
115 Int_t nphibins = GetHistoNPhiBins();
\r
116 Int_t netabins = GetHistoNEtaBins();
\r
117 Float_t ptmax = GetHistoPtMax();
\r
118 Float_t phimax = GetHistoPhiMax();
\r
119 Float_t etamax = GetHistoEtaMax();
\r
120 Float_t ptmin = GetHistoPtMin();
\r
121 Float_t phimin = GetHistoPhiMin();
\r
122 Float_t etamin = GetHistoEtaMin();
\r
124 fhPtPi0 = new TH1F("hPtPi0","Number of identified #pi^{0} decay",nptbins,ptmin,ptmax);
\r
125 fhPtPi0->SetYTitle("N");
\r
126 fhPtPi0->SetXTitle("p_{T #pi^{0}}(GeV/c)");
\r
127 outputContainer->Add(fhPtPi0) ;
\r
129 fhPhiPi0 = new TH2F
\r
130 ("hPhiPi0","#phi_{#pi^{0}}",nptbins,ptmin,ptmax,nphibins,phimin,phimax);
\r
131 fhPhiPi0->SetYTitle("#phi");
\r
132 fhPhiPi0->SetXTitle("p_{T #pi^{0}} (GeV/c)");
\r
133 outputContainer->Add(fhPhiPi0) ;
\r
135 fhEtaPi0 = new TH2F
\r
136 ("hEtaPi0","#phi_{#pi^{0}}",nptbins,ptmin,ptmax,netabins,etamin,etamax);
\r
137 fhEtaPi0->SetYTitle("#eta");
\r
138 fhEtaPi0->SetXTitle("p_{T #pi^{0}} (GeV/c)");
\r
139 outputContainer->Add(fhEtaPi0) ;
\r
143 fhPtMCPi0 = new TH1F("hPtMCPi0","Identified pi0 from pi0",nptbins,ptmin,ptmax);
\r
144 fhPtMCPi0->SetYTitle("N");
\r
145 fhPtMCPi0->SetXTitle("p_{T #pi^{0}}(GeV/c)");
\r
146 outputContainer->Add(fhPtMCPi0) ;
\r
148 fhPhiMCPi0 = new TH2F
\r
149 ("hPhiMCPi0","Identified pi0 from pi0",nptbins,ptmin,ptmax,nphibins,phimin,phimax);
\r
150 fhPhiMCPi0->SetYTitle("#phi");
\r
151 fhPhiMCPi0->SetXTitle("p_{T #pi^{0}} (GeV/c)");
\r
152 outputContainer->Add(fhPhiMCPi0) ;
\r
154 fhEtaMCPi0 = new TH2F
\r
155 ("hEtaMCPi0","Identified pi0 from pi0",nptbins,ptmin,ptmax,netabins,etamin,etamax);
\r
156 fhEtaMCPi0->SetYTitle("#eta");
\r
157 fhEtaMCPi0->SetXTitle("p_{T #pi^{0}} (GeV/c)");
\r
158 outputContainer->Add(fhEtaMCPi0) ;
\r
160 fhPtMCNoPi0 = new TH1F("hPtMCNoPi0","Identified pi0 not from pi0",nptbins,ptmin,ptmax);
\r
161 fhPtMCNoPi0->SetYTitle("N");
\r
162 fhPtMCNoPi0->SetXTitle("p_{T #pi^{0}}(GeV/c)");
\r
163 outputContainer->Add(fhPtMCNoPi0) ;
\r
165 fhPhiMCNoPi0 = new TH2F
\r
166 ("hPhiMCNoPi0","Identified pi0 not from pi0",nptbins,ptmin,ptmax,nphibins,phimin,phimax);
\r
167 fhPhiMCNoPi0->SetYTitle("#phi");
\r
168 fhPhiMCNoPi0->SetXTitle("p_{T #pi^{0}} (GeV/c)");
\r
169 outputContainer->Add(fhPhiMCNoPi0) ;
\r
171 fhEtaMCNoPi0 = new TH2F
\r
172 ("hEtaMCNoPi0","Identified pi0 not from pi0",nptbins,ptmin,ptmax,netabins,etamin,etamax);
\r
173 fhEtaMCNoPi0->SetYTitle("#eta");
\r
174 fhEtaMCNoPi0->SetXTitle("p_{T #pi^{0}} (GeV/c)");
\r
175 outputContainer->Add(fhEtaMCNoPi0) ;
\r
180 //Keep neutral meson selection histograms if requiered
\r
181 //Setting done in AliNeutralMesonSelection
\r
183 if(fAnaType!=kSSCalo && GetNeutralMesonSelection()){
\r
185 TList * nmsHistos = GetNeutralMesonSelection()->GetCreateOutputObjects() ;
\r
186 if(GetNeutralMesonSelection()->AreNeutralMesonSelectionHistosKept())
\r
187 for(Int_t i = 0; i < nmsHistos->GetEntries(); i++) outputContainer->Add(nmsHistos->At(i)) ;
\r
190 //Save parameters used for analysis
\r
191 TString parList ; //this will be list of parameters used for this analysis.
\r
194 sprintf(onePar,"--- AliAnaPi0EbE ---\n") ;
\r
196 sprintf(onePar,"fAnaType=%d (Pi0 selection type) \n",fAnaType) ;
\r
199 if(fAnaType == kSSCalo){
\r
200 sprintf(onePar,"Calorimeter: %s\n",fCalorimeter.Data()) ;
\r
202 sprintf(onePar,"fMinDist =%2.2f (Minimal distance to bad channel to accept cluster) \n",fMinDist) ;
\r
204 sprintf(onePar,"fMinDist2=%2.2f (Cuts on Minimal distance to study acceptance evaluation) \n",fMinDist2) ;
\r
206 sprintf(onePar,"fMinDist3=%2.2f (One more cut on distance used for acceptance-efficiency study) \n",fMinDist3) ;
\r
210 //Get parameters set in base class.
\r
211 parList += GetBaseParametersList() ;
\r
213 //Get parameters set in PID class.
\r
214 if(fAnaType == kSSCalo) parList += GetCaloPID()->GetPIDParametersList() ;
\r
216 TObjString *oString= new TObjString(parList) ;
\r
217 outputContainer->Add(oString);
\r
219 return outputContainer ;
\r
223 //__________________________________________________________________
\r
224 void AliAnaPi0EbE::MakeAnalysisFillAOD()
\r
226 //Do analysis and fill aods
\r
231 MakeInvMassInCalorimeter();
\r
235 MakeShowerShapeIdentification();
\r
238 case kIMCaloTracks:
\r
239 MakeInvMassInCalorimeterAndCTS();
\r
245 //__________________________________________________________________
\r
246 void AliAnaPi0EbE::MakeInvMassInCalorimeter()
\r
248 //Do analysis and fill aods
\r
249 //Search for the photon decay in calorimeters
\r
250 //Read photon list from AOD, produced in class AliAnaPhoton
\r
251 //Check if 2 photons have the mass of the pi0.
\r
253 TLorentzVector mom1;
\r
254 TLorentzVector mom2;
\r
255 TLorentzVector mom ;
\r
260 if(!GetInputAODBranch())
\r
261 AliFatal(Form("MakeInvMassInCalo: No input calo photons in AOD with name branch < %s > \n",GetInputAODName().Data()));
\r
263 for(Int_t iphoton = 0; iphoton < GetInputAODBranch()->GetEntriesFast(); iphoton++){
\r
264 AliAODPWG4Particle * photon1 = (AliAODPWG4Particle*) (GetInputAODBranch()->At(iphoton));
\r
265 mom1 = *(photon1->Momentum());
\r
267 //Play with the MC stack if available
\r
269 for(Int_t jphoton = iphoton+1; jphoton < GetInputAODBranch()->GetEntriesFast()-1; jphoton++){
\r
271 AliAODPWG4ParticleCorrelation * photon2 = (AliAODPWG4ParticleCorrelation*) (GetInputAODBranch()->At(jphoton));
\r
272 mom2 = *(photon2->Momentum());
\r
273 //Select good pair (good phi, pt cuts, aperture and invariant mass)
\r
274 if(GetNeutralMesonSelection()->SelectPair(mom1, mom2))
\r
276 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
280 //Check origin of the candidates
\r
281 tag1 = GetCaloPID()->CheckOrigin(photon1->GetLabel(), GetMCStack());
\r
282 tag2 = GetCaloPID()->CheckOrigin(photon2->GetLabel(), GetMCStack());
\r
283 if(GetDebug() > 0) printf("Origin of: photon1 %d; photon2 %d \n",tag1, tag2);
\r
284 if(tag1 == AliCaloPID::kMCPi0Decay && tag2 == AliCaloPID::kMCPi0Decay) tag = AliCaloPID::kMCPi0;
\r
285 }//Work with stack also
\r
287 //Create AOD for analysis
\r
289 AliAODPWG4Particle pi0 = AliAODPWG4Particle(mom);
\r
290 //pi0.SetLabel(calo->GetLabel(0));
\r
291 pi0.SetPdg(AliCaloPID::kPi0);
\r
292 pi0.SetDetector(photon1->GetDetector());
\r
294 //Set the indeces of the original caloclusters
\r
295 pi0.SetCaloLabel(photon1->GetCaloLabel(0), photon2->GetCaloLabel(0));
\r
296 AddAODParticle(pi0);
\r
302 if(GetDebug() > 1) printf("End fill AODs \n");
\r
306 //__________________________________________________________________
\r
307 void AliAnaPi0EbE::MakeInvMassInCalorimeterAndCTS()
\r
309 //Do analysis and fill aods
\r
310 //Search for the photon decay in calorimeters
\r
311 //Read photon list from AOD, produced in class AliAnaPhoton and AliGammaConversion
\r
312 //Check if 2 photons have the mass of the pi0.
\r
314 TLorentzVector mom1;
\r
315 TLorentzVector mom2;
\r
316 TLorentzVector mom ;
\r
321 if(!GetInputAODBranch())
\r
322 AliFatal(Form("MakeInvMassInCalo: No input calo photons in AOD branch with name < %s > \n",GetInputAODName().Data()));
\r
324 for(Int_t iphoton = 0; iphoton < GetInputAODBranch()->GetEntriesFast(); iphoton++){
\r
325 AliAODPWG4Particle * photon1 = (AliAODPWG4Particle*) (GetInputAODBranch()->At(iphoton));
\r
326 mom1 = *(photon1->Momentum());
\r
328 //Play with the MC stack if available
\r
329 fInputAODGammaConv = (TClonesArray *) GetReader()->GetAOD()->FindListObject(fInputAODGammaConvName);
\r
330 if(!fInputAODGammaConv)
\r
331 AliFatal(Form("MakeInvMassInCaloAndCTS: No input gamma conversions in AOD branch with name < %s >",fInputAODGammaConvName.Data()));
\r
333 for(Int_t jphoton = iphoton+1; jphoton < fInputAODGammaConv->GetEntriesFast()-1; jphoton++){
\r
334 AliAODPWG4ParticleCorrelation * photon2 = (AliAODPWG4ParticleCorrelation*) (fInputAODGammaConv->At(jphoton));
\r
335 mom2 = *(photon2->Momentum());
\r
336 //Select good pair (good phi, pt cuts, aperture and invariant mass)
\r
337 if(GetNeutralMesonSelection()->SelectPair(mom1, mom2)){
\r
338 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
341 //Check origin of the candidates
\r
342 tag1 = GetCaloPID()->CheckOrigin(photon1->GetLabel(), GetMCStack());
\r
343 tag2 = GetCaloPID()->CheckOrigin(photon2->GetLabel(), GetMCStack());
\r
344 if(GetDebug() > 0) printf("Origin of: photon1 %d; photon2 %d \n",tag1, tag2);
\r
345 if(tag1 == AliCaloPID::kMCPi0Decay && tag2 == AliCaloPID::kMCPi0Decay) tag = AliCaloPID::kMCPi0;
\r
346 }//Work with stack also
\r
348 //Create AOD for analysis
\r
350 AliAODPWG4ParticleCorrelation pi0 = AliAODPWG4ParticleCorrelation(mom);
\r
351 //pi0.SetLabel(calo->GetLabel(0));
\r
352 pi0.SetPdg(AliCaloPID::kPi0);
\r
353 pi0.SetDetector(photon1->GetDetector());
\r
355 //Set the indeces of the original tracks or caloclusters
\r
356 pi0.SetCaloLabel(photon1->GetCaloLabel(0), -1);
\r
357 pi0.SetTrackLabel(photon2->GetTrackLabel(0), photon2->GetTrackLabel(1));
\r
358 AddAODParticle(pi0);
\r
364 if(GetDebug() > 1) printf("End fill AODs \n");
\r
369 //__________________________________________________________________
\r
370 void AliAnaPi0EbE::MakeShowerShapeIdentification()
\r
372 //Search for pi0 in fCalorimeter with shower shape analysis
\r
374 TClonesArray * pl = new TClonesArray;
\r
376 //Get vertex for photon momentum calculation
\r
377 Double_t vertex[]={0,0,0} ; //vertex ;
\r
378 if(!GetReader()->GetDataType()== AliCaloTrackReader::kMC) GetReader()->GetVertex(vertex);
\r
380 //Select the Calorimeter of the photon
\r
381 if(fCalorimeter == "PHOS")
\r
383 else if (fCalorimeter == "EMCAL")
\r
384 pl = GetAODEMCAL();
\r
385 //Fill AODCaloClusters and AODParticle with PHOS aods
\r
386 TLorentzVector mom ;
\r
387 for(Int_t icalo = 0; icalo < pl->GetEntriesFast(); icalo++){
\r
388 AliAODCaloCluster * calo = (AliAODCaloCluster*) (pl->At(icalo));
\r
390 //Cluster selection, not charged, with photon id and in fidutial cut
\r
391 //Get Momentum vector,
\r
392 calo->GetMomentum(mom,vertex);//Assume that come from vertex in straight line
\r
393 //If too small or big pt, skip it
\r
394 if(mom.Pt() < GetMinPt() || mom.Pt() > GetMaxPt() ) continue ;
\r
395 //Check acceptance selection
\r
396 if(IsFidutialCutOn()){
\r
397 Bool_t in = GetFidutialCut()->IsInFidutialCut(mom,fCalorimeter) ;
\r
398 if(! in ) continue ;
\r
401 //Create AOD for analysis
\r
402 AliAODPWG4Particle aodph = AliAODPWG4Particle(mom);
\r
403 aodph.SetLabel(calo->GetLabel(0));
\r
404 aodph.SetDetector(fCalorimeter);
\r
405 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
407 //Check Distance to Bad channel, set bit.
\r
408 Double_t distBad=calo->GetDistToBadChannel() ; //Distance to bad channel
\r
409 if(distBad < 0.) distBad=9999. ; //workout strange convension dist = -1. ;
\r
410 if(distBad < fMinDist) //In bad channel (PHOS cristal size 2.2x2.2 cm)
\r
413 if(GetDebug() > 1) printf("AliAnaPhoton::FillAOD: Bad channel cut passed %4.2f\n",distBad);
\r
415 if(distBad > fMinDist3) aodph.SetDistToBad(2) ;
\r
416 else if(distBad > fMinDist2) aodph.SetDistToBad(1) ;
\r
417 else aodph.SetDistToBad(0) ;
\r
420 //PID selection or bit setting
\r
421 if(GetReader()->GetDataType() == AliCaloTrackReader::kMC){
\r
422 //Get most probable PID, check PID weights (in MC this option is mandatory)
\r
423 aodph.SetPdg(GetCaloPID()->GetPdg(fCalorimeter,calo->PID(),mom.E()));//PID with weights
\r
424 if(GetDebug() > 1) printf("AliAnaPhoton::FillAOD: PDG of identified particle %d\n",aodph.GetPdg());
\r
425 //If primary is not photon, skip it.
\r
426 if(aodph.GetPdg() != AliCaloPID::kPhoton) continue ;
\r
428 else if(IsCaloPIDOn()){
\r
429 //Skip matched clusters with tracks
\r
430 if(calo->GetNTracksMatched() > 0) continue ;
\r
432 //Get most probable PID, 2 options check PID weights
\r
433 //or redo PID, recommended option for EMCal.
\r
434 if(!IsCaloPIDRecalculationOn())
\r
435 aodph.SetPdg(GetCaloPID()->GetPdg(fCalorimeter,calo->PID(),mom.E()));//PID with weights
\r
437 aodph.SetPdg(GetCaloPID()->GetPdg(fCalorimeter,mom,calo));//PID recalculated
\r
439 if(GetDebug() > 1) printf("AliAnaPhoton::FillAOD: PDG of identified particle %d\n",aodph.GetPdg());
\r
441 //If cluster does not pass pid, not photon, skip it.
\r
442 if(aodph.GetPdg() != AliCaloPID::kPhoton) continue ;
\r
446 //Set PID bits for later selection (AliAnaPi0 for example)
\r
447 //GetPDG already called in SetPIDBits.
\r
448 GetCaloPID()->SetPIDBits(fCalorimeter,calo,&aodph);
\r
449 if(GetDebug() > 1) printf("AliAnaPhoton::FillAOD: PID Bits set \n");
\r
452 if(GetDebug() > 1) printf("AliAnaPhoton::FillAOD: Photon selection cuts passed: pT %3.2f, pdg %d\n",aodph.Pt(), aodph.GetPdg());
\r
454 //Play with the MC stack if available
\r
455 //Check origin of the candidates
\r
457 aodph.SetTag(GetCaloPID()->CheckOrigin(calo->GetLabel(0),GetMCStack()));
\r
458 if(GetDebug() > 0) printf("AliAnaPi0EbE::FillAOD: Origin of candidate %d\n",aodph.GetTag());
\r
459 }//Work with stack also
\r
461 //Add AOD with photon object to aod branch
\r
462 AddAODParticle(aodph);
\r
466 if(GetDebug() > 1) printf("AliAnaPhoton::FillAOD: End fill AODs \n");
\r
469 //__________________________________________________________________
\r
470 void AliAnaPi0EbE::MakeAnalysisFillHistograms()
\r
472 //Do analysis and fill histograms
\r
474 if(!GetOutputAODBranch())
\r
475 AliFatal(Form("MakeInvMassInCalo: No output pi0 in AOD branch with name < %s > \n",GetOutputAODName().Data()));
\r
477 //Loop on stored AOD pi0
\r
478 Int_t naod = GetOutputAODBranch()->GetEntriesFast();
\r
479 if(GetDebug() > 0) printf("pi0 aod branch entries %d\n", naod);
\r
481 for(Int_t iaod = 0; iaod < naod ; iaod++){
\r
483 AliAODPWG4Particle* pi0 = (AliAODPWG4Particle*) (GetOutputAODBranch()->At(iaod));
\r
484 Int_t pdg = pi0->GetPdg();
\r
486 if(pdg != AliCaloPID::kPi0) continue;
\r
488 //Fill pi0 histograms
\r
489 Float_t pt = pi0->Pt();
\r
490 Float_t phi = pi0->Phi();
\r
491 Float_t eta = pi0->Eta();
\r
494 fhPtPi0 ->Fill(pt);
\r
495 fhPhiPi0 ->Fill(pt,phi);
\r
496 fhEtaPi0 ->Fill(pt,eta);
\r
499 if(pi0->GetTag()== AliCaloPID::kMCPi0Decay){
\r
500 fhPtMCPi0 ->Fill(pt);
\r
501 fhPhiMCPi0 ->Fill(pt,phi);
\r
502 fhEtaMCPi0 ->Fill(pt,eta);
\r
505 fhPtMCNoPi0 ->Fill(pt);
\r
506 fhPhiMCNoPi0 ->Fill(pt,phi);
\r
507 fhEtaMCNoPi0 ->Fill(pt,eta);
\r
509 }//Histograms with MC
\r
516 //____________________________________________________________________________
\r
517 void AliAnaPi0EbE::InitParameters()
\r
520 //Initialize the parameters of the analysis.
\r
521 SetOutputAODClassName("AliAODPWG4Particle");
\r
522 SetOutputAODName("pi0s");
\r
523 fInputAODGammaConvName = "gammaconv" ;
\r
524 fAnaType = kIMCalo ;
\r
525 fCalorimeter = "PHOS" ;
\r
532 //__________________________________________________________________
\r
533 void AliAnaPi0EbE::Print(const Option_t * opt) const
\r
536 //Print some relevant parameters set for the analysis
\r
540 printf("**** Print %s %s ****\n", GetName(), GetTitle() ) ;
\r
541 AliAnaPartCorrBaseClass::Print("");
\r
542 printf("Analysis Type = %d \n", fAnaType) ;
\r
543 if(fAnaType == kSSCalo){
\r
544 printf("Calorimeter = %s\n", fCalorimeter.Data()) ;
\r
545 printf("Min Distance to Bad Channel = %2.1f\n",fMinDist);
\r
546 printf("Min Distance to Bad Channel 2 = %2.1f\n",fMinDist2);
\r
547 printf("Min Distance to Bad Channel 3 = %2.1f\n",fMinDist3);
\r