]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWG4/PartCorrDep/AliAnaPi0.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Fixing error on gcc 4.5.1
[u/mrichter/AliRoot.git] / PWG4 / PartCorrDep / AliAnaPi0.cxx
CommitLineData
1c5acb87 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 *
8 * documentation strictly for non-commercial purposes is hereby granted *
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: $ */
16
17//_________________________________________________________________________
18// Class to collect two-photon invariant mass distributions for
19// extractin raw pi0 yield.
20//
21//-- Author: Dmitri Peressounko (RRC "KI")
22//-- Adapted to PartCorr frame by Lamia Benhabib (SUBATECH)
23//-- and Gustavo Conesa (INFN-Frascati)
24//_________________________________________________________________________
25
26
27// --- ROOT system ---
28#include "TH3.h"
50f39b97 29#include "TH2D.h"
1c5acb87 30//#include "Riostream.h"
6639984f 31#include "TCanvas.h"
32#include "TPad.h"
33#include "TROOT.h"
477d6cee 34#include "TClonesArray.h"
0c1383b5 35#include "TObjString.h"
1c5acb87 36
37//---- AliRoot system ----
38#include "AliAnaPi0.h"
39#include "AliCaloTrackReader.h"
40#include "AliCaloPID.h"
6639984f 41#include "AliStack.h"
ff45398a 42#include "AliFiducialCut.h"
477d6cee 43#include "TParticle.h"
477d6cee 44#include "AliVEvent.h"
6921fa00 45#include "AliESDCaloCluster.h"
46#include "AliESDEvent.h"
47#include "AliAODEvent.h"
50f39b97 48#include "AliNeutralMesonSelection.h"
c8fe2783 49#include "AliMixedEvent.h"
50
591cc579 51
1c5acb87 52ClassImp(AliAnaPi0)
53
54//________________________________________________________________________________________________________________________________________________
55AliAnaPi0::AliAnaPi0() : AliAnaPartCorrBaseClass(),
7e7694bb 56fDoOwnMix(kFALSE),fNCentrBin(0),fNZvertBin(0),fNrpBin(0),
1c5acb87 57fNPID(0),fNmaxMixEv(0), fZvtxCut(0.),fCalorimeter(""),
50f39b97 58fNModules(12), fUseAngleCut(kFALSE), fEventsList(0x0), //fhEtalon(0x0),
5a2dbc3c 59fhReMod(0x0), fhRe1(0x0),fhMi1(0x0),fhRe2(0x0),fhMi2(0x0),fhRe3(0x0),fhMi3(0x0),fhEvents(0x0),
50f39b97 60fhRealOpeningAngle(0x0),fhRealCosOpeningAngle(0x0),
61fhPrimPt(0x0), fhPrimAccPt(0x0), fhPrimY(0x0), fhPrimAccY(0x0), fhPrimPhi(0x0), fhPrimAccPhi(0x0),
62fhPrimOpeningAngle(0x0),fhPrimCosOpeningAngle(0x0)
1c5acb87 63{
64//Default Ctor
65 InitParameters();
66
67}
78219bac 68/*
1c5acb87 69//________________________________________________________________________________________________________________________________________________
70AliAnaPi0::AliAnaPi0(const AliAnaPi0 & ex) : AliAnaPartCorrBaseClass(ex),
71fNCentrBin(ex.fNCentrBin),fNZvertBin(ex.fNZvertBin),fNrpBin(ex.fNrpBin),
c0649ccf 72fNPID(ex.fNPID),fNmaxMixEv(ex.fNmaxMixEv),fZvtxCut(ex.fZvtxCut), fCalorimeter(ex.fCalorimeter),
614701c6 73fNModules(ex.fNModules), fUseAngleCut(ex.fUseAngleCut), fEventsList(0x0), //fhEtalon(ex.fhEtalon),
5a2dbc3c 74fhReMod(ex.fhReMod), fhRe1(ex.fhRe1),fhMi1(ex.fhMi1),fhRe2(ex.fhRe2),fhMi2(ex.fhMi2),
75fhRe3(ex.fhRe3),fhMi3(ex.fhMi3),fhEvents(ex.fhEvents),
50f39b97 76fhRealOpeningAngle(ex.fhRealOpeningAngle),fhRealCosOpeningAngle(ex.fhRealCosOpeningAngle),
1c5acb87 77fhPrimPt(ex.fhPrimPt), fhPrimAccPt(ex.fhPrimAccPt), fhPrimY(ex.fhPrimY),
50f39b97 78fhPrimAccY(ex.fhPrimAccY), fhPrimPhi(ex.fhPrimPhi), fhPrimAccPhi(ex.fhPrimAccPhi),
79fhPrimOpeningAngle(ex.fhPrimOpeningAngle),fhPrimCosOpeningAngle(ex.fhPrimCosOpeningAngle)
1c5acb87 80{
477d6cee 81 // cpy ctor
82 //Do not need it
1c5acb87 83}
78219bac 84*/
614701c6 85//
86////________________________________________________________________________________________________________________________________________________
87//AliAnaPi0 & AliAnaPi0::operator = (const AliAnaPi0 & ex)
88//{
89// // assignment operator
90//
91// if(this == &ex)return *this;
92// ((AliAnaPartCorrBaseClass *)this)->operator=(ex);
93//
94// fNCentrBin = ex.fNCentrBin ; fNZvertBin = ex.fNZvertBin ; fNrpBin = ex.fNrpBin ;
95// fNPID = ex.fNPID ; fNmaxMixEv = ex.fNmaxMixEv ; fZvtxCut = ex.fZvtxCut ; fCalorimeter = ex.fCalorimeter ;
96// fNModules = ex.fNModules; fEventsList = ex.fEventsList ; //fhEtalon = ex.fhEtalon ;
97// fhRe1 = ex.fhRe1 ; fhMi1 = ex.fhMi1 ; fhRe2 = ex.fhRe2 ; fhMi2 = ex.fhMi2 ; fhReMod = ex.fhReMod;
98// fhRe3 = ex.fhRe3 ; fhMi3 = ex.fhMi3 ; fhEvents = ex.fhEvents ; fUseAngleCut = ex.fUseAngleCut;
99// fhPrimPt = ex.fhPrimPt ; fhPrimAccPt = ex.fhPrimAccPt ; fhPrimY = ex.fhPrimY ;
100// fhPrimAccY = ex.fhPrimAccY ; fhPrimPhi = ex.fhPrimPhi ; fhPrimAccPhi = ex.fhPrimAccPhi ;
101// fhRealOpeningAngle = ex.fhRealOpeningAngle; fhRealCosOpeningAngle = ex.fhRealCosOpeningAngle;
102// fhPrimOpeningAngle = ex.fhPrimOpeningAngle; fhPrimCosOpeningAngle = ex.fhPrimCosOpeningAngle;
103//
104// return *this;
105//
106//}
1c5acb87 107
108//________________________________________________________________________________________________________________________________________________
109AliAnaPi0::~AliAnaPi0() {
477d6cee 110 // Remove event containers
7e7694bb 111
112 if(fDoOwnMix && fEventsList){
477d6cee 113 for(Int_t ic=0; ic<fNCentrBin; ic++){
114 for(Int_t iz=0; iz<fNZvertBin; iz++){
7e7694bb 115 for(Int_t irp=0; irp<fNrpBin; irp++){
116 fEventsList[ic*fNZvertBin*fNrpBin+iz*fNrpBin+irp]->Delete() ;
117 delete fEventsList[ic*fNZvertBin*fNrpBin+iz*fNrpBin+irp] ;
118 }
477d6cee 119 }
120 }
121 delete[] fEventsList;
122 fEventsList=0 ;
123 }
591cc579 124
1c5acb87 125}
126
127//________________________________________________________________________________________________________________________________________________
128void AliAnaPi0::InitParameters()
129{
130//Init parameters when first called the analysis
131//Set default parameters
a3aebfff 132 SetInputAODName("PWG4Particle");
133
134 AddToHistogramsName("AnaPi0_");
6921fa00 135 fNModules = 12; // set maximum to maximum number of EMCAL modules
477d6cee 136 fNCentrBin = 1;
137 fNZvertBin = 1;
138 fNrpBin = 1;
139 fNPID = 9;
140 fNmaxMixEv = 10;
141 fZvtxCut = 40;
142 fCalorimeter = "PHOS";
50f39b97 143 fUseAngleCut = kFALSE;
144
1c5acb87 145}
146//________________________________________________________________________________________________________________________________________________
5a2dbc3c 147//void AliAnaPi0::Init()
148//{
477d6cee 149 //Init some data members needed in analysis
150
151 //Histograms binning and range
5a2dbc3c 152 // if(!fhEtalon){ // p_T alpha d m_gg
153 // fhEtalon = new TH3D("hEtalon","Histo with binning parameters",50,0.,25.,10,0.,1.,200,0.,1.) ;
154 // fhEtalon->SetXTitle("P_{T} (GeV)") ;
155 // fhEtalon->SetYTitle("#alpha") ;
156 // fhEtalon->SetZTitle("m_{#gamma#gamma} (GeV)") ;
157 // }
477d6cee 158
5a2dbc3c 159//}
1c5acb87 160
0c1383b5 161
162//________________________________________________________________________________________________________________________________________________
163TObjString * AliAnaPi0::GetAnalysisCuts()
164{
165 //Save parameters used for analysis
166 TString parList ; //this will be list of parameters used for this analysis.
167 char onePar[255] ;
168 sprintf(onePar,"--- AliAnaPi0 ---\n") ;
169 parList+=onePar ;
170 sprintf(onePar,"Number of bins in Centrality: %d \n",fNCentrBin) ;
171 parList+=onePar ;
172 sprintf(onePar,"Number of bins in Z vert. pos: %d \n",fNZvertBin) ;
173 parList+=onePar ;
174 sprintf(onePar,"Number of bins in Reac. Plain: %d \n",fNrpBin) ;
175 parList+=onePar ;
176 sprintf(onePar,"Depth of event buffer: %d \n",fNmaxMixEv) ;
177 parList+=onePar ;
178 sprintf(onePar,"Number of different PID used: %d \n",fNPID) ;
179 parList+=onePar ;
180 sprintf(onePar,"Cuts: \n") ;
181 parList+=onePar ;
182 sprintf(onePar,"Z vertex position: -%f < z < %f \n",fZvtxCut,fZvtxCut) ;
183 parList+=onePar ;
184 sprintf(onePar,"Calorimeter: %s \n",fCalorimeter.Data()) ;
185 parList+=onePar ;
186 sprintf(onePar,"Number of modules: %d \n",fNModules) ;
187 parList+=onePar ;
188
189 return new TObjString(parList) ;
190}
191
2e557d1c 192//________________________________________________________________________________________________________________________________________________
193TList * AliAnaPi0::GetCreateOutputObjects()
194{
477d6cee 195 // Create histograms to be saved in output file and
196 // store them in fOutputContainer
197
198 //create event containers
199 fEventsList = new TList*[fNCentrBin*fNZvertBin*fNrpBin] ;
1c5acb87 200
477d6cee 201 for(Int_t ic=0; ic<fNCentrBin; ic++){
202 for(Int_t iz=0; iz<fNZvertBin; iz++){
203 for(Int_t irp=0; irp<fNrpBin; irp++){
7e7694bb 204 fEventsList[ic*fNZvertBin*fNrpBin+iz*fNrpBin+irp] = new TList() ;
477d6cee 205 }
206 }
207 }
7e7694bb 208
477d6cee 209 TList * outputContainer = new TList() ;
210 outputContainer->SetName(GetName());
6921fa00 211
212 fhReMod = new TH3D*[fNModules] ;
50305ae9 213 fhRe1 = new TH3D*[fNCentrBin*fNPID] ;
214 fhRe2 = new TH3D*[fNCentrBin*fNPID] ;
215 fhRe3 = new TH3D*[fNCentrBin*fNPID] ;
216 fhMi1 = new TH3D*[fNCentrBin*fNPID] ;
217 fhMi2 = new TH3D*[fNCentrBin*fNPID] ;
218 fhMi3 = new TH3D*[fNCentrBin*fNPID] ;
477d6cee 219
220 char key[255] ;
221 char title[255] ;
5a2dbc3c 222 Int_t nptbins = GetHistoPtBins();
223 Int_t nphibins = GetHistoPhiBins();
224 Int_t netabins = GetHistoEtaBins();
225 Float_t ptmax = GetHistoPtMax();
226 Float_t phimax = GetHistoPhiMax();
227 Float_t etamax = GetHistoEtaMax();
228 Float_t ptmin = GetHistoPtMin();
229 Float_t phimin = GetHistoPhiMin();
230 Float_t etamin = GetHistoEtaMin();
231
232 Int_t nmassbins = GetHistoMassBins();
233 Int_t nasymbins = GetHistoAsymmetryBins();
234 Float_t massmax = GetHistoMassMax();
235 Float_t asymmax = GetHistoAsymmetryMax();
236 Float_t massmin = GetHistoMassMin();
237 Float_t asymmin = GetHistoAsymmetryMin();
238
477d6cee 239 for(Int_t ic=0; ic<fNCentrBin; ic++){
240 for(Int_t ipid=0; ipid<fNPID; ipid++){
7e7694bb 241
477d6cee 242 //Distance to bad module 1
243 sprintf(key,"hRe_cen%d_pid%d_dist1",ic,ipid) ;
244 sprintf(title,"Real m_{#gamma#gamma} distr. for centrality=%d and PID=%d",ic,ipid) ;
5a2dbc3c 245 //fhEtalon->Clone(key);
246 //fhRe1[ic*fNPID+ipid]=(TH3D*)fhEtalon->Clone(key) ;
247 //fhRe1[ic*fNPID+ipid]->SetName(key) ;
248 //fhRe1[ic*fNPID+ipid]->SetTitle(title) ;
7e7694bb 249 fhRe1[ic*fNPID+ipid] = new TH3D(key,title,nptbins,ptmin,ptmax,nasymbins,asymmin,asymmax,nmassbins,massmin,massmax) ;
477d6cee 250 outputContainer->Add(fhRe1[ic*fNPID+ipid]) ;
7e7694bb 251
477d6cee 252 //Distance to bad module 2
253 sprintf(key,"hRe_cen%d_pid%d_dist2",ic,ipid) ;
254 sprintf(title,"Real m_{#gamma#gamma} distr. for centrality=%d and PID=%d",ic,ipid) ;
5a2dbc3c 255 //fhRe2[ic*fNPID+ipid]=(TH3D*)fhEtalon->Clone(key) ;
256 //fhRe2[ic*fNPID+ipid]->SetName(key) ;
257 //fhRe2[ic*fNPID+ipid]->SetTitle(title) ;
7e7694bb 258 fhRe2[ic*fNPID+ipid] = new TH3D(key,title,nptbins,ptmin,ptmax,nasymbins,asymmin,asymmax,nmassbins,massmin,massmax) ;
477d6cee 259 outputContainer->Add(fhRe2[ic*fNPID+ipid]) ;
260
477d6cee 261 //Distance to bad module 3
262 sprintf(key,"hRe_cen%d_pid%d_dist3",ic,ipid) ;
263 sprintf(title,"Real m_{#gamma#gamma} distr. for centrality=%d and PID=%d",ic,ipid) ;
5a2dbc3c 264 //fhRe3[ic*fNPID+ipid]=(TH3D*)fhEtalon->Clone(key) ;
265 //fhRe3[ic*fNPID+ipid]->SetName(key) ;
266 //fhRe3[ic*fNPID+ipid]->SetTitle(title) ;
7e7694bb 267 fhRe3[ic*fNPID+ipid] = new TH3D(key,title,nptbins,ptmin,ptmax,nasymbins,asymmin,asymmax,nmassbins,massmin,massmax) ;
477d6cee 268 outputContainer->Add(fhRe3[ic*fNPID+ipid]) ;
269
7e7694bb 270 if(fDoOwnMix){
271 //Distance to bad module 1
272 sprintf(key,"hMi_cen%d_pid%d_dist1",ic,ipid) ;
273 sprintf(title,"Mixed m_{#gamma#gamma} distr. for centrality=%d and PID=%d",ic,ipid) ;
274 //fhMi1[ic*fNPID+ipid]=(TH3D*)fhEtalon->Clone(key) ;
275 //fhMi1[ic*fNPID+ipid]->SetName(key) ;
276 //fhMi1[ic*fNPID+ipid]->SetTitle(title) ;
277 fhMi1[ic*fNPID+ipid] = new TH3D(key,title,nptbins,ptmin,ptmax,nasymbins,asymmin,asymmax,nmassbins,massmin,massmax) ;
278 outputContainer->Add(fhMi1[ic*fNPID+ipid]) ;
279
280 //Distance to bad module 2
281 sprintf(key,"hMi_cen%d_pid%d_dist2",ic,ipid) ;
282 sprintf(title,"Mixed m_{#gamma#gamma} distr. for centrality=%d and PID=%d",ic,ipid) ;
283 //fhMi2[ic*fNPID+ipid]=(TH3D*)fhEtalon->Clone(key) ;
284 //fhMi2[ic*fNPID+ipid]->SetName(key) ;
285 //fhMi2[ic*fNPID+ipid]->SetTitle(title) ;
286 fhMi2[ic*fNPID+ipid] = new TH3D(key,title,nptbins,ptmin,ptmax,nasymbins,asymmin,asymmax,nmassbins,massmin,massmax) ;
287 outputContainer->Add(fhMi2[ic*fNPID+ipid]) ;
288
289 //Distance to bad module 3
290 sprintf(key,"hMi_cen%d_pid%d_dist3",ic,ipid) ;
291 sprintf(title,"Mixed m_{#gamma#gamma} distr. for centrality=%d and PID=%d",ic,ipid) ;
292 //fhMi3[ic*fNPID+ipid]=(TH3D*)fhEtalon->Clone(key) ;
293 //fhMi3[ic*fNPID+ipid]->SetName(key) ;
294 //fhMi3[ic*fNPID+ipid]->SetTitle(title) ;
295 fhMi3[ic*fNPID+ipid] = new TH3D(key,title,nptbins,ptmin,ptmax,nasymbins,asymmin,asymmax,nmassbins,massmin,massmax) ;
296 outputContainer->Add(fhMi3[ic*fNPID+ipid]) ;
297 }
1c5acb87 298 }
477d6cee 299 }
300
477d6cee 301 fhEvents=new TH3D("hEvents","Number of events",fNCentrBin,0.,1.*fNCentrBin,
7e7694bb 302 fNZvertBin,0.,1.*fNZvertBin,fNrpBin,0.,1.*fNrpBin) ;
477d6cee 303 outputContainer->Add(fhEvents) ;
50f39b97 304
305 fhRealOpeningAngle = new TH2D
306 ("hRealOpeningAngle","Angle between all #gamma pair vs E_{#pi^{0}}",nptbins,ptmin,ptmax,200,0,0.5);
307 fhRealOpeningAngle->SetYTitle("#theta(rad)");
308 fhRealOpeningAngle->SetXTitle("E_{ #pi^{0}} (GeV)");
309 outputContainer->Add(fhRealOpeningAngle) ;
7e7694bb 310
50f39b97 311 fhRealCosOpeningAngle = new TH2D
312 ("hRealCosOpeningAngle","Cosinus of angle between all #gamma pair vs E_{#pi^{0}}",nptbins,ptmin,ptmax,200,-1,1);
313 fhRealCosOpeningAngle->SetYTitle("cos (#theta) ");
314 fhRealCosOpeningAngle->SetXTitle("E_{ #pi^{0}} (GeV)");
315 outputContainer->Add(fhRealCosOpeningAngle) ;
316
317
477d6cee 318 //Histograms filled only if MC data is requested
0ae57829 319 if(IsDataMC()){
7e7694bb 320 // if(fhEtalon->GetXaxis()->GetXbins() && fhEtalon->GetXaxis()->GetXbins()->GetSize()){ //Variable bin size
321 // fhPrimPt = new TH1D("hPrimPt","Primary pi0 pt",fhEtalon->GetXaxis()->GetNbins(),fhEtalon->GetXaxis()->GetXbins()->GetArray()) ;
322 // fhPrimAccPt = new TH1D("hPrimAccPt","Primary pi0 pt with both photons in acceptance",fhEtalon->GetXaxis()->GetNbins(),
323 // fhEtalon->GetXaxis()->GetXbins()->GetArray()) ;
324 // }
325 // else{
326 // fhPrimPt = new TH1D("hPrimPt","Primary pi0 pt",fhEtalon->GetXaxis()->GetNbins(),fhEtalon->GetXaxis()->GetXmin(),fhEtalon->GetXaxis()->GetXmax()) ;
327 // fhPrimAccPt = new TH1D("hPrimAccPt","Primary pi0 pt with both photons in acceptance",
328 // fhEtalon->GetXaxis()->GetNbins(),fhEtalon->GetXaxis()->GetXmin(),fhEtalon->GetXaxis()->GetXmax()) ;
329 // }
330
331 fhPrimPt = new TH1D("hPrimPt","Primary pi0 pt",nptbins,ptmin,ptmax) ;
5a2dbc3c 332 fhPrimAccPt = new TH1D("hPrimAccPt","Primary pi0 pt with both photons in acceptance",nptbins,ptmin,ptmax) ;
477d6cee 333 outputContainer->Add(fhPrimPt) ;
334 outputContainer->Add(fhPrimAccPt) ;
335
5a2dbc3c 336 fhPrimY = new TH1D("hPrimaryRapidity","Rapidity of primary pi0",netabins,etamin,etamax) ;
477d6cee 337 outputContainer->Add(fhPrimY) ;
338
5a2dbc3c 339 fhPrimAccY = new TH1D("hPrimAccRapidity","Rapidity of primary pi0",netabins,etamin,etamax) ;
477d6cee 340 outputContainer->Add(fhPrimAccY) ;
341
7e7694bb 342 fhPrimPhi = new TH1D("hPrimaryPhi","Azimithal of primary pi0",nphibins,phimin*TMath::RadToDeg(),phimax*TMath::RadToDeg()) ;
477d6cee 343 outputContainer->Add(fhPrimPhi) ;
344
5a2dbc3c 345 fhPrimAccPhi = new TH1D("hPrimAccPhi","Azimithal of primary pi0 with accepted daughters",nphibins,phimin*TMath::RadToDeg(),phimax*TMath::RadToDeg()) ;
477d6cee 346 outputContainer->Add(fhPrimAccPhi) ;
50f39b97 347
348
349 fhPrimOpeningAngle = new TH2D
7e7694bb 350 ("hPrimOpeningAngle","Angle between all primary #gamma pair vs E_{#pi^{0}}",nptbins,ptmin,ptmax,100,0,0.5);
50f39b97 351 fhPrimOpeningAngle->SetYTitle("#theta(rad)");
352 fhPrimOpeningAngle->SetXTitle("E_{ #pi^{0}} (GeV)");
353 outputContainer->Add(fhPrimOpeningAngle) ;
354
355 fhPrimCosOpeningAngle = new TH2D
7e7694bb 356 ("hPrimCosOpeningAngle","Cosinus of angle between all primary #gamma pair vs E_{#pi^{0}}",nptbins,ptmin,ptmax,100,-1,1);
50f39b97 357 fhPrimCosOpeningAngle->SetYTitle("cos (#theta) ");
358 fhPrimCosOpeningAngle->SetXTitle("E_{ #pi^{0}} (GeV)");
359 outputContainer->Add(fhPrimCosOpeningAngle) ;
360
477d6cee 361 }
50f39b97 362
6921fa00 363 for(Int_t imod=0; imod<fNModules; imod++){
50f39b97 364 //Module dependent invariant mass
365 sprintf(key,"hReMod_%d",imod) ;
366 sprintf(title,"Real m_{#gamma#gamma} distr. for Module %d",imod) ;
367 //fhEtalon->Clone(key);
368 //fhReMod[imod]=(TH3D*)fhEtalon->Clone(key) ;
369 //fhReMod[imod]->SetName(key) ;
370 //fhReMod[imod]->SetTitle(title) ;
371 fhReMod[imod] = new TH3D(key,title,nptbins,ptmin,ptmax,nasymbins,asymmin,asymmax,nmassbins,massmin,massmax) ;
372 outputContainer->Add(fhReMod[imod]) ;
6921fa00 373 }
50f39b97 374
477d6cee 375 return outputContainer;
1c5acb87 376}
377
378//_________________________________________________________________________________________________________________________________________________
379void AliAnaPi0::Print(const Option_t * /*opt*/) const
380{
477d6cee 381 //Print some relevant parameters set for the analysis
a3aebfff 382 printf("**** Print %s %s ****\n", GetName(), GetTitle() ) ;
477d6cee 383 AliAnaPartCorrBaseClass::Print(" ");
a3aebfff 384
477d6cee 385 printf("Number of bins in Centrality: %d \n",fNCentrBin) ;
386 printf("Number of bins in Z vert. pos: %d \n",fNZvertBin) ;
387 printf("Number of bins in Reac. Plain: %d \n",fNrpBin) ;
388 printf("Depth of event buffer: %d \n",fNmaxMixEv) ;
389 printf("Number of different PID used: %d \n",fNPID) ;
390 printf("Cuts: \n") ;
391 printf("Z vertex position: -%2.3f < z < %2.3f \n",fZvtxCut,fZvtxCut) ;
50f39b97 392 printf("Number of modules: %d \n",fNModules) ;
393 printf("Select pairs with their angle: %d \n",fUseAngleCut) ;
477d6cee 394 printf("------------------------------------------------------\n") ;
1c5acb87 395}
396
c8fe2783 397
1c5acb87 398//____________________________________________________________________________________________________________________________________________________
6639984f 399void AliAnaPi0::MakeAnalysisFillHistograms()
1c5acb87 400{
477d6cee 401 //Process one event and extract photons from AOD branch
402 // filled with AliAnaPhoton and fill histos with invariant mass
403
404 //Apply some cuts on event: vertex position and centrality range
405 Int_t iRun=(GetReader()->GetInputEvent())->GetRunNumber() ;
406 if(IsBadRun(iRun)) return ;
407
477d6cee 408 Int_t nPhot = GetInputAODBranch()->GetEntriesFast() ;
c8fe2783 409 if(GetDebug() > 1)
410 printf("AliAnaPi0::MakeAnalysisFillHistograms() - Photon entries %d\n", nPhot);
411 if(nPhot < 2 )
412 return ;
6921fa00 413 Int_t module1 = -1;
414 Int_t module2 = -1;
7e7694bb 415 Double_t vert[] = {0.0, 0.0, 0.0} ; //vertex
416 Int_t evtIndex1 = 0 ;
417 Int_t currentEvtIndex = -1 ;
418 Int_t curCentrBin = 0 ;
419 Int_t curRPBin = 0 ;
420 Int_t curZvertBin = 0 ;
421
477d6cee 422 for(Int_t i1=0; i1<nPhot-1; i1++){
423 AliAODPWG4Particle * p1 = (AliAODPWG4Particle*) (GetInputAODBranch()->At(i1)) ;
7e7694bb 424 // get the event index in the mixed buffer where the photon comes from
425 // in case of mixing with analysis frame, not own mixing
c8fe2783 426 evtIndex1 = GetEventIndex(p1, vert) ;
7e7694bb 427 if(vert[2]<-fZvtxCut || vert[2]> fZvtxCut) return ; //Event can not be used (vertex, centrality,... cuts not fulfilled)
c8fe2783 428 if ( evtIndex1 == -1 )
429 return ;
430 if ( evtIndex1 == -2 )
431 continue ;
432 if (evtIndex1 != currentEvtIndex) {
7e7694bb 433 //Get Reaction Plan position and calculate RP bin
434 //does not exist in ESD yet????
c8fe2783 435 curCentrBin = 0 ;
436 curRPBin = 0 ;
437 curZvertBin = (Int_t)(0.5*fNZvertBin*(vert[2]+fZvtxCut)/fZvtxCut) ;
438 fhEvents->Fill(curCentrBin+0.5,curZvertBin+0.5,curRPBin+0.5) ;
439 currentEvtIndex = evtIndex1 ;
440 }
7e7694bb 441
477d6cee 442 TLorentzVector photon1(p1->Px(),p1->Py(),p1->Pz(),p1->E());
59b6bd99 443 //Get Module number
444 module1 = GetModuleNumber(p1);
477d6cee 445 for(Int_t i2=i1+1; i2<nPhot; i2++){
446 AliAODPWG4Particle * p2 = (AliAODPWG4Particle*) (GetInputAODBranch()->At(i2)) ;
c8fe2783 447 Int_t evtIndex2 = GetEventIndex(p2, vert) ;
448 if ( evtIndex2 == -1 )
449 return ;
450 if ( evtIndex2 == -2 )
451 continue ;
452 if (GetMixedEvent() && (evtIndex1 == evtIndex2))
7e7694bb 453 continue ;
477d6cee 454 TLorentzVector photon2(p2->Px(),p2->Py(),p2->Pz(),p2->E());
59b6bd99 455 //Get module number
456 module2 = GetModuleNumber(p2);
477d6cee 457 Double_t m = (photon1 + photon2).M() ;
458 Double_t pt = (photon1 + photon2).Pt();
459 Double_t a = TMath::Abs(p1->E()-p2->E())/(p1->E()+p2->E()) ;
460 if(GetDebug() > 2)
c8fe2783 461 printf("AliAnaPi0::MakeAnalysisFillHistograms() - Current Event: pT: photon1 %2.2f, photon2 %2.2f; Pair: pT %2.2f, mass %2.3f, a %f2.3\n",
462 p1->Pt(), p2->Pt(), pt,m,a);
50f39b97 463 //Check if opening angle is too large or too small compared to what is expected
464 Double_t angle = photon1.Angle(photon2.Vect());
465 //if(fUseAngleCut && !GetNeutralMesonSelection()->IsAngleInWindow((photon1+photon2).E(),angle)) continue;
466 //printf("angle %f\n",angle);
c8fe2783 467 if(fUseAngleCut && angle < 0.1)
468 continue;
50f39b97 469 fhRealOpeningAngle ->Fill(pt,angle);
470 fhRealCosOpeningAngle->Fill(pt,TMath::Cos(angle));
59b6bd99 471 //Fill module dependent histograms
472 //if(module1==module2) printf("mod1 %d\n",module1);
473 if(module1==module2 && module1 >=0 && module1<fNModules)
c8fe2783 474 fhReMod[module1]->Fill(pt,a,m) ;
7e7694bb 475
c8fe2783 476 for(Int_t ipid=0; ipid<fNPID; ipid++){
477 if((p1->IsPIDOK(ipid,AliCaloPID::kPhoton)) && (p2->IsPIDOK(ipid,AliCaloPID::kPhoton))){
478 fhRe1[curCentrBin*fNPID+ipid]->Fill(pt,a,m) ;
479 if(p1->DistToBad()>0 && p2->DistToBad()>0){
480 fhRe2[curCentrBin*fNPID+ipid]->Fill(pt,a,m) ;
481 if(p1->DistToBad()>1 && p2->DistToBad()>1){
482 fhRe3[curCentrBin*fNPID+ipid]->Fill(pt,a,m) ;
483 }
484 }
485 }
486 }
7e7694bb 487 }// second same event particle
488 }// first cluster
489
490 if(fDoOwnMix){
491 //Fill mixed
492 TList * evMixList=fEventsList[curCentrBin*fNZvertBin*fNrpBin+curZvertBin*fNrpBin+curRPBin] ;
493 Int_t nMixed = evMixList->GetSize() ;
494 for(Int_t ii=0; ii<nMixed; ii++){
495 TClonesArray* ev2= (TClonesArray*) (evMixList->At(ii));
496 Int_t nPhot2=ev2->GetEntriesFast() ;
497 Double_t m = -999;
498 if(GetDebug() > 1) printf("AliAnaPi0::MakeAnalysisFillHistograms() - Mixed event %d photon entries %d\n", ii, nPhot);
499
500 for(Int_t i1=0; i1<nPhot; i1++){
501 AliAODPWG4Particle * p1 = (AliAODPWG4Particle*) (GetInputAODBranch()->At(i1)) ;
502 TLorentzVector photon1(p1->Px(),p1->Py(),p1->Pz(),p1->E());
503 for(Int_t i2=0; i2<nPhot2; i2++){
504 AliAODPWG4Particle * p2 = (AliAODPWG4Particle*) (ev2->At(i2)) ;
505
506 TLorentzVector photon2(p2->Px(),p2->Py(),p2->Pz(),p2->E());
507 m = (photon1+photon2).M() ;
508 Double_t pt = (photon1 + photon2).Pt();
509 Double_t a = TMath::Abs(p1->E()-p2->E())/(p1->E()+p2->E()) ;
510
511 //Check if opening angle is too large or too small compared to what is expected
512 Double_t angle = photon1.Angle(photon2.Vect());
513 //if(fUseAngleCut && !GetNeutralMesonSelection()->IsAngleInWindow((photon1+photon2).E(),angle)) continue;
514 if(fUseAngleCut && angle < 0.1) continue;
515
516 if(GetDebug() > 2)
517 printf("AliAnaPi0::MakeAnalysisFillHistograms() - Mixed Event: pT: photon1 %2.2f, photon2 %2.2f; Pair: pT %2.2f, mass %2.3f, a %f2.3\n",
518 p1->Pt(), p2->Pt(), pt,m,a);
519 for(Int_t ipid=0; ipid<fNPID; ipid++){
520 if((p1->IsPIDOK(ipid,AliCaloPID::kPhoton)) && (p2->IsPIDOK(ipid,AliCaloPID::kPhoton))){
521 fhMi1[curCentrBin*fNPID+ipid]->Fill(pt,a,m) ;
522 if(p1->DistToBad()>0 && p2->DistToBad()>0){
523 fhMi2[curCentrBin*fNPID+ipid]->Fill(pt,a,m) ;
524 if(p1->DistToBad()>1 && p2->DistToBad()>1){
525 fhMi3[curCentrBin*fNPID+ipid]->Fill(pt,a,m) ;
526 }
527
528 }
529 }
530 }//loop for histograms
531 }// second cluster loop
532 }//first cluster loop
533 }//loop on mixed events
534
535 TClonesArray *currentEvent = new TClonesArray(*GetInputAODBranch());
536 //Ad d current event to buffer and Remove redundant events
537 if(currentEvent->GetEntriesFast()>0){
538 evMixList->AddFirst(currentEvent) ;
539 currentEvent=0 ; //Now list of particles belongs to buffer and it will be deleted with buffer
540 if(evMixList->GetSize()>=fNmaxMixEv)
541 {
542 TClonesArray * tmp = (TClonesArray*) (evMixList->Last()) ;
543 evMixList->RemoveLast() ;
544 delete tmp ;
545 }
546 }
547 else{ //empty event
548 delete currentEvent ;
549 currentEvent=0 ;
477d6cee 550 }
7e7694bb 551 }// DoOwnMix
c8fe2783 552
7e7694bb 553 //Acceptance
0ae57829 554
de8a210f 555 if(IsDataMC() && GetReader()->ReadStack()){
50f39b97 556 AliStack * stack = GetMCStack();
557 if(stack && (IsDataMC() || (GetReader()->GetDataType() == AliCaloTrackReader::kMC)) ){
558 for(Int_t i=0 ; i<stack->GetNprimary(); i++){
c8fe2783 559 TParticle * prim = stack->Particle(i) ;
560 if(prim->GetPdgCode() == 111){
561 Double_t pi0Pt = prim->Pt() ;
7e7694bb 562 //printf("pi0, pt %2.2f\n",pi0Pt);
c8fe2783 563 if(prim->Energy() == TMath::Abs(prim->Pz())) continue ; //Protection against floating point exception
564 Double_t pi0Y = 0.5*TMath::Log((prim->Energy()-prim->Pz())/(prim->Energy()+prim->Pz())) ;
565 Double_t phi = TMath::RadToDeg()*prim->Phi() ;
566 if(TMath::Abs(pi0Y) < 0.5){
567 fhPrimPt->Fill(pi0Pt) ;
568 }
569 fhPrimY ->Fill(pi0Y) ;
570 fhPrimPhi->Fill(phi) ;
571
7e7694bb 572 //Check if both photons hit Calorimeter
c8fe2783 573 Int_t iphot1=prim->GetFirstDaughter() ;
574 Int_t iphot2=prim->GetLastDaughter() ;
575 if(iphot1>-1 && iphot1<stack->GetNtrack() && iphot2>-1 && iphot2<stack->GetNtrack()){
576 TParticle * phot1 = stack->Particle(iphot1) ;
577 TParticle * phot2 = stack->Particle(iphot2) ;
578 if(phot1 && phot2 && phot1->GetPdgCode()==22 && phot2->GetPdgCode()==22){
7e7694bb 579 //printf("2 photons: photon 1: pt %2.2f, phi %3.2f, eta %1.2f; photon 2: pt %2.2f, phi %3.2f, eta %1.2f\n",
580 // phot1->Pt(), phot1->Phi()*180./3.1415, phot1->Eta(), phot2->Pt(), phot2->Phi()*180./3.1415, phot2->Eta());
c8fe2783 581
582 TLorentzVector lv1, lv2;
583 phot1->Momentum(lv1);
584 phot2->Momentum(lv2);
585
586 Bool_t inacceptance = kFALSE;
587 if(fCalorimeter == "PHOS"){
588 if(GetPHOSGeometry() && GetCaloUtils()->IsPHOSGeoMatrixSet()){
589 Int_t mod ;
590 Double_t x,z ;
591 if(GetPHOSGeometry()->ImpactOnEmc(phot1,mod,z,x) && GetPHOSGeometry()->ImpactOnEmc(phot2,mod,z,x))
592 inacceptance = kTRUE;
593 if(GetDebug() > 2) printf("In %s Real acceptance? %d\n",fCalorimeter.Data(),inacceptance);
594 }
595 else{
596
597 if(GetFiducialCut()->IsInFiducialCut(lv1,fCalorimeter) && GetFiducialCut()->IsInFiducialCut(lv2,fCalorimeter))
598 inacceptance = kTRUE ;
599 if(GetDebug() > 2) printf("In %s fiducial cut acceptance? %d\n",fCalorimeter.Data(),inacceptance);
600 }
601
602 }
603 else if(fCalorimeter == "EMCAL" && GetCaloUtils()->IsEMCALGeoMatrixSet()){
604 if(GetEMCALGeometry()){
605 if(GetEMCALGeometry()->Impact(phot1) && GetEMCALGeometry()->Impact(phot2))
606 inacceptance = kTRUE;
607 if(GetDebug() > 2) printf("In %s Real acceptance? %d\n",fCalorimeter.Data(),inacceptance);
608 }
609 else{
610 if(GetFiducialCut()->IsInFiducialCut(lv1,fCalorimeter) && GetFiducialCut()->IsInFiducialCut(lv2,fCalorimeter))
611 inacceptance = kTRUE ;
612 if(GetDebug() > 2) printf("In %s fiducial cut acceptance? %d\n",fCalorimeter.Data(),inacceptance);
613 }
614 }
615
616 if(inacceptance){
617
618 fhPrimAccPt->Fill(pi0Pt) ;
619 fhPrimAccPhi->Fill(phi) ;
620 fhPrimAccY->Fill(pi0Y) ;
621 Double_t angle = lv1.Angle(lv2.Vect());
622 fhPrimOpeningAngle ->Fill(pi0Pt,angle);
623 fhPrimCosOpeningAngle->Fill(pi0Pt,TMath::Cos(angle));
624
625 }//Accepted
626 }// 2 photons
627 }//Check daughters exist
628 }// Primary pi0
50f39b97 629 }//loop on primaries
630 }//stack exists and data is MC
591cc579 631 }//read stack
632 else if(GetReader()->ReadAODMCParticles()){
50f39b97 633 if(GetDebug() >= 0) printf("AliAnaPi0::MakeAnalysisFillHistograms() - Acceptance calculation with MCParticles not implemented yet\n");
591cc579 634 }
1c5acb87 635}
636
a5cc4f03 637//________________________________________________________________________
638void AliAnaPi0::ReadHistograms(TList* outputList)
639{
50f39b97 640 // Needed when Terminate is executed in distributed environment
641 // Refill analysis histograms of this class with corresponding histograms in output list.
642
643 // Histograms of this analsys are kept in the same list as other analysis, recover the position of
644 // the first one and then add the next.
645 Int_t index = outputList->IndexOf(outputList->FindObject(GetAddedHistogramsStringToName()+"hRe_cen0_pid0_dist1"));
646
647 if(!fhRe1) fhRe1 = new TH3D*[fNCentrBin*fNPID] ;
648 if(!fhRe2) fhRe2 = new TH3D*[fNCentrBin*fNPID] ;
649 if(!fhRe3) fhRe3 = new TH3D*[fNCentrBin*fNPID] ;
650 if(!fhMi1) fhMi1 = new TH3D*[fNCentrBin*fNPID] ;
651 if(!fhMi2) fhMi2 = new TH3D*[fNCentrBin*fNPID] ;
652 if(!fhMi3) fhMi3 = new TH3D*[fNCentrBin*fNPID] ;
653 if(!fhReMod) fhReMod = new TH3D*[fNModules] ;
654
655 for(Int_t ic=0; ic<fNCentrBin; ic++){
656 for(Int_t ipid=0; ipid<fNPID; ipid++){
657 fhRe1[ic*fNPID+ipid] = (TH3D*) outputList->At(index++);
658 fhMi1[ic*fNPID+ipid] = (TH3D*) outputList->At(index++);
659 fhRe2[ic*fNPID+ipid] = (TH3D*) outputList->At(index++);
660 fhMi2[ic*fNPID+ipid] = (TH3D*) outputList->At(index++);
661 fhRe3[ic*fNPID+ipid] = (TH3D*) outputList->At(index++);
662 fhMi3[ic*fNPID+ipid] = (TH3D*) outputList->At(index++);
663 }
664 }
665
666 fhEvents = (TH3D *) outputList->At(index++);
667
668 //Histograms filled only if MC data is requested
669 if(IsDataMC() || (GetReader()->GetDataType() == AliCaloTrackReader::kMC) ){
670 fhPrimPt = (TH1D*) outputList->At(index++);
671 fhPrimAccPt = (TH1D*) outputList->At(index++);
672 fhPrimY = (TH1D*) outputList->At(index++);
673 fhPrimAccY = (TH1D*) outputList->At(index++);
674 fhPrimPhi = (TH1D*) outputList->At(index++);
675 fhPrimAccPhi = (TH1D*) outputList->At(index++);
676 }
677
678 for(Int_t imod=0; imod < fNModules; imod++)
679 fhReMod[imod] = (TH3D*) outputList->At(index++);
680
a5cc4f03 681}
682
683
6639984f 684//____________________________________________________________________________________________________________________________________________________
a5cc4f03 685void AliAnaPi0::Terminate(TList* outputList)
6639984f 686{
687 //Do some calculations and plots from the final histograms.
477d6cee 688
fbeaf916 689 printf(" *** %s Terminate:\n", GetName()) ;
50f39b97 690
a5cc4f03 691 //Recover histograms from output histograms list, needed for distributed analysis.
692 ReadHistograms(outputList);
50f39b97 693
2e557d1c 694 if (!fhRe1) {
50f39b97 695 printf("AliAnaPi0::Terminate() - Error: Remote output histograms not imported in AliAnaPi0 object");
696 return;
2e557d1c 697 }
50f39b97 698
a3aebfff 699 printf("AliAnaPi0::Terminate() Mgg Real : %5.3f , RMS : %5.3f \n", fhRe1[0]->GetMean(), fhRe1[0]->GetRMS() ) ;
50f39b97 700
71dd883b 701 char nameIM[128];
702 sprintf(nameIM,"AliAnaPi0_%s_cPt",fCalorimeter.Data());
703 TCanvas * cIM = new TCanvas(nameIM, "", 400, 10, 600, 700) ;
6639984f 704 cIM->Divide(2, 2);
50f39b97 705
6639984f 706 cIM->cd(1) ;
707 //gPad->SetLogy();
583c48f1 708 TH1D * hIMAllPt = (TH1D*) fhRe1[0]->ProjectionZ(Form("IMPtAll_%s",fCalorimeter.Data()));
6639984f 709 hIMAllPt->SetLineColor(2);
710 hIMAllPt->SetTitle("No cut on p_{T, #gamma#gamma} ");
711 hIMAllPt->Draw();
712
713 cIM->cd(2) ;
583c48f1 714 TH3F * hRe1Pt5 = (TH3F*)fhRe1[0]->Clone(Form("IMPt5_%s",fCalorimeter.Data()));
6639984f 715 hRe1Pt5->GetXaxis()->SetRangeUser(0,5);
583c48f1 716 TH1D * hIMPt5 = (TH1D*) hRe1Pt5->Project3D(Form("IMPt5_%s_pz",fCalorimeter.Data()));
6639984f 717 hIMPt5->SetLineColor(2);
718 hIMPt5->SetTitle("0 < p_{T, #gamma#gamma} < 5 GeV/c");
719 hIMPt5->Draw();
720
721 cIM->cd(3) ;
583c48f1 722 TH3F * hRe1Pt10 = (TH3F*)fhRe1[0]->Clone(Form("IMPt10_%s",fCalorimeter.Data()));
6639984f 723 hRe1Pt10->GetXaxis()->SetRangeUser(5,10);
583c48f1 724 TH1D * hIMPt10 = (TH1D*) hRe1Pt10->Project3D(Form("IMPt10_%s_pz",fCalorimeter.Data()));
6639984f 725 hIMPt10->SetLineColor(2);
726 hIMPt10->SetTitle("5 < p_{T, #gamma#gamma} < 10 GeV/c");
727 hIMPt10->Draw();
728
729 cIM->cd(4) ;
583c48f1 730 TH3F * hRe1Pt20 = (TH3F*)fhRe1[0]->Clone(Form("IMPt20_%s",fCalorimeter.Data()));
6639984f 731 hRe1Pt20->GetXaxis()->SetRangeUser(10,20);
583c48f1 732 TH1D * hIMPt20 = (TH1D*) hRe1Pt20->Project3D(Form("IMPt20_%s_pz",fCalorimeter.Data()));
6639984f 733 hIMPt20->SetLineColor(2);
734 hIMPt20->SetTitle("10 < p_{T, #gamma#gamma} < 20 GeV/c");
735 hIMPt20->Draw();
736
71dd883b 737 char nameIMF[128];
738 sprintf(nameIMF,"AliAnaPi0_%s_Mgg.eps",fCalorimeter.Data());
739 cIM->Print(nameIMF);
6639984f 740
71dd883b 741 char namePt[128];
742 sprintf(namePt,"AliAnaPi0_%s_cPt",fCalorimeter.Data());
743 TCanvas * cPt = new TCanvas(namePt, "", 400, 10, 600, 700) ;
6639984f 744 cPt->Divide(2, 2);
745
746 cPt->cd(1) ;
747 //gPad->SetLogy();
748 TH1D * hPt = (TH1D*) fhRe1[0]->Project3D("x");
749 hPt->SetLineColor(2);
750 hPt->SetTitle("No cut on M_{#gamma#gamma} ");
751 hPt->Draw();
752
753 cPt->cd(2) ;
583c48f1 754 TH3F * hRe1IM1 = (TH3F*)fhRe1[0]->Clone(Form("Pt1_%s",fCalorimeter.Data()));
6639984f 755 hRe1IM1->GetZaxis()->SetRangeUser(0.05,0.21);
756 TH1D * hPtIM1 = (TH1D*) hRe1IM1->Project3D("x");
757 hPtIM1->SetLineColor(2);
758 hPtIM1->SetTitle("0.05 < M_{#gamma#gamma} < 0.21 GeV/c^{2}");
759 hPtIM1->Draw();
760
761 cPt->cd(3) ;
583c48f1 762 TH3F * hRe1IM2 = (TH3F*)fhRe1[0]->Clone(Form("Pt2_%s",fCalorimeter.Data()));
6639984f 763 hRe1IM2->GetZaxis()->SetRangeUser(0.09,0.17);
764 TH1D * hPtIM2 = (TH1D*) hRe1IM2->Project3D("x");
765 hPtIM2->SetLineColor(2);
766 hPtIM2->SetTitle("0.09 < M_{#gamma#gamma} < 0.17 GeV/c^{2}");
767 hPtIM2->Draw();
768
769 cPt->cd(4) ;
583c48f1 770 TH3F * hRe1IM3 = (TH3F*)fhRe1[0]->Clone(Form("Pt3_%s",fCalorimeter.Data()));
6639984f 771 hRe1IM3->GetZaxis()->SetRangeUser(0.11,0.15);
772 TH1D * hPtIM3 = (TH1D*) hRe1IM1->Project3D("x");
773 hPtIM3->SetLineColor(2);
774 hPtIM3->SetTitle("0.11 < M_{#gamma#gamma} < 0.15 GeV/c^{2}");
775 hPtIM3->Draw();
776
71dd883b 777 char namePtF[128];
778 sprintf(namePtF,"AliAnaPi0_%s_Pt.eps",fCalorimeter.Data());
779 cPt->Print(namePtF);
1c5acb87 780
6639984f 781
782 char line[1024] ;
71dd883b 783 sprintf(line, ".!tar -zcf %s_%s.tar.gz *.eps", GetName(),fCalorimeter.Data()) ;
6639984f 784 gROOT->ProcessLine(line);
71dd883b 785 sprintf(line, ".!rm -fR AliAnaPi0_%s*.eps",fCalorimeter.Data());
6639984f 786 gROOT->ProcessLine(line);
787
71dd883b 788 printf(" AliAnaPi0::Terminate() - !! All the eps files are in %s_%s.tar.gz !!!\n", GetName(), fCalorimeter.Data());
1c5acb87 789
6639984f 790}
c8fe2783 791 //____________________________________________________________________________________________________________________________________________________
792Int_t AliAnaPi0::GetEventIndex(AliAODPWG4Particle * part, Double_t * vert)
793{
794 // retieves the event index and checks the vertex
795 // in the mixed buffer returns -2 if vertex NOK
796 // for normal events returns 0 if vertex OK and -1 if vertex NOK
797
798 Int_t rv = -1 ;
799 if (GetMixedEvent()){
800 TObjArray * pl = 0x0;
801 if (part->GetDetector().Contains("PHOS")) {
802 pl = GetAODPHOS();
803 } else if (part->GetDetector().Contains("EMCAL")) {
804 pl = GetAODEMCAL();
805 } else {
806 AliFatal(Form("%s is an unknown calorimeter", part->GetDetector().Data())) ;
807 }
808 rv = GetMixedEvent()->EventIndexForCaloCluster(part->GetCaloLabel(0)) ;
809 GetMixedEvent()->GetVertexOfEvent(rv)->GetXYZ(vert);
810 if(vert[2] < -fZvtxCut || vert[2] > fZvtxCut)
811 rv = -2 ; //Event can not be used (vertex, centrality,... cuts not fulfilled)
0ae57829 812 } else if(GetReader()->GetDataType()!=AliCaloTrackReader::kMC){
c8fe2783 813 Double_t * tempo = GetReader()->GetVertex() ;
814 vert[0] = tempo[0] ;
815 vert[1] = tempo[1] ;
816 vert[2] = tempo[2] ;
817 if(vert[2] < -fZvtxCut || vert[2] > fZvtxCut)
818 rv = -1 ; //Event can not be used (vertex, centrality,... cuts not fulfilled)
819 else
820 rv = 0 ;
0ae57829 821 }//No MC reader
822 else rv = 0;
823
c8fe2783 824 return rv ;
825}