]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWG4/PartCorrDep/AliAnaPi0EbE.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Possilibity to change the binning and range of time histogras
[u/mrichter/AliRoot.git] / PWG4 / PartCorrDep / AliAnaPi0EbE.cxx
CommitLineData
477d6cee 1/**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
6 * *
7 * Permission to use, copy, modify and distribute this software and its *
21a4b1c0 8 * documentation strictly for non-commercial purposes is hereby granted *
477d6cee 9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
15/* $Id: AliAnaPi0EbE.cxx 28688 2008-09-11 15:04:07Z gconesab $ */
16
17//_________________________________________________________________________
18// Class for the analysis of high pT pi0 event by event
19// Pi0 identified by one of the following:
20// -Invariant mass of 2 cluster in calorimeter
21// -Shower shape analysis in calorimeter
21a4b1c0 22// -Invariant mass of one cluster in calorimeter and one photon reconstructed in CTS
477d6cee 23//
24// -- Author: Gustavo Conesa (LNF-INFN) & Raphaelle Ichou (SUBATECH)
25//////////////////////////////////////////////////////////////////////////////
26
27
28// --- ROOT system ---
29#include <TList.h>
30#include <TClonesArray.h>
0c1383b5 31#include <TObjString.h>
57b97dc6 32#include <TH3F.h>
477d6cee 33//#include "Riostream.h"
34
35// --- Analysis system ---
36#include "AliAnaPi0EbE.h"
37#include "AliCaloTrackReader.h"
38#include "AliIsolationCut.h"
39#include "AliNeutralMesonSelection.h"
40#include "AliCaloPID.h"
41#include "AliMCAnalysisUtils.h"
477d6cee 42#include "AliStack.h"
ff45398a 43#include "AliFiducialCut.h"
477d6cee 44#include "TParticle.h"
0ae57829 45#include "AliVCluster.h"
477d6cee 46#include "AliAODEvent.h"
591cc579 47#include "AliAODMCParticle.h"
477d6cee 48
49ClassImp(AliAnaPi0EbE)
50
691bdd02 51//____________________________________________________________________________
477d6cee 52AliAnaPi0EbE::AliAnaPi0EbE() :
53AliAnaPartCorrBaseClass(), fAnaType(kIMCalo),fCalorimeter(""),
54fMinDist(0.),fMinDist2(0.),fMinDist3(0.),
55fInputAODGammaConv(0x0),fInputAODGammaConvName(""),
57b97dc6 56fHistoSSBins(100), fHistoSSMax(5), fHistoSSMin(0),
ddc0a8a5 57fHistoDiffTimeBins(800), fHistoDiffTimeMax(400), fHistoDiffTimeMin(-400),
c4a7d28a 58fhPtPi0(0), fhEPi0(0), fhEEtaPhiPi0(0),
59fhEDispPi0(0), fhEDispBkg(0),
60fhELambda0Pi0(0), fhELambda0Bkg(0),
61fhELambda1Pi0(0), fhELambda1Bkg(0),
62fhELambdaPi0EtaCen(0),fhELambdaPi0EtaMid(0),fhELambdaPi0EtaBor(0),
63fhClusterPairDiffTimeE(0),fhClusterPairDiffTimeAsy(0),
64//MC histos
65fhELambdaPhNoConv(0),fhELambdaConvPhotons(0),fhELambdaElectrons(0),
66fhELambdaPi0NoPh(0),fhELambdaOtherParts(0),
477d6cee 67fhPtMCNoPi0(0),fhPhiMCNoPi0(0),fhEtaMCNoPi0(0),
68fhPtMCPi0(0),fhPhiMCPi0(0),fhEtaMCPi0(0)
69{
70 //default ctor
71
72 //Initialize parameters
73 InitParameters();
74
75}
477d6cee 76
477d6cee 77//____________________________________________________________________________
78AliAnaPi0EbE::~AliAnaPi0EbE()
79{
80 //dtor
81 if(fInputAODGammaConv){
82 fInputAODGammaConv->Clear() ;
83 delete fInputAODGammaConv ;
84 }
85}
86
0c1383b5 87//________________________________________________________________________
88TObjString * AliAnaPi0EbE::GetAnalysisCuts()
89{
90 //Save parameters used for analysis
91 TString parList ; //this will be list of parameters used for this analysis.
5ae09196 92 const Int_t buffersize = 255;
93 char onePar[buffersize] ;
0c1383b5 94
5ae09196 95 snprintf(onePar,buffersize,"--- AliAnaPi0EbE ---\n") ;
0c1383b5 96 parList+=onePar ;
5ae09196 97 snprintf(onePar,buffersize,"fAnaType=%d (Pi0 selection type) \n",fAnaType) ;
0c1383b5 98 parList+=onePar ;
99
100 if(fAnaType == kSSCalo){
5ae09196 101 snprintf(onePar,buffersize,"Calorimeter: %s\n",fCalorimeter.Data()) ;
0c1383b5 102 parList+=onePar ;
5ae09196 103 snprintf(onePar,buffersize,"fMinDist =%2.2f (Minimal distance to bad channel to accept cluster) \n",fMinDist) ;
0c1383b5 104 parList+=onePar ;
5ae09196 105 snprintf(onePar,buffersize,"fMinDist2=%2.2f (Cuts on Minimal distance to study acceptance evaluation) \n",fMinDist2) ;
0c1383b5 106 parList+=onePar ;
5ae09196 107 snprintf(onePar,buffersize,"fMinDist3=%2.2f (One more cut on distance used for acceptance-efficiency study) \n",fMinDist3) ;
0c1383b5 108 parList+=onePar ;
109 }
110
111 //Get parameters set in base class.
112 parList += GetBaseParametersList() ;
113
114 //Get parameters set in PID class.
115 if(fAnaType == kSSCalo) parList += GetCaloPID()->GetPIDParametersList() ;
116
117 return new TObjString(parList) ;
118}
119
477d6cee 120//________________________________________________________________________
121TList * AliAnaPi0EbE::GetCreateOutputObjects()
122{
123 // Create histograms to be saved in output file and
124 // store them in outputContainer
125 TList * outputContainer = new TList() ;
126 outputContainer->SetName("Pi0EbEHistos") ;
127
57b97dc6 128 Int_t nptbins = GetHistoPtBins(); Float_t ptmax = GetHistoPtMax(); Float_t ptmin = GetHistoPtMin();
129 Int_t nphibins = GetHistoPhiBins(); Float_t phimax = GetHistoPhiMax(); Float_t phimin = GetHistoPhiMin();
130 Int_t netabins = GetHistoEtaBins(); Float_t etamax = GetHistoEtaMax(); Float_t etamin = GetHistoEtaMin();
131 Int_t ssbins = GetHistoShowerShapeBins(); Float_t ssmax = GetHistoShowerShapeMax(); Float_t ssmin = GetHistoShowerShapeMin();
c4a7d28a 132
477d6cee 133 fhPtPi0 = new TH1F("hPtPi0","Number of identified #pi^{0} decay",nptbins,ptmin,ptmax);
134 fhPtPi0->SetYTitle("N");
135 fhPtPi0->SetXTitle("p_{T #pi^{0}}(GeV/c)");
136 outputContainer->Add(fhPtPi0) ;
c4a7d28a 137
138 fhEPi0 = new TH1F("hEPi0","Number of identified #pi^{0} decay",nptbins,ptmin,ptmax);
139 fhEPi0->SetYTitle("N");
140 fhEPi0->SetXTitle("E #pi^{0}(GeV)");
141 outputContainer->Add(fhEPi0) ;
477d6cee 142
c4a7d28a 143 fhEEtaPhiPi0 = new TH3F
144 ("hEEtaPhiPi0","Selected #pi^{0} pairs: E vs #eta vs #phi",nptbins,ptmin,ptmax,netabins,etamin,etamax, nphibins,phimin,phimax);
145 fhEEtaPhiPi0->SetZTitle("#phi");
146 fhEEtaPhiPi0->SetYTitle("#eta");
147 fhEEtaPhiPi0->SetXTitle("E (GeV)");
148 outputContainer->Add(fhEEtaPhiPi0) ;
57b97dc6 149
c4a7d28a 150 ////////
57b97dc6 151
152 if(fAnaType == kIMCalo){
c4a7d28a 153
154 fhEDispPi0 = new TH2F
155 ("hEDispPi0","Selected #pi^{0} pairs: E vs dispersion",nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
156 fhEDispPi0->SetYTitle("dispersion");
157 fhEDispPi0->SetXTitle("E (GeV)");
158 outputContainer->Add(fhEDispPi0) ;
159
160 fhEDispBkg = new TH2F
161 ("hEDispBkg","Rejected #pi^{0} pairs: E vs dispersion",nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
162 fhEDispBkg->SetYTitle("dispersion");
163 fhEDispBkg->SetXTitle("E (GeV)");
164 outputContainer->Add(fhEDispBkg) ;
165
166 fhELambda0Pi0 = new TH2F
167 ("hELambda0Pi0","Selected #pi^{0} pairs: E vs #lambda_{0}",nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
168 fhELambda0Pi0->SetYTitle("#lambda_{0}");
169 fhELambda0Pi0->SetXTitle("E (GeV)");
170 outputContainer->Add(fhELambda0Pi0) ;
171
172 fhELambda0Bkg = new TH2F
173 ("hELambda0Bkg","Rejected #pi^{0} pairs: E vs #lambda_{0}",nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
174 fhELambda0Bkg->SetYTitle("#lambda_{0}");
175 fhELambda0Bkg->SetXTitle("E (GeV)");
176 outputContainer->Add(fhELambda0Bkg) ;
177
178 fhELambda1Pi0 = new TH2F
179 ("hELambda1Pi0","Selected #pi^{0} pairs: E vs #lambda_{1}",nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
180 fhELambda1Pi0->SetYTitle("#lambda_{1}");
181 fhELambda1Pi0->SetXTitle("E (GeV)");
182 outputContainer->Add(fhELambda1Pi0) ;
183
184 fhELambda1Bkg = new TH2F
185 ("hELambda1Bkg","Rejected #pi^{0} pairs: E vs #lambda_{1}",nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
186 fhELambda1Bkg->SetYTitle("#lambda_{1}");
187 fhELambda1Bkg->SetXTitle("E (GeV)");
188 outputContainer->Add(fhELambda1Bkg) ;
189
190 fhELambdaPi0EtaCen = new TH2F
191 ("hELambdaPi0EtaCen","Selected #pi^{0} pairs: E vs #lambda_{0}, |#eta < 0.2|",nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
192 fhELambdaPi0EtaCen->SetYTitle("#lambda_{0}");
193 fhELambdaPi0EtaCen->SetXTitle("E (GeV)");
194 outputContainer->Add(fhELambdaPi0EtaCen) ;
195
196 fhELambdaPi0EtaMid = new TH2F
197 ("hELambdaPi0EtaMid","Selected #pi^{0} pairs: E vs #lambda_{0}, |0.2< #eta < 0.5|",nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
198 fhELambdaPi0EtaMid->SetYTitle("#lambda_{0}");
199 fhELambdaPi0EtaMid->SetXTitle("E (GeV)");
200 outputContainer->Add(fhELambdaPi0EtaMid) ;
201
202 fhELambdaPi0EtaBor = new TH2F
203 ("hELambdaPi0EtaBor","Selected #pi^{0} pairs: E vs #lambda_{0}, |#eta > 0.5|",nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
204 fhELambdaPi0EtaBor->SetYTitle("#lambda_{0}");
205 fhELambdaPi0EtaBor->SetXTitle("E (GeV)");
206 outputContainer->Add(fhELambdaPi0EtaBor) ;
207
ddc0a8a5 208 fhClusterPairDiffTimeE = new TH2F("hClusterPairDiffTimeE","cluster pair time difference vs E",nptbins,ptmin,ptmax, fHistoDiffTimeBins,fHistoDiffTimeMin,fHistoDiffTimeMax);
c4a7d28a 209 fhClusterPairDiffTimeE->SetXTitle("E_{pair} (GeV)");
210 fhClusterPairDiffTimeE->SetYTitle("#Delta t (ns)");
211 outputContainer->Add(fhClusterPairDiffTimeE);
212
ddc0a8a5 213 fhClusterPairDiffTimeAsy = new TH2F("hClusterPairDiffTime","cluster pair time difference vs pair asymmetry",100,0,1, fHistoDiffTimeBins,fHistoDiffTimeMin,fHistoDiffTimeMax);
c4a7d28a 214 fhClusterPairDiffTimeAsy->SetXTitle("Asymmetry");
215 fhClusterPairDiffTimeAsy->SetYTitle("#Delta t (ns)");
216 outputContainer->Add(fhClusterPairDiffTimeAsy);
57b97dc6 217
218 }// Invariant mass analysis in calorimeters only
477d6cee 219
220 if(IsDataMC()) {
221 if((GetReader()->GetDataType() == AliCaloTrackReader::kMC && fAnaType!=kSSCalo) ||
222 GetReader()->GetDataType() != AliCaloTrackReader::kMC){
223
224 fhPtMCPi0 = new TH1F("hPtMCPi0","Identified pi0 from pi0",nptbins,ptmin,ptmax);
225 fhPtMCPi0->SetYTitle("N");
226 fhPtMCPi0->SetXTitle("p_{T #pi^{0}}(GeV/c)");
227 outputContainer->Add(fhPtMCPi0) ;
228
229 fhPhiMCPi0 = new TH2F
5ae09196 230 ("hPhiMCPi0","Identified pi0 from pi0",nptbins,ptmin,ptmax,nphibins,phimin,phimax);
477d6cee 231 fhPhiMCPi0->SetYTitle("#phi");
232 fhPhiMCPi0->SetXTitle("p_{T #pi^{0}} (GeV/c)");
233 outputContainer->Add(fhPhiMCPi0) ;
234
235 fhEtaMCPi0 = new TH2F
5ae09196 236 ("hEtaMCPi0","Identified pi0 from pi0",nptbins,ptmin,ptmax,netabins,etamin,etamax);
477d6cee 237 fhEtaMCPi0->SetYTitle("#eta");
238 fhEtaMCPi0->SetXTitle("p_{T #pi^{0}} (GeV/c)");
239 outputContainer->Add(fhEtaMCPi0) ;
240
241 fhPtMCNoPi0 = new TH1F("hPtMCNoPi0","Identified pi0 not from pi0",nptbins,ptmin,ptmax);
242 fhPtMCNoPi0->SetYTitle("N");
243 fhPtMCNoPi0->SetXTitle("p_{T #pi^{0}}(GeV/c)");
244 outputContainer->Add(fhPtMCNoPi0) ;
245
246 fhPhiMCNoPi0 = new TH2F
5ae09196 247 ("hPhiMCNoPi0","Identified pi0 not from pi0",nptbins,ptmin,ptmax,nphibins,phimin,phimax);
477d6cee 248 fhPhiMCNoPi0->SetYTitle("#phi");
249 fhPhiMCNoPi0->SetXTitle("p_{T #pi^{0}} (GeV/c)");
250 outputContainer->Add(fhPhiMCNoPi0) ;
251
252 fhEtaMCNoPi0 = new TH2F
5ae09196 253 ("hEtaMCNoPi0","Identified pi0 not from pi0",nptbins,ptmin,ptmax,netabins,etamin,etamax);
477d6cee 254 fhEtaMCNoPi0->SetYTitle("#eta");
255 fhEtaMCNoPi0->SetXTitle("p_{T #pi^{0}} (GeV/c)");
256 outputContainer->Add(fhEtaMCNoPi0) ;
257
c4a7d28a 258 if(fAnaType == kIMCalo){
259
260 fhELambdaPhNoConv = new TH2F
261 ("hELambdaPhNoConv","Selected #pi^{0} pairs (Really photons and not conversion): E vs #lambda_{0} vs #lambda_{1}",nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
262 fhELambdaPhNoConv->SetZTitle("#lambda_{1}");
263 fhELambdaPhNoConv->SetYTitle("#lambda_{0}");
264 fhELambdaPhNoConv->SetXTitle("E (GeV)");
265 outputContainer->Add(fhELambdaPhNoConv) ;
266
267
268 fhELambdaConvPhotons = new TH2F
269 ("hELambdaConvPhotons","Selected #pi^{0} pairs (Converted photons): E vs #lambda_{0} vs #lambda_{1}",nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
270 fhELambdaConvPhotons->SetZTitle("#lambda_{1}");
271 fhELambdaConvPhotons->SetYTitle("#lambda_{0}");
272 fhELambdaConvPhotons->SetXTitle("E (GeV)");
273 outputContainer->Add(fhELambdaConvPhotons) ;
274
275 fhELambdaElectrons = new TH2F
276 ("hELambdaElectrons","Selected #pi^{0} pairs(Electrons ): E vs #lambda_{0} vs #lambda_{1}",nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
277 fhELambdaElectrons->SetZTitle("#lambda_{1}");
278 fhELambdaElectrons->SetYTitle("#lambda_{0}");
279 fhELambdaElectrons->SetXTitle("E (GeV)");
280 outputContainer->Add(fhELambdaElectrons) ;
281
282 fhELambdaPi0NoPh = new TH2F
283 ("hELambdaPi0NoPh","Selected #pi^{0} pairs(Pi0s instead of photons): E vs #lambda_{0} vs #lambda_{1}",nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
284 fhELambdaPi0NoPh->SetZTitle("#lambda_{1}");
285 fhELambdaPi0NoPh->SetYTitle("#lambda_{0}");
286 fhELambdaPi0NoPh->SetXTitle("E (GeV)");
287 outputContainer->Add(fhELambdaPi0NoPh) ;
288
289 fhELambdaOtherParts = new TH2F
290 ("hELambdaOtherParts","Selected #pi^{0} pairs (Other parts): E vs #lambda_{0} vs #lambda_{1}",nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
291 fhELambdaOtherParts->SetZTitle("#lambda_{1}");
292 fhELambdaOtherParts->SetYTitle("#lambda_{0}");
293 fhELambdaOtherParts->SetXTitle("E (GeV)");
294 outputContainer->Add(fhELambdaOtherParts) ;
295
296 }//kIMCalo
297
477d6cee 298 }
299 }//Histos with MC
300
301
302 //Keep neutral meson selection histograms if requiered
303 //Setting done in AliNeutralMesonSelection
304
305 if(fAnaType!=kSSCalo && GetNeutralMesonSelection()){
306
307 TList * nmsHistos = GetNeutralMesonSelection()->GetCreateOutputObjects() ;
308 if(GetNeutralMesonSelection()->AreNeutralMesonSelectionHistosKept())
309 for(Int_t i = 0; i < nmsHistos->GetEntries(); i++) outputContainer->Add(nmsHistos->At(i)) ;
5ae09196 310 delete nmsHistos;
a14fd526 311
477d6cee 312 }
313
477d6cee 314 return outputContainer ;
315
316}
317
318//__________________________________________________________________
319void AliAnaPi0EbE::MakeAnalysisFillAOD()
320{
321 //Do analysis and fill aods
322
323 switch(fAnaType)
324 {
325 case kIMCalo:
326 MakeInvMassInCalorimeter();
327 break;
328
329 case kSSCalo:
330 MakeShowerShapeIdentification();
331 break;
332
333 case kIMCaloTracks:
334 MakeInvMassInCalorimeterAndCTS();
335 break;
336
337 }
338}
339
340//__________________________________________________________________
341void AliAnaPi0EbE::MakeInvMassInCalorimeter()
342{
57b97dc6 343 //Do analysis and fill aods
344 //Search for the photon decay in calorimeters
345 //Read photon list from AOD, produced in class AliAnaPhoton
346 //Check if 2 photons have the mass of the pi0.
477d6cee 347
348 TLorentzVector mom1;
349 TLorentzVector mom2;
350 TLorentzVector mom ;
591cc579 351 Int_t tag1 = 0;
352 Int_t tag2 = 0;
353 Int_t tag = 0;
477d6cee 354
355 if(!GetInputAODBranch()){
a3aebfff 356 printf("AliAnaPi0EbE::MakeInvMassInCalorimeter() - No input calo photons in AOD with name branch < %s >, STOP \n",GetInputAODName().Data());
477d6cee 357 abort();
358 }
f8006433 359
c4a7d28a 360 for(Int_t iphoton = 0; iphoton < GetInputAODBranch()->GetEntriesFast()-1; iphoton++){
477d6cee 361 AliAODPWG4Particle * photon1 = (AliAODPWG4Particle*) (GetInputAODBranch()->At(iphoton));
c8fe2783 362
c4a7d28a 363 //Vertex cut in case of mixed events
c8fe2783 364 Int_t evtIndex1 = 0 ;
365 if(GetMixedEvent())
366 evtIndex1 = GetMixedEvent()->EventIndexForCaloCluster(photon1->GetCaloLabel(0)) ;
5025c139 367 if(TMath::Abs(GetVertex(evtIndex1)[2]) > GetZvertexCut()) continue ; //vertex cut
477d6cee 368 mom1 = *(photon1->Momentum());
369
c4a7d28a 370 //Get shower shape information of clusters
371 TObjArray *clusters = 0;
372 if (fCalorimeter="EMCAL") clusters = GetEMCALClusters();
373 else if(fCalorimeter="PHOS" ) clusters = GetPHOSClusters() ;
9ab9e937 374
c4a7d28a 375 Bool_t bFound1 = kFALSE;
9ab9e937 376 Int_t caloLabel1 = photon1->GetCaloLabel(0);
377 Bool_t iclus1 = -1;
378 AliVCluster *cluster1 = 0;
379 if(clusters) {
380 //Get original cluster, to recover some information
381 for(Int_t iclus = 0; iclus < clusters->GetEntriesFast(); iclus++){
382 AliVCluster *cluster= dynamic_cast<AliVCluster*> (clusters->At(iclus));
383 if(cluster){
384 if (cluster->GetID()==caloLabel1) {
385 bFound1 = kTRUE ;
386 cluster1 = cluster;
387 iclus1 = iclus;
388 }
389 }
390 if(bFound1) break;
391 }// calorimeter clusters loop
392 }
c4a7d28a 393
394 for(Int_t jphoton = iphoton+1; jphoton < GetInputAODBranch()->GetEntriesFast(); jphoton++){
477d6cee 395
a3aebfff 396 AliAODPWG4Particle * photon2 = (AliAODPWG4Particle*) (GetInputAODBranch()->At(jphoton));
c8fe2783 397 Int_t evtIndex2 = 0 ;
398 if(GetMixedEvent())
399 evtIndex2 = GetMixedEvent()->EventIndexForCaloCluster(photon2->GetCaloLabel(0)) ;
400 if(GetMixedEvent() && (evtIndex1 == evtIndex2))
401 continue ;
5025c139 402 if(TMath::Abs(GetVertex(evtIndex2)[2]) > GetZvertexCut()) continue ; //vertex cut
477d6cee 403 mom2 = *(photon2->Momentum());
c4a7d28a 404
9ab9e937 405 //Photon1
c4a7d28a 406 Float_t e1 = photon1->E();
407 Float_t eta1 = photon1->Eta();
c4a7d28a 408 Float_t tof1 = -1;
409 Float_t disp1 = -1;
410 Float_t l01 = -1;
411 Float_t l11 = -1;
c4a7d28a 412
9ab9e937 413 //Photon2
c4a7d28a 414 Float_t e2 = photon2->E();
415 Float_t eta2 = photon2->Eta();
416 Float_t disp2 = -1;
417 Float_t tof2 = -1;
418 Float_t l02 = -1;
419 Float_t l12 = -1;
c4a7d28a 420
9ab9e937 421 Bool_t bFound2 = kFALSE;
422 Int_t caloLabel2 = photon2->GetCaloLabel(0);
423 AliVCluster *cluster2 = 0;
424 if(clusters){
425 for(Int_t iclus = iclus1+1; iclus < clusters->GetEntriesFast(); iclus++){
426 AliVCluster *cluster= dynamic_cast<AliVCluster*> (clusters->At(iclus));
427 if(cluster){
428 if(cluster->GetID()==caloLabel2) {
429 bFound2 = kTRUE ;
430 cluster2 = cluster;
431 }
432 }
433 if(bFound2) break;
434 }// calorimeter clusters loop
435
436 //Photon/Cluster 1
437 if(cluster1 && bFound1){
438 disp1 = cluster1->GetDispersion();
439 l11 = cluster1->GetM20();
440 l01 = cluster1->GetM02();
441 tof1 = cluster1->GetTOF()*1e9;
442 }
443// else printf("cluster1 not available: calo label %d / %d, cluster ID %d\n",
444// photon2->GetCaloLabel(0),(GetReader()->GetInputEvent())->GetNumberOfCaloClusters()-1,cluster1->GetID());
445
446 //Photon/Cluster 2
447 if(cluster2 && bFound2){
448 disp2 = cluster2->GetDispersion();
449 l12 = cluster2->GetM20();
450 l02 = cluster2->GetM02();
451 tof2 = cluster2->GetTOF()*1e9;
452 }
453// else printf("cluster2 not available: calo label %d / %d, cluster ID %d\n",
454// photon2->GetCaloLabel(0),(GetReader()->GetInputEvent())->GetNumberOfCaloClusters()-1,cluster2->GetID());
455
456 //Select clusters with good time window difference
457 Double_t t12diff = tof1-tof2;
458 Float_t asymmetry = TMath::Abs(e1-e2)/(e1+e2);
459 fhClusterPairDiffTimeE ->Fill(e1+e2, t12diff);
460 fhClusterPairDiffTimeAsy->Fill(asymmetry,t12diff);
461 if(TMath::Abs(t12diff) > GetPairTimeCut()) continue;
462 }
c4a7d28a 463
57b97dc6 464 //Select good pair (good phi, pt cuts, aperture and invariant mass)
477d6cee 465 if(GetNeutralMesonSelection()->SelectPair(mom1, mom2))
c8fe2783 466 {
467 if(GetDebug()>1)
468 printf("AliAnaPi0EbE::MakeInvMassInCalorimeter() - Selected gamma pair: pt %f, phi %f, eta%f \n",(mom1+mom2).Pt(), (mom1+mom2).Phi()*180./3.1416, (mom1+mom2).Eta());
469
57b97dc6 470 //Play with the MC stack if available
c8fe2783 471 if(IsDataMC()){
57b97dc6 472 //Check origin of the candidates
c8fe2783 473 Int_t label1 = photon1->GetLabel();
474 Int_t label2 = photon2->GetLabel();
0db79843 475 if(label1>=0)tag1 = GetMCAnalysisUtils()->CheckOrigin(label1, GetReader(), photon1->GetInputFileIndex());
476 if(label2>=0)tag2 = GetMCAnalysisUtils()->CheckOrigin(label2, GetReader(), photon2->GetInputFileIndex());
c8fe2783 477
478 if(GetDebug() > 0) printf("AliAnaPi0EbE::MakeInvMassInCalorimeter() - Origin of: photon1 %d; photon2 %d \n",tag1, tag2);
0db79843 479 if(GetMCAnalysisUtils()->CheckTagBit(tag1,AliMCAnalysisUtils::kMCPi0Decay) &&
480 GetMCAnalysisUtils()->CheckTagBit(tag2,AliMCAnalysisUtils::kMCPi0Decay)){
c8fe2783 481
57b97dc6 482 //Check if pi0 mother is the same
c8fe2783 483 if(GetReader()->ReadStack()){
0db79843 484 if(label1>=0){
485 TParticle * mother1 = GetMCStack()->Particle(label1);//photon in kine tree
486 label1 = mother1->GetFirstMother();
487 //mother1 = GetMCStack()->Particle(label1);//pi0
488 }
489 if(label2>=0){
490 TParticle * mother2 = GetMCStack()->Particle(label2);//photon in kine tree
491 label2 = mother2->GetFirstMother();
492 //mother2 = GetMCStack()->Particle(label2);//pi0
493 }
c8fe2783 494 }
f8006433 495 else if(GetReader()->ReadAODMCParticles()){//&& (input > -1)){
0db79843 496 if(label1>=0){
497 AliAODMCParticle * mother1 = (AliAODMCParticle *) (GetReader()->GetAODMCParticles(photon1->GetInputFileIndex()))->At(label1);//photon in kine tree
498 label1 = mother1->GetMother();
499 //mother1 = GetMCStack()->Particle(label1);//pi0
500 }
501 if(label2>=0){
502 AliAODMCParticle * mother2 = (AliAODMCParticle *) (GetReader()->GetAODMCParticles(photon2->GetInputFileIndex()))->At(label2);//photon in kine tree
503 label2 = mother2->GetMother();
504 //mother2 = GetMCStack()->Particle(label2);//pi0
505 }
c8fe2783 506 }
507
57b97dc6 508 //printf("mother1 %d, mother2 %d\n",label1,label2);
0db79843 509 if(label1 == label2 && label1>=0)
c8fe2783 510 GetMCAnalysisUtils()->SetTagBit(tag,AliMCAnalysisUtils::kMCPi0);
511 }
512 }//Work with stack also
513
57b97dc6 514 //Fill some histograms about shower shape
9ab9e937 515 if(clusters && GetReader()->GetDataType()!=AliCaloTrackReader::kMC){
edffc439 516 //Photon1
c4a7d28a 517
518 //printf("Signal Cl1: e %f, pt %f, disp %f, l1 %f, l0 %f, eta %f, phi %f \n",
519 // e,pt,disp,l1,l0,photon2->Eta(),photon2->Phi());
520
521 fhEDispPi0 ->Fill(e1, disp1);
522 fhELambda0Pi0->Fill(e1, l01 );
523 fhELambda1Pi0->Fill(e1, l11 );
524
525 if(TMath::Abs(eta1)<0.2){
526 fhELambdaPi0EtaCen->Fill(e1, l01);
527 }
528 else if (TMath::Abs(eta1)>0.2 && TMath::Abs(eta1)<0.5){
529 fhELambdaPi0EtaMid->Fill(e1, l01);
530 }
531 else {
532 fhELambdaPi0EtaBor->Fill(e1, l01);
533 }
534
edffc439 535 //Photon2
c4a7d28a 536 //printf("Signal Cl2: e %f, pt %f, disp %f, l1 %f, l0 %f, eta %f, phi %f \n",e
537 // ,pt,disp,l1,l0,photon2->Eta(),photon2->Phi());
538
539 fhEDispPi0 ->Fill(e2, disp2);
540 fhELambda0Pi0->Fill(e2, l02 );
541 fhELambda1Pi0->Fill(e2, l12 );
542
543 if(TMath::Abs(eta2)<0.2){
544 fhELambdaPi0EtaCen->Fill(e2, l02);
545 }
546 else if (TMath::Abs(eta2)>0.2 && TMath::Abs(eta2)<0.5){
547 fhELambdaPi0EtaMid->Fill(e2, l02);
548 }
549 else {
550 fhELambdaPi0EtaBor->Fill(e2, l02);
551 }
552
553 if(IsDataMC()) {
554 //Photon1
555 if( GetMCAnalysisUtils()->CheckTagBit(tag1,AliMCAnalysisUtils::kMCPhoton) && !GetMCAnalysisUtils()->CheckTagBit(tag1,AliMCAnalysisUtils::kMCConversion) ){
556 fhELambdaPhNoConv->Fill(e1, l01);
557 }//photon no conversion
558 else if ( GetMCAnalysisUtils()->CheckTagBit(tag1,AliMCAnalysisUtils::kMCElectron)){
559 fhELambdaElectrons->Fill(e1, l01);
560 }//electron
561 else if ( GetMCAnalysisUtils()->CheckTagBit(tag1,AliMCAnalysisUtils::kMCConversion) && GetMCAnalysisUtils()->CheckTagBit(tag1,AliMCAnalysisUtils::kMCConversion) ){
562 fhELambdaConvPhotons->Fill(e1, l01);
563 }//convesion photon
564
565 else if ( GetMCAnalysisUtils()->CheckTagBit(tag1,AliMCAnalysisUtils::kMCPi0) ){
566 fhELambdaPi0NoPh->Fill(e1, l01);
567 }//pi0
568 else {
569 fhELambdaOtherParts->Fill(e1, l01);
570 }//other particules
571
572 //Photon 2
573 if( GetMCAnalysisUtils()->CheckTagBit(tag2,AliMCAnalysisUtils::kMCPhoton) && !GetMCAnalysisUtils()->CheckTagBit(tag2,AliMCAnalysisUtils::kMCConversion) ){
574 fhELambdaPhNoConv->Fill(e2, l02);
575 }//photon no conversion
576 else if ( GetMCAnalysisUtils()->CheckTagBit(tag2,AliMCAnalysisUtils::kMCElectron)){
577 fhELambdaElectrons->Fill(e2, l02);
578 }//electron
579 else if ( GetMCAnalysisUtils()->CheckTagBit(tag2,AliMCAnalysisUtils::kMCConversion) && GetMCAnalysisUtils()->CheckTagBit(tag2,AliMCAnalysisUtils::kMCConversion) ){
580 fhELambdaConvPhotons->Fill(e2, l02);
581 }//convesion photon
582
583 else if ( GetMCAnalysisUtils()->CheckTagBit(tag2,AliMCAnalysisUtils::kMCPi0) ){
584 fhELambdaPi0NoPh->Fill(e2, l02);
585 }//pi0
586 else {
587 fhELambdaOtherParts->Fill(e2, l02);
588 }//other particules
589 }//is datamc
590
591 }//MC histograms
592
57b97dc6 593 //Create AOD for analysis
c8fe2783 594 mom = mom1+mom2;
595 AliAODPWG4Particle pi0 = AliAODPWG4Particle(mom);
57b97dc6 596 //pi0.SetLabel(calo->GetLabel());
21a4b1c0 597 pi0.SetIdentifiedParticleType(AliCaloPID::kPi0);
c8fe2783 598 pi0.SetDetector(photon1->GetDetector());
599 pi0.SetTag(tag);
57b97dc6 600 //Set the indeces of the original caloclusters
c8fe2783 601 pi0.SetCaloLabel(photon1->GetCaloLabel(0), photon2->GetCaloLabel(0));
f8006433 602 //pi0.SetInputFileIndex(input);
c8fe2783 603 AddAODParticle(pi0);
604 }//pi0
57b97dc6 605 else{
606 Float_t phi = (mom1+mom2).Phi();
607 if(phi < 0) phi+=TMath::TwoPi();
57b97dc6 608
609 //Fill some histograms about shower shape
edffc439 610 if(GetReader()->GetDataType()!=AliCaloTrackReader::kMC){
611 //Photon1
c4a7d28a 612
613 //printf("Bkg Cl1: e %f, pt %f, disp %f, l1 %f, l0 %f, eta %f, phi %f \n",
614 // e,pt,disp,l1,l0,photon2->Eta(),photon2->Phi());
615
616 fhEDispBkg ->Fill(e1, disp1);
617 fhELambda0Bkg->Fill(e1, l01 );
618 fhELambda1Bkg->Fill(e1, l11 );
619
edffc439 620 //Photon2
c4a7d28a 621
622 //printf("Bkg Cl2: e %f, pt %f, disp %f, l1 %f, l0 %f, eta %f, phi %f \n",
623 //e,pt,disp,l1,l0,photon2->Eta(),photon2->Phi());
624
625 fhEDispBkg ->Fill(e2, disp2);
626 fhELambda0Bkg->Fill(e2, l02 );
627 fhELambda1Bkg->Fill(e2, l12 );
628
edffc439 629 }
57b97dc6 630
631 }//bkg pair
632
477d6cee 633 }//2n photon loop
634
635 }//1st photon loop
636
a3aebfff 637 if(GetDebug() > 1) printf("AliAnaPi0EbE::MakeInvMassInCalorimeter() - End fill AODs \n");
477d6cee 638
639}
640
641//__________________________________________________________________
642void AliAnaPi0EbE::MakeInvMassInCalorimeterAndCTS()
643{
644 //Do analysis and fill aods
645 //Search for the photon decay in calorimeters
646 //Read photon list from AOD, produced in class AliAnaPhoton and AliGammaConversion
647 //Check if 2 photons have the mass of the pi0.
648
649 TLorentzVector mom1;
650 TLorentzVector mom2;
651 TLorentzVector mom ;
591cc579 652 Int_t tag1 = 0;
653 Int_t tag2 = 0;
654 Int_t tag = 0;
5025c139 655 Int_t evtIndex = 0;
477d6cee 656 if(!GetInputAODBranch()){
a3aebfff 657 printf("AliAnaPi0EbE::MakeInvMassInCalorimeterAndCTS() - No input calo photons in AOD branch with name < %s > , STOP\n",GetInputAODName().Data());
477d6cee 658 abort();
659 }
57b97dc6 660
477d6cee 661 for(Int_t iphoton = 0; iphoton < GetInputAODBranch()->GetEntriesFast(); iphoton++){
662 AliAODPWG4Particle * photon1 = (AliAODPWG4Particle*) (GetInputAODBranch()->At(iphoton));
663 mom1 = *(photon1->Momentum());
664
665 //Play with the MC stack if available
666 fInputAODGammaConv = (TClonesArray *) GetReader()->GetOutputEvent()->FindListObject(fInputAODGammaConvName);
667 if(!fInputAODGammaConv) {
a3aebfff 668 printf("AliAnaPi0EbE::MakeInvMassInCalorimeterAndCTS() - No input gamma conversions in AOD branch with name < %s >, STOP\n",fInputAODGammaConvName.Data());
477d6cee 669 abort();
670 }
c4a7d28a 671 for(Int_t jphoton = 0; jphoton < fInputAODGammaConv->GetEntriesFast(); jphoton++){
a3aebfff 672 AliAODPWG4Particle * photon2 = (AliAODPWG4Particle*) (fInputAODGammaConv->At(jphoton));
5025c139 673 if(GetMixedEvent())
674 evtIndex = GetMixedEvent()->EventIndexForCaloCluster(photon2->GetCaloLabel(0)) ;
675 if(TMath::Abs(GetVertex(evtIndex)[2]) > GetZvertexCut()) continue ; //vertex cut
676
477d6cee 677 mom2 = *(photon2->Momentum());
57b97dc6 678
f8006433 679 //Int_t input = -1; //if -1 photons come from different files, not a pi0
680 //if(photon1->GetInputFileIndex() == photon2->GetInputFileIndex()) input = photon1->GetInputFileIndex();
57b97dc6 681
477d6cee 682 //Select good pair (good phi, pt cuts, aperture and invariant mass)
683 if(GetNeutralMesonSelection()->SelectPair(mom1, mom2)){
57b97dc6 684 if(GetDebug() > 1) printf("AliAnaPi0EbE::MakeInvMassInCalorimeterAndCTS() - Selected gamma pair: pt %f, phi %f, eta%f\n",(mom1+mom2).Pt(), (mom1+mom2).Phi()*180./3.1416, (mom1+mom2).Eta());
685
686 if(IsDataMC()){
687 Int_t label1 = photon1->GetLabel();
688 Int_t label2 = photon2->GetLabel();
0db79843 689 if(label1>=0)tag1 = GetMCAnalysisUtils()->CheckOrigin(label1, GetReader(), photon1->GetInputFileIndex());
690 if(label2>=0)tag2 = GetMCAnalysisUtils()->CheckOrigin(label2, GetReader(), photon2->GetInputFileIndex());
57b97dc6 691 if(GetDebug() > 0) printf("AliAnaPi0EbE::MakeInvMassInCalorimeterAndCTS() - Origin of: photon1 %d; photon2 %d \n",tag1, tag2);
0db79843 692 if(GetMCAnalysisUtils()->CheckTagBit(tag1,AliMCAnalysisUtils::kMCPi0Decay) &&
693 GetMCAnalysisUtils()->CheckTagBit(tag2,AliMCAnalysisUtils::kMCPi0Decay)){
57b97dc6 694 //Check if pi0 mother is the same
695
696 if(GetReader()->ReadStack()){
0db79843 697 if(label1>=0){
698 TParticle * mother1 = GetMCStack()->Particle(label1);//photon in kine tree
699 label1 = mother1->GetFirstMother();
700 //mother1 = GetMCStack()->Particle(label1);//pi0
701 }
702 if(label2>=0){
703 TParticle * mother2 = GetMCStack()->Particle(label2);//photon in kine tree
704 label2 = mother2->GetFirstMother();
705 //mother2 = GetMCStack()->Particle(label2);//pi0
706 }
57b97dc6 707 }
0db79843 708 else if(GetReader()->ReadAODMCParticles()&& label1>=0 && label2>=0){ //&& (input > -1)){
709 if(label1>=0){
710 AliAODMCParticle * mother1 = (AliAODMCParticle *) (GetReader()->GetAODMCParticles(photon1->GetInputFileIndex()))->At(label1);//photon in kine tree
711 label1 = mother1->GetMother();
712 //mother1 = GetMCStack()->Particle(label1);//pi0
713 }
714 if(label2>=0){
715 AliAODMCParticle * mother2 = (AliAODMCParticle *) (GetReader()->GetAODMCParticles(photon2->GetInputFileIndex()))->At(label2);//photon in kine tree
716 label2 = mother2->GetMother();
717 //mother2 = GetMCStack()->Particle(label2);//pi0
718 }
57b97dc6 719 }
720
721 //printf("mother1 %d, mother2 %d\n",label1,label2);
0db79843 722 if(label1 == label2 && label1>=0)
57b97dc6 723 GetMCAnalysisUtils()->SetTagBit(tag,AliMCAnalysisUtils::kMCPi0);
724 }
725 }//Work with stack also
726
727 //Create AOD for analysis
728 mom = mom1+mom2;
729 AliAODPWG4Particle pi0 = AliAODPWG4Particle(mom);
730 //pi0.SetLabel(calo->GetLabel());
21a4b1c0 731 pi0.SetIdentifiedParticleType(AliCaloPID::kPi0);
57b97dc6 732 pi0.SetDetector(photon1->GetDetector());
733 pi0.SetTag(tag);
734 //Set the indeces of the original tracks or caloclusters
735 pi0.SetCaloLabel(photon1->GetCaloLabel(0), -1);
736 pi0.SetTrackLabel(photon2->GetTrackLabel(0), photon2->GetTrackLabel(1));
f8006433 737 //pi0.SetInputFileIndex(input);
57b97dc6 738 AddAODParticle(pi0);
477d6cee 739 }//pi0
740 }//2n photon loop
741
742 }//1st photon loop
743
a3aebfff 744 if(GetDebug() > 1) printf("AliAnaPi0EbE::MakeInvMassInCalorimeterAndCTS() - End fill AODs \n");
477d6cee 745
746}
747
748
749//__________________________________________________________________
750void AliAnaPi0EbE::MakeShowerShapeIdentification()
751{
752 //Search for pi0 in fCalorimeter with shower shape analysis
753
4a745797 754 TObjArray * pl = 0x0;
5ae09196 755 //Select the Calorimeter of the photon
756 if(fCalorimeter == "PHOS")
be518ab0 757 pl = GetPHOSClusters();
5ae09196 758 else if (fCalorimeter == "EMCAL")
be518ab0 759 pl = GetEMCALClusters();
57b97dc6 760
5ae09196 761 if(!pl) {
762 Info("MakeShowerShapeIdentification","TObjArray with %s clusters is NULL!\n",fCalorimeter.Data());
763 return;
764 }
477d6cee 765
766 //Get vertex for photon momentum calculation
f37fa8d2 767 //Double_t vertex2[] = {0,0,0} ; //vertex from second aod input
f8006433 768 //if(GetReader()->GetDataType() != AliCaloTrackReader::kMC)
769 //{
770 //if(GetReader()->GetSecondInputAODTree()) GetReader()->GetSecondInputAODVertex(vertex2);
771 //}
233e0df8 772
57b97dc6 773
477d6cee 774 TLorentzVector mom ;
775 for(Int_t icalo = 0; icalo < pl->GetEntriesFast(); icalo++){
0ae57829 776 AliVCluster * calo = (AliVCluster*) (pl->At(icalo));
477d6cee 777
f8006433 778 Int_t evtIndex = 0 ;
779 if (GetMixedEvent()) {
780 evtIndex=GetMixedEvent()->EventIndexForCaloCluster(calo->GetID()) ;
781 }
5025c139 782 if(TMath::Abs(GetVertex(evtIndex)[2]) > GetZvertexCut()) continue ; //vertex cut
783
ff45398a 784 //Cluster selection, not charged, with pi0 id and in fiducial cut
233e0df8 785
57b97dc6 786 //Input from second AOD?
f8006433 787 //Int_t input = 0;
be518ab0 788 // if (fCalorimeter == "EMCAL" && GetReader()->GetEMCALClustersNormalInputEntries() <= icalo) input = 1 ;
789 // else if(fCalorimeter == "PHOS" && GetReader()->GetPHOSClustersNormalInputEntries() <= icalo) input = 1;
233e0df8 790
57b97dc6 791 //Get Momentum vector,
f8006433 792 //if (input == 0)
793 if(GetReader()->GetDataType() != AliCaloTrackReader::kMC){
794 calo->GetMomentum(mom,GetVertex(evtIndex)) ;}//Assume that come from vertex in straight line
795 else{
796 Double_t vertex[]={0,0,0};
797 calo->GetMomentum(mom,vertex) ;
798 }
57b97dc6 799 //else if(input == 1) calo->GetMomentum(mom,vertex2);//Assume that come from vertex in straight line
233e0df8 800
57b97dc6 801 //If too small or big pt, skip it
477d6cee 802 if(mom.Pt() < GetMinPt() || mom.Pt() > GetMaxPt() ) continue ;
803 //Check acceptance selection
ff45398a 804 if(IsFiducialCutOn()){
805 Bool_t in = GetFiducialCut()->IsInFiducialCut(mom,fCalorimeter) ;
477d6cee 806 if(! in ) continue ;
807 }
808
809 //Create AOD for analysis
810 AliAODPWG4Particle aodpi0 = AliAODPWG4Particle(mom);
c8fe2783 811 aodpi0.SetLabel(calo->GetLabel());
477d6cee 812 //Set the indeces of the original caloclusters
813 aodpi0.SetCaloLabel(calo->GetID(),-1);
814 aodpi0.SetDetector(fCalorimeter);
815 if(GetDebug() > 1)
ff45398a 816 printf("AliAnaPi0EbE::MakeShowerShapeIdentification() - FillAOD: Min pt cut and fiducial cut passed: pt %3.2f, phi %2.2f, eta %1.2f\n",aodpi0.Pt(),aodpi0.Phi(),aodpi0.Eta());
477d6cee 817
818 //Check Distance to Bad channel, set bit.
c8fe2783 819 Double_t distBad=calo->GetDistanceToBadChannel() ; //Distance to bad channel
477d6cee 820 if(distBad < 0.) distBad=9999. ; //workout strange convension dist = -1. ;
821 if(distBad < fMinDist) //In bad channel (PHOS cristal size 2.2x2.2 cm)
822 continue ;
823
a3aebfff 824 if(GetDebug() > 1) printf("AliAnaPi0EbE::MakeShowerShapeIdentification() - FillAOD: Bad channel cut passed %4.2f\n",distBad);
477d6cee 825
826 if(distBad > fMinDist3) aodpi0.SetDistToBad(2) ;
827 else if(distBad > fMinDist2) aodpi0.SetDistToBad(1) ;
828 else aodpi0.SetDistToBad(0) ;
829
830 //Check PID
831 //PID selection or bit setting
832 if(GetReader()->GetDataType() == AliCaloTrackReader::kMC){
833 //Get most probable PID, check PID weights (in MC this option is mandatory)
21a4b1c0 834 aodpi0.SetIdentifiedParticleType(GetCaloPID()->GetIdentifiedParticleType(fCalorimeter,calo->GetPID(),mom.E()));//PID with weights
477d6cee 835 if(GetDebug() > 1)
21a4b1c0 836 printf("AliAnaPi0EbE::MakeShowerShapeIdentification() - FillAOD: PDG of identified particle %d\n",aodpi0.GetIdentifiedParticleType());
477d6cee 837 //If primary is not pi0, skip it.
21a4b1c0 838 if(aodpi0.GetIdentifiedParticleType() != AliCaloPID::kPi0) continue ;
477d6cee 839 }
840 else if(IsCaloPIDOn()){
841 //Skip matched clusters with tracks
f37fa8d2 842 if(IsTrackMatched(calo)) continue ;
477d6cee 843
844 //Get most probable PID, 2 options check PID weights
845 //or redo PID, recommended option for EMCal.
846 if(!IsCaloPIDRecalculationOn())
21a4b1c0 847 aodpi0.SetIdentifiedParticleType(GetCaloPID()->GetIdentifiedParticleType(fCalorimeter,calo->GetPID(),mom.E()));//PID with weights
477d6cee 848 else
21a4b1c0 849 aodpi0.SetIdentifiedParticleType(GetCaloPID()->GetIdentifiedParticleType(fCalorimeter,mom,calo));//PID recalculated
477d6cee 850
21a4b1c0 851 if(GetDebug() > 1) printf("AliAnaPi0EbE::MakeShowerShapeIdentification() - PDG of identified particle %d\n",aodpi0.GetIdentifiedParticleType());
477d6cee 852
853 //If cluster does not pass pid, not pi0, skip it.
21a4b1c0 854 if(aodpi0.GetIdentifiedParticleType() != AliCaloPID::kPi0) continue ;
477d6cee 855
856 }
857 else{
858 //Set PID bits for later selection (AliAnaPi0 for example)
859 //GetPDG already called in SetPIDBits.
f2ccb5b8 860 GetCaloPID()->SetPIDBits(fCalorimeter,calo,&aodpi0, GetCaloUtils());
a3aebfff 861 if(GetDebug() > 1) printf("AliAnaPi0EbE::MakeShowerShapeIdentification() - PID Bits set \n");
477d6cee 862 }
863
21a4b1c0 864 if(GetDebug() > 1) printf("AliAnaPi0EbE::MakeShowerShapeIdentification() - Pi0 selection cuts passed: pT %3.2f, pdg %d\n",aodpi0.Pt(), aodpi0.GetIdentifiedParticleType());
477d6cee 865
866 //Play with the MC stack if available
867 //Check origin of the candidates
868 if(IsDataMC()){
869 if((GetReader()->GetDataType() == AliCaloTrackReader::kMC && fAnaType!=kSSCalo) ||
57b97dc6 870 GetReader()->GetDataType() != AliCaloTrackReader::kMC){
f8006433 871 //aodpi0.SetInputFileIndex(input);
57b97dc6 872 Int_t tag =0;
873 tag = GetMCAnalysisUtils()->CheckOrigin(calo->GetLabel(),GetReader(), aodpi0.GetInputFileIndex());
874 //GetMCAnalysisUtils()->CheckMultipleOrigin(calo->GetLabels(),calo->GetNLabels(), GetReader(), aodpi0.GetInputFileIndex(), tag);
875 aodpi0.SetTag(tag);
876 if(GetDebug() > 0) printf("AliAnaPi0EbE::MakeShowerShapeIdentification() - Origin of candidate %d\n",aodpi0.GetTag());
477d6cee 877 }
878 }//Work with stack also
879
880 //Add AOD with pi0 object to aod branch
881 AddAODParticle(aodpi0);
882
883 }//loop
884
a3aebfff 885 if(GetDebug() > 1) printf("AliAnaPi0EbE::MakeShowerShapeIdentification() - End fill AODs \n");
477d6cee 886
887}
888//__________________________________________________________________
889void AliAnaPi0EbE::MakeAnalysisFillHistograms()
691bdd02 890{
477d6cee 891 //Do analysis and fill histograms
691bdd02 892
477d6cee 893 if(!GetOutputAODBranch()){
a3aebfff 894 printf("AliAnaPi0EbE::MakeAnalysisFillHistograms() - No output pi0 in AOD branch with name < %s >,STOP \n",GetOutputAODName().Data());
477d6cee 895 abort();
896 }
897 //Loop on stored AOD pi0
898 Int_t naod = GetOutputAODBranch()->GetEntriesFast();
a3aebfff 899 if(GetDebug() > 0) printf("AliAnaPi0EbE::MakeAnalysisFillHistograms() - aod branch entries %d\n", naod);
477d6cee 900
901 for(Int_t iaod = 0; iaod < naod ; iaod++){
902
903 AliAODPWG4Particle* pi0 = (AliAODPWG4Particle*) (GetOutputAODBranch()->At(iaod));
21a4b1c0 904 Int_t pdg = pi0->GetIdentifiedParticleType();
9415d854 905
906 if(IsCaloPIDOn() && pdg != AliCaloPID::kPi0) continue;
477d6cee 907
908 //Fill pi0 histograms
c4a7d28a 909 Float_t ener = pi0->E();
910 Float_t pt = pi0->Pt();
911 Float_t phi = pi0->Phi();
57b97dc6 912 if(phi < 0) phi+=TMath::TwoPi();
477d6cee 913 Float_t eta = pi0->Eta();
914
c4a7d28a 915 fhPtPi0 ->Fill(pt);
916 fhEPi0 ->Fill(ener);
917 fhEEtaPhiPi0 ->Fill(ener,eta,phi);
477d6cee 918
919 if(IsDataMC()){
920 if((GetReader()->GetDataType() == AliCaloTrackReader::kMC && fAnaType!=kSSCalo) ||
57b97dc6 921 GetReader()->GetDataType() != AliCaloTrackReader::kMC){
922 if(GetMCAnalysisUtils()->CheckTagBit(pi0->GetTag(), AliMCAnalysisUtils::kMCPi0)){
923 fhPtMCPi0 ->Fill(pt);
924 fhPhiMCPi0 ->Fill(pt,phi);
925 fhEtaMCPi0 ->Fill(pt,eta);
926 }
927 else{
928 fhPtMCNoPi0 ->Fill(pt);
929 fhPhiMCNoPi0 ->Fill(pt,phi);
930 fhEtaMCNoPi0 ->Fill(pt,eta);
931 }
477d6cee 932 }
933 }//Histograms with MC
934
935 }// aod loop
936
937}
938
939
940//____________________________________________________________________________
941void AliAnaPi0EbE::Init()
942{
691bdd02 943 //Init
944 //Do some checks
1e86c71e 945 if(fCalorimeter == "PHOS" && !GetReader()->IsPHOSSwitchedOn() && NewOutputAOD()){
591cc579 946 printf("AliAnaPi0EbE::Init() - !!STOP: You want to use PHOS in analysis but it is not read!! \n!!Check the configuration file!!\n");
691bdd02 947 abort();
948 }
1e86c71e 949 else if(fCalorimeter == "EMCAL" && !GetReader()->IsEMCALSwitchedOn() && NewOutputAOD()){
591cc579 950 printf("AliAnaPi0EbE::Init() - !!STOP: You want to use EMCAL in analysis but it is not read!! \n!!Check the configuration file!!\n");
691bdd02 951 abort();
952 }
953
954}
955
477d6cee 956//____________________________________________________________________________
957void AliAnaPi0EbE::InitParameters()
958{
1e86c71e 959 //Initialize the parameters of the analysis.
a3aebfff 960 AddToHistogramsName("AnaPi0EbE_");
961
477d6cee 962 fInputAODGammaConvName = "gammaconv" ;
963 fAnaType = kIMCalo ;
964 fCalorimeter = "EMCAL" ;
965 fMinDist = 2.;
966 fMinDist2 = 4.;
967 fMinDist3 = 5.;
968
969}
970
971//__________________________________________________________________
972void AliAnaPi0EbE::Print(const Option_t * opt) const
973{
974 //Print some relevant parameters set for the analysis
975 if(! opt)
976 return;
977
978 printf("**** Print %s %s ****\n", GetName(), GetTitle() ) ;
979 AliAnaPartCorrBaseClass::Print("");
980 printf("Analysis Type = %d \n", fAnaType) ;
981 if(fAnaType == kSSCalo){
982 printf("Calorimeter = %s\n", fCalorimeter.Data()) ;
983 printf("Min Distance to Bad Channel = %2.1f\n",fMinDist);
984 printf("Min Distance to Bad Channel 2 = %2.1f\n",fMinDist2);
985 printf("Min Distance to Bad Channel 3 = %2.1f\n",fMinDist3);
986 }
987 printf(" \n") ;
988
989}