AliPhysics  v5-06-21-01 (1eac791)
 All Classes Namespaces Files Functions Variables Enumerations Enumerator Macros
AliAnaPi0.cxx
Go to the documentation of this file.
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 
16 // --- ROOT system ---
17 #include "TH3.h"
18 #include "TH2F.h"
19 //#include "Riostream.h"
20 #include "TCanvas.h"
21 #include "TPad.h"
22 #include "TROOT.h"
23 #include "TClonesArray.h"
24 #include "TObjString.h"
25 #include "TDatabasePDG.h"
26 
27 //---- AliRoot system ----
28 #include "AliAnaPi0.h"
29 #include "AliCaloTrackReader.h"
30 #include "AliCaloPID.h"
31 #include "AliStack.h"
32 #include "AliFiducialCut.h"
33 #include "TParticle.h"
34 #include "AliVEvent.h"
35 #include "AliESDCaloCluster.h"
36 #include "AliESDEvent.h"
37 #include "AliAODEvent.h"
38 #include "AliNeutralMesonSelection.h"
39 #include "AliMixedEvent.h"
40 #include "AliAODMCParticle.h"
41 
42 // --- Detectors ---
43 #include "AliPHOSGeoUtils.h"
44 #include "AliEMCALGeometry.h"
45 
47 ClassImp(AliAnaPi0) ;
49 
50 //______________________________________________________
52 //______________________________________________________
53 AliAnaPi0::AliAnaPi0() : AliAnaCaloTrackCorrBaseClass(),
54 fEventsList(0x0),
55 fNModules(22),
56 fUseAngleCut(kFALSE), fUseAngleEDepCut(kFALSE), fAngleCut(0), fAngleMaxCut(7.),
57 fMultiCutAna(kFALSE), fMultiCutAnaSim(kFALSE),
58 fNPtCuts(0), fNAsymCuts(0), fNCellNCuts(0), fNPIDBits(0),
59 fMakeInvPtPlots(kFALSE), fSameSM(kFALSE),
60 fFillSMCombinations(kFALSE), fCheckConversion(kFALSE),
61 fFillBadDistHisto(kFALSE), fFillSSCombinations(kFALSE),
62 fFillAngleHisto(kFALSE), fFillAsymmetryHisto(kFALSE), fFillOriginHisto(0),
63 fFillArmenterosThetaStar(0), fFillOnlyMCAcceptanceHisto(0),
64 fCheckAccInSector(kFALSE),
65 fPhotonMom1(), fPhotonMom1Boost(), fPhotonMom2(), fPi0Mom(),
66 fProdVertex(),
67 
68 // Histograms
69 fhAverTotECluster(0), fhAverTotECell(0), fhAverTotECellvsCluster(0),
70 fhEDensityCluster(0), fhEDensityCell(0), fhEDensityCellvsCluster(0),
71 fhReMod(0x0), fhReSameSideEMCALMod(0x0), fhReSameSectorEMCALMod(0x0), fhReDiffPHOSMod(0x0),
72 fhMiMod(0x0), fhMiSameSideEMCALMod(0x0), fhMiSameSectorEMCALMod(0x0), fhMiDiffPHOSMod(0x0),
73 fhReConv(0x0), fhMiConv(0x0), fhReConv2(0x0), fhMiConv2(0x0),
74 fhRe1(0x0), fhMi1(0x0), fhRe2(0x0), fhMi2(0x0),
75 fhRe3(0x0), fhMi3(0x0), fhReInvPt1(0x0), fhMiInvPt1(0x0),
76 fhReInvPt2(0x0), fhMiInvPt2(0x0), fhReInvPt3(0x0), fhMiInvPt3(0x0),
77 fhRePtNCellAsymCuts(0x0), fhMiPtNCellAsymCuts(0x0), fhRePtNCellAsymCutsSM(),
78 fhRePIDBits(0x0), fhRePtMult(0x0), fhReSS(),
79 fhRePtAsym(0x0), fhRePtAsymPi0(0x0), fhRePtAsymEta(0x0),
80 fhEventBin(0), fhEventMixBin(0),
81 fhCentrality(0x0), fhCentralityNoPair(0x0),
82 fhEventPlaneResolution(0x0),
83 fhRealOpeningAngle(0x0), fhRealCosOpeningAngle(0x0), fhMixedOpeningAngle(0x0), fhMixedCosOpeningAngle(0x0),
84 
85 // MC histograms
86 fhPrimPi0E(0x0), fhPrimPi0Pt(0x0),
87 fhPrimPi0AccE(0x0), fhPrimPi0AccPt(0x0),
88 fhPrimPi0Y(0x0), fhPrimPi0AccY(0x0),
89 fhPrimPi0Yeta(0x0), fhPrimPi0YetaYcut(0x0), fhPrimPi0AccYeta(0x0),
90 fhPrimPi0Phi(0x0), fhPrimPi0AccPhi(0x0),
91 fhPrimPi0OpeningAngle(0x0), fhPrimPi0OpeningAngleAsym(0x0),fhPrimPi0CosOpeningAngle(0x0),
92 fhPrimPi0PtCentrality(0), fhPrimPi0PtEventPlane(0),
93 fhPrimPi0AccPtCentrality(0), fhPrimPi0AccPtEventPlane(0),
94 fhPrimEtaE(0x0), fhPrimEtaPt(0x0),
95 fhPrimEtaAccE(0x0), fhPrimEtaAccPt(0x0),
96 fhPrimEtaY(0x0), fhPrimEtaAccY(0x0),
97 fhPrimEtaYeta(0x0), fhPrimEtaYetaYcut(0x0), fhPrimEtaAccYeta(0x0),
98 fhPrimEtaPhi(0x0), fhPrimEtaAccPhi(0x0),
99 fhPrimEtaOpeningAngle(0x0), fhPrimEtaOpeningAngleAsym(0x0),fhPrimEtaCosOpeningAngle(0x0),
100 fhPrimEtaPtCentrality(0), fhPrimEtaPtEventPlane(0),
101 fhPrimEtaAccPtCentrality(0), fhPrimEtaAccPtEventPlane(0),
102 fhPrimPi0PtOrigin(0x0), fhPrimEtaPtOrigin(0x0),
103 fhPrimNotResonancePi0PtOrigin(0x0), fhPrimPi0PtStatus(0x0),
104 fhMCOrgMass(), fhMCOrgAsym(), fhMCOrgDeltaEta(), fhMCOrgDeltaPhi(),
105 fhMCPi0MassPtRec(), fhMCPi0MassPtTrue(), fhMCPi0PtTruePtRec(),
106 fhMCEtaMassPtRec(), fhMCEtaMassPtTrue(), fhMCEtaPtTruePtRec(),
107 fhMCPi0PtOrigin(0x0), fhMCEtaPtOrigin(0x0),
108 fhMCNotResonancePi0PtOrigin(0),fhMCPi0PtStatus(0x0),
109 fhMCPi0ProdVertex(0), fhMCEtaProdVertex(0),
110 fhPrimPi0ProdVertex(0), fhPrimEtaProdVertex(0),
111 fhReMCFromConversion(0), fhReMCFromNotConversion(0), fhReMCFromMixConversion(0),
112 fhCosThStarPrimPi0(0), fhCosThStarPrimEta(0),
113 fhEPairDiffTime(0)
114 {
115  InitParameters();
116 
117  for(Int_t i = 0; i < 4; i++)
118  {
119  fhArmPrimEta[i] = 0;
120  fhArmPrimPi0[i] = 0;
121  }
122 }
123 
124 //_____________________
126 //_____________________
127 AliAnaPi0::~AliAnaPi0()
128 {
129  // Remove event containers
130 
131  if(DoOwnMix() && fEventsList)
132  {
133  for(Int_t ic=0; ic<GetNCentrBin(); ic++)
134  {
135  for(Int_t iz=0; iz<GetNZvertBin(); iz++)
136  {
137  for(Int_t irp=0; irp<GetNRPBin(); irp++)
138  {
139  Int_t bin = GetEventMixBin(ic,iz,irp);
140  fEventsList[bin]->Delete() ;
141  delete fEventsList[bin] ;
142  }
143  }
144  }
145  delete[] fEventsList;
146  }
147 }
148 
149 //______________________________
152 //______________________________
153 void AliAnaPi0::InitParameters()
154 {
155  SetInputAODName("PWG4Particle");
156 
157  AddToHistogramsName("AnaPi0_");
158 
159  fUseAngleCut = kFALSE;
160  fUseAngleEDepCut = kFALSE;
161  fAngleCut = 0.;
162  fAngleMaxCut = TMath::Pi();
163 
164  fMultiCutAna = kFALSE;
165 
166  fNPtCuts = 3;
167  fPtCuts[0] = 0.; fPtCuts[1] = 0.3; fPtCuts[2] = 0.5;
168  for(Int_t i = fNPtCuts; i < 10; i++)fPtCuts[i] = 0.;
169 
170  fNAsymCuts = 2;
171  fAsymCuts[0] = 1.; fAsymCuts[1] = 0.7; //fAsymCuts[2] = 0.6; // fAsymCuts[3] = 0.1;
172  for(Int_t i = fNAsymCuts; i < 10; i++)fAsymCuts[i] = 0.;
173 
174  fNCellNCuts = 3;
175  fCellNCuts[0] = 0; fCellNCuts[1] = 1; fCellNCuts[2] = 2;
176  for(Int_t i = fNCellNCuts; i < 10; i++)fCellNCuts[i] = 0;
177 
178  fNPIDBits = 2;
179  fPIDBits[0] = 0; fPIDBits[1] = 2; // fPIDBits[2] = 4; fPIDBits[3] = 6;// check, no cut, dispersion, neutral, dispersion&&neutral
180  for(Int_t i = fNPIDBits; i < 10; i++)fPIDBits[i] = 0;
181 }
182 
183 //_______________________________________
185 //_______________________________________
186 TObjString * AliAnaPi0::GetAnalysisCuts()
187 {
188  TString parList ; //this will be list of parameters used for this analysis.
189  const Int_t buffersize = 255;
190  char onePar[buffersize] ;
191  snprintf(onePar,buffersize,"--- AliAnaPi0 ---:") ;
192  parList+=onePar ;
193  snprintf(onePar,buffersize,"Number of bins in Centrality: %d;",GetNCentrBin()) ;
194  parList+=onePar ;
195  snprintf(onePar,buffersize,"Number of bins in Z vert. pos: %d;",GetNZvertBin()) ;
196  parList+=onePar ;
197  snprintf(onePar,buffersize,"Number of bins in Reac. Plain: %d;",GetNRPBin()) ;
198  parList+=onePar ;
199  snprintf(onePar,buffersize,"Depth of event buffer: %d;",GetNMaxEvMix()) ;
200  parList+=onePar ;
201  snprintf(onePar,buffersize,"Select pairs with their angle: %d, edep %d, min angle %2.3f, max angle %2.3f;",fUseAngleCut, fUseAngleEDepCut,fAngleCut,fAngleMaxCut) ;
202  parList+=onePar ;
203  snprintf(onePar,buffersize," Asymmetry cuts: n = %d, asymmetry < ",fNAsymCuts) ;
204  for(Int_t i = 0; i < fNAsymCuts; i++) snprintf(onePar,buffersize,"%s %2.2f;",onePar,fAsymCuts[i]);
205  parList+=onePar ;
206  snprintf(onePar,buffersize," PID selection bits: n = %d, PID bit =",fNPIDBits) ;
207  for(Int_t i = 0; i < fNPIDBits; i++) snprintf(onePar,buffersize,"%s %d;",onePar,fPIDBits[i]);
208  parList+=onePar ;
209  snprintf(onePar,buffersize,"Cuts:") ;
210  parList+=onePar ;
211  snprintf(onePar,buffersize,"Z vertex position: -%f < z < %f;",GetZvertexCut(),GetZvertexCut()) ;
212  parList+=onePar ;
213  snprintf(onePar,buffersize,"Calorimeter: %s;",GetCalorimeterString().Data()) ;
214  parList+=onePar ;
215  snprintf(onePar,buffersize,"Number of modules: %d:",fNModules) ;
216  parList+=onePar ;
217  if(fMultiCutAna)
218  {
219  snprintf(onePar, buffersize," pT cuts: n = %d, pt > ",fNPtCuts) ;
220  for(Int_t i = 0; i < fNPtCuts; i++) snprintf(onePar,buffersize,"%s %2.2f;",onePar,fPtCuts[i]);
221  parList+=onePar ;
222  snprintf(onePar,buffersize, " N cell in cluster cuts: n = %d, nCell > ",fNCellNCuts) ;
223  for(Int_t i = 0; i < fNCellNCuts; i++) snprintf(onePar,buffersize,"%s %d;",onePar,fCellNCuts[i]);
224  parList+=onePar ;
225  }
226 
227  return new TObjString(parList) ;
228 }
229 
230 //_________________________________________
233 //_________________________________________
234 TList * AliAnaPi0::GetCreateOutputObjects()
235 {
236  TList * outputContainer = new TList() ;
237  outputContainer->SetName(GetName());
238 
239  // Set sizes and ranges
240  const Int_t buffersize = 255;
241  char key[buffersize] ;
242  char title[buffersize] ;
243 
244  Int_t nptbins = GetHistogramRanges()->GetHistoPtBins();
245  Int_t nphibins = GetHistogramRanges()->GetHistoPhiBins();
246  Int_t netabins = GetHistogramRanges()->GetHistoEtaBins();
247  Float_t ptmax = GetHistogramRanges()->GetHistoPtMax();
248  Float_t phimax = GetHistogramRanges()->GetHistoPhiMax();
249  Float_t etamax = GetHistogramRanges()->GetHistoEtaMax();
250  Float_t ptmin = GetHistogramRanges()->GetHistoPtMin();
251  Float_t phimin = GetHistogramRanges()->GetHistoPhiMin();
252  Float_t etamin = GetHistogramRanges()->GetHistoEtaMin();
253 
254  Int_t nmassbins = GetHistogramRanges()->GetHistoMassBins();
255  Int_t nasymbins = GetHistogramRanges()->GetHistoAsymmetryBins();
256  Float_t massmax = GetHistogramRanges()->GetHistoMassMax();
257  Float_t asymmax = GetHistogramRanges()->GetHistoAsymmetryMax();
258  Float_t massmin = GetHistogramRanges()->GetHistoMassMin();
259  Float_t asymmin = GetHistogramRanges()->GetHistoAsymmetryMin();
260  Int_t ntrmbins = GetHistogramRanges()->GetHistoTrackMultiplicityBins();
261  Int_t ntrmmax = GetHistogramRanges()->GetHistoTrackMultiplicityMax();
262  Int_t ntrmmin = GetHistogramRanges()->GetHistoTrackMultiplicityMin();
263  Int_t tdbins = GetHistogramRanges()->GetHistoDiffTimeBins() ;
264  Float_t tdmax = GetHistogramRanges()->GetHistoDiffTimeMax();
265  Float_t tdmin = GetHistogramRanges()->GetHistoDiffTimeMin();
266 
267  // Start with pure MC kinematics histograms
268  // In case other tasks just need this info like AliAnaPi0EbE
269  if(IsDataMC())
270  {
271  // Pi0
272 
273  fhPrimPi0E = new TH1F("hPrimPi0E","Primary #pi^{0} E, |#it{Y}|<1",nptbins,ptmin,ptmax) ;
274  fhPrimPi0AccE = new TH1F("hPrimPi0AccE","Primary #pi^{0} #it{E} with both photons in acceptance",nptbins,ptmin,ptmax) ;
275  fhPrimPi0E ->SetXTitle("#it{E} (GeV)");
276  fhPrimPi0AccE->SetXTitle("#it{E} (GeV)");
277  outputContainer->Add(fhPrimPi0E) ;
278  outputContainer->Add(fhPrimPi0AccE) ;
279 
280  fhPrimPi0Pt = new TH1F("hPrimPi0Pt","Primary #pi^{0} #it{p}_{T} , |#it{Y}|<1",nptbins,ptmin,ptmax) ;
281  fhPrimPi0AccPt = new TH1F("hPrimPi0AccPt","Primary #pi^{0} #it{p}_{T} with both photons in acceptance",nptbins,ptmin,ptmax) ;
282  fhPrimPi0Pt ->SetXTitle("#it{p}_{T} (GeV/#it{c})");
283  fhPrimPi0AccPt->SetXTitle("#it{p}_{T} (GeV/#it{c})");
284  outputContainer->Add(fhPrimPi0Pt) ;
285  outputContainer->Add(fhPrimPi0AccPt) ;
286 
287  Int_t netabinsopen = TMath::Nint(netabins*4/(etamax-etamin));
288  fhPrimPi0Y = new TH2F("hPrimPi0Rapidity","Rapidity of primary #pi^{0}",nptbins,ptmin,ptmax,netabinsopen,-2, 2) ;
289  fhPrimPi0Y ->SetYTitle("#it{Rapidity}");
290  fhPrimPi0Y ->SetXTitle("#it{p}_{T} (GeV/#it{c})");
291  outputContainer->Add(fhPrimPi0Y) ;
292 
293  fhPrimPi0Yeta = new TH2F("hPrimPi0PseudoRapidity","PseudoRapidity of primary #pi^{0}",nptbins,ptmin,ptmax,netabinsopen,-2, 2) ;
294  fhPrimPi0Yeta ->SetYTitle("#eta");
295  fhPrimPi0Yeta ->SetXTitle("#it{p}_{T} (GeV/#it{c})");
296  outputContainer->Add(fhPrimPi0Yeta) ;
297 
298  fhPrimPi0YetaYcut = new TH2F("hPrimPi0PseudoRapidityYcut","PseudoRapidity of primary #pi^{0}, |#it{Y}|<1",nptbins,ptmin,ptmax,netabinsopen,-2, 2) ;
299  fhPrimPi0YetaYcut ->SetYTitle("#eta");
300  fhPrimPi0YetaYcut ->SetXTitle("#it{p}_{T} (GeV/#it{c})");
301  outputContainer->Add(fhPrimPi0YetaYcut) ;
302 
303  fhPrimPi0AccY = new TH2F("hPrimPi0AccRapidity","Rapidity of primary #pi^{0} with accepted daughters",nptbins,ptmin,ptmax,netabins,etamin,etamax) ;
304  fhPrimPi0AccY->SetYTitle("Rapidity");
305  fhPrimPi0AccY->SetXTitle("#it{p}_{T} (GeV/#it{c})");
306  outputContainer->Add(fhPrimPi0AccY) ;
307 
308  fhPrimPi0AccYeta = new TH2F("hPrimPi0AccPseudoRapidity","PseudoRapidity of primary #pi^{0} with accepted daughters",nptbins,ptmin,ptmax,netabins,etamin,etamax) ;
309  fhPrimPi0AccYeta ->SetYTitle("#eta");
310  fhPrimPi0AccYeta ->SetXTitle("#it{p}_{T} (GeV/#it{c})");
311  outputContainer->Add(fhPrimPi0AccYeta) ;
312 
313  Int_t nphibinsopen = TMath::Nint(nphibins*TMath::TwoPi()/(phimax-phimin));
314  fhPrimPi0Phi = new TH2F("hPrimPi0Phi","#phi of primary #pi^{0}, |#it{Y}|<1",nptbins,ptmin,ptmax,nphibinsopen,0,360) ;
315  fhPrimPi0Phi->SetYTitle("#phi (deg)");
316  fhPrimPi0Phi->SetXTitle("#it{p}_{T} (GeV/#it{c})");
317  outputContainer->Add(fhPrimPi0Phi) ;
318 
319  fhPrimPi0AccPhi = new TH2F("hPrimPi0AccPhi","#phi of primary #pi^{0} with accepted daughters",nptbins,ptmin,ptmax,
320  nphibins,phimin*TMath::RadToDeg(),phimax*TMath::RadToDeg()) ;
321  fhPrimPi0AccPhi->SetYTitle("#phi (deg)");
322  fhPrimPi0AccPhi->SetXTitle("#it{p}_{T} (GeV/#it{c})");
323  outputContainer->Add(fhPrimPi0AccPhi) ;
324 
325  // Eta
326 
327  fhPrimEtaE = new TH1F("hPrimEtaE","Primary eta E",nptbins,ptmin,ptmax) ;
328  fhPrimEtaAccE = new TH1F("hPrimEtaAccE","Primary #eta #it{E} with both photons in acceptance",nptbins,ptmin,ptmax) ;
329  fhPrimEtaE ->SetXTitle("#it{E} (GeV)");
330  fhPrimEtaAccE->SetXTitle("#it{E} (GeV)");
331  outputContainer->Add(fhPrimEtaE) ;
332  outputContainer->Add(fhPrimEtaAccE) ;
333 
334  fhPrimEtaPt = new TH1F("hPrimEtaPt","Primary #eta #it{p}_{T}",nptbins,ptmin,ptmax) ;
335  fhPrimEtaAccPt = new TH1F("hPrimEtaAccPt","Primary eta #it{p}_{T} with both photons in acceptance",nptbins,ptmin,ptmax) ;
336  fhPrimEtaPt ->SetXTitle("#it{p}_{T} (GeV/#it{c})");
337  fhPrimEtaAccPt->SetXTitle("#it{p}_{T} (GeV/#it{c})");
338  outputContainer->Add(fhPrimEtaPt) ;
339  outputContainer->Add(fhPrimEtaAccPt) ;
340 
341  fhPrimEtaY = new TH2F("hPrimEtaRapidity","Rapidity of primary #eta",nptbins,ptmin,ptmax,netabinsopen,-2, 2) ;
342  fhPrimEtaY->SetYTitle("#it{Rapidity}");
343  fhPrimEtaY->SetXTitle("#it{p}_{T} (GeV/#it{c})");
344  outputContainer->Add(fhPrimEtaY) ;
345 
346  fhPrimEtaYeta = new TH2F("hPrimEtaPseudoRapidityEta","PseudoRapidity of primary #eta",nptbins,ptmin,ptmax,netabinsopen,-2, 2) ;
347  fhPrimEtaYeta->SetYTitle("#it{Rapidity}");
348  fhPrimEtaYeta->SetXTitle("#it{p}_{T} (GeV/#it{c})");
349  outputContainer->Add(fhPrimEtaYeta) ;
350 
351  fhPrimEtaYetaYcut = new TH2F("hPrimEtaPseudoRapidityEtaYcut","PseudoRapidity of primary #eta, |#it{Y}|<1",nptbins,ptmin,ptmax,netabinsopen,-2, 2) ;
352  fhPrimEtaYetaYcut->SetYTitle("#it{Pseudorapidity}");
353  fhPrimEtaYetaYcut->SetXTitle("#it{p}_{T} (GeV/#it{c})");
354  outputContainer->Add(fhPrimEtaYetaYcut) ;
355 
356  fhPrimEtaAccY = new TH2F("hPrimEtaAccRapidity","Rapidity of primary #eta",nptbins,ptmin,ptmax, netabins,etamin,etamax) ;
357  fhPrimEtaAccY->SetYTitle("#it{Rapidity}");
358  fhPrimEtaAccY->SetXTitle("#it{p}_{T} (GeV/#it{c})");
359  outputContainer->Add(fhPrimEtaAccY) ;
360 
361  fhPrimEtaAccYeta = new TH2F("hPrimEtaAccPseudoRapidity","PseudoRapidity of primary #eta",nptbins,ptmin,ptmax, netabins,etamin,etamax) ;
362  fhPrimEtaAccYeta->SetYTitle("#it{Pseudorapidity}");
363  fhPrimEtaAccYeta->SetXTitle("#it{p}_{T} (GeV/#it{c})");
364  outputContainer->Add(fhPrimEtaAccYeta) ;
365 
366  fhPrimEtaPhi = new TH2F("hPrimEtaPhi","Azimuthal of primary #eta",nptbins,ptmin,ptmax, nphibinsopen,0,360) ;
367  fhPrimEtaPhi->SetYTitle("#phi (deg)");
368  fhPrimEtaPhi->SetXTitle("#it{p}_{T} (GeV/#it{c})");
369  outputContainer->Add(fhPrimEtaPhi) ;
370 
371  fhPrimEtaAccPhi = new TH2F("hPrimEtaAccPhi","Azimuthal of primary #eta with accepted daughters",nptbins,ptmin,ptmax, nphibins,phimin*TMath::RadToDeg(),phimax*TMath::RadToDeg()) ;
372  fhPrimEtaAccPhi->SetYTitle("#phi (deg)");
373  fhPrimEtaAccPhi->SetXTitle("#it{p}_{T} (GeV/#it{c})");
374  outputContainer->Add(fhPrimEtaAccPhi) ;
375 
376  // Create histograms only for PbPb or high multiplicity analysis analysis
377  if( IsHighMultiplicityAnalysisOn() )
378  {
379  fhPrimPi0PtCentrality = new TH2F("hPrimPi0PtCentrality","Primary #pi^{0} #it{p}_{T} vs reco centrality, |#it{Y}|<1",
380  nptbins,ptmin,ptmax, 100, 0, 100) ;
381  fhPrimPi0AccPtCentrality = new TH2F("hPrimPi0AccPtCentrality","Primary #pi^{0} with both photons in acceptance #it{p}_{T} vs reco centrality",
382  nptbins,ptmin,ptmax, 100, 0, 100) ;
383  fhPrimPi0PtCentrality ->SetXTitle("#it{p}_{T} (GeV/#it{c})");
384  fhPrimPi0AccPtCentrality->SetXTitle("#it{p}_{T} (GeV/#it{c})");
385  fhPrimPi0PtCentrality ->SetYTitle("Centrality");
386  fhPrimPi0AccPtCentrality->SetYTitle("Centrality");
387  outputContainer->Add(fhPrimPi0PtCentrality) ;
388  outputContainer->Add(fhPrimPi0AccPtCentrality) ;
389 
390  fhPrimEtaPtCentrality = new TH2F("hPrimEtaPtCentrality","Primary #eta #it{p}_{T} vs reco centrality, |#it{Y}|<1",nptbins,ptmin,ptmax, 100, 0, 100) ;
391  fhPrimEtaAccPtCentrality = new TH2F("hPrimEtaAccPtCentrality","Primary #eta with both photons in acceptance #it{p}_{T} vs reco centrality",nptbins,ptmin,ptmax, 100, 0, 100) ;
392  fhPrimEtaPtCentrality ->SetXTitle("#it{p}_{T} (GeV/#it{c})");
393  fhPrimEtaAccPtCentrality->SetXTitle("#it{p}_{T} (GeV/#it{c})");
394  fhPrimEtaPtCentrality ->SetYTitle("Centrality");
395  fhPrimEtaAccPtCentrality->SetYTitle("Centrality");
396  outputContainer->Add(fhPrimEtaPtCentrality) ;
397  outputContainer->Add(fhPrimEtaAccPtCentrality) ;
398 
399  fhPrimPi0PtEventPlane = new TH2F("hPrimPi0PtEventPlane","Primary #pi^{0} #it{p}_{T} vs reco event plane angle, |#it{Y}|<1",
400  nptbins,ptmin,ptmax, 100, 0, TMath::Pi()) ;
401  fhPrimPi0AccPtEventPlane = new TH2F("hPrimPi0AccPtEventPlane","Primary #pi^{0} with both photons in acceptance #it{p}_{T} vs reco event plane angle",
402  nptbins,ptmin,ptmax, 100, 0, TMath::Pi()) ;
403  fhPrimPi0PtEventPlane ->SetXTitle("#it{p}_{T} (GeV/#it{c})");
404  fhPrimPi0AccPtEventPlane->SetXTitle("#it{p}_{T} (GeV/#it{c})");
405  fhPrimPi0PtEventPlane ->SetYTitle("Event Plane Angle (rad)");
406  fhPrimPi0AccPtEventPlane->SetYTitle("Event Plane Angle (rad)");
407  outputContainer->Add(fhPrimPi0PtEventPlane) ;
408  outputContainer->Add(fhPrimPi0AccPtEventPlane) ;
409 
410  fhPrimEtaPtEventPlane = new TH2F("hPrimEtaPtEventPlane","Primary #eta #it{p}_{T} vs reco event plane angle, |#it{Y}|<1",nptbins,ptmin,ptmax, 100, 0, TMath::Pi()) ;
411  fhPrimEtaAccPtEventPlane = new TH2F("hPrimEtaAccPtEventPlane","Primary #eta with both #gamma_{decay} in acceptance #it{p}_{T} vs reco event plane angle",nptbins,ptmin,ptmax, 100, 0, TMath::Pi()) ;
412  fhPrimEtaPtEventPlane ->SetXTitle("#it{p}_{T} (GeV/#it{c})");
413  fhPrimEtaAccPtEventPlane->SetXTitle("#it{p}_{T} (GeV/#it{c})");
414  fhPrimEtaPtEventPlane ->SetYTitle("Event Plane Angle (rad)");
415  fhPrimEtaAccPtEventPlane->SetYTitle("Event Plane Angle (rad)");
416  outputContainer->Add(fhPrimEtaPtEventPlane) ;
417  outputContainer->Add(fhPrimEtaAccPtEventPlane) ;
418  }
419 
420  if(fFillAngleHisto)
421  {
422  fhPrimPi0OpeningAngle = new TH2F
423  ("hPrimPi0OpeningAngle","Angle between all primary #gamma pair vs E_{#pi^{0}}",nptbins,ptmin,ptmax,100,0,0.5);
424  fhPrimPi0OpeningAngle->SetYTitle("#theta(rad)");
425  fhPrimPi0OpeningAngle->SetXTitle("E_{ #pi^{0}} (GeV)");
426  outputContainer->Add(fhPrimPi0OpeningAngle) ;
427 
428  fhPrimPi0OpeningAngleAsym = new TH2F
429  ("hPrimPi0OpeningAngleAsym","Angle between all primary #gamma pair vs #it{Asymmetry}, #it{p}_{T}>5 GeV/#it{c}",100,0,1,100,0,0.5);
430  fhPrimPi0OpeningAngleAsym->SetXTitle("|A|=| (E_{1}-E_{2}) / (E_{1}+E_{2}) |");
431  fhPrimPi0OpeningAngleAsym->SetYTitle("#theta(rad)");
432  outputContainer->Add(fhPrimPi0OpeningAngleAsym) ;
433 
434  fhPrimPi0CosOpeningAngle = new TH2F
435  ("hPrimPi0CosOpeningAngle","Cosinus of angle between all primary #gamma pair vs E_{#pi^{0}}",nptbins,ptmin,ptmax,100,-1,1);
436  fhPrimPi0CosOpeningAngle->SetYTitle("cos (#theta) ");
437  fhPrimPi0CosOpeningAngle->SetXTitle("E_{ #pi^{0}} (GeV)");
438  outputContainer->Add(fhPrimPi0CosOpeningAngle) ;
439 
440  fhPrimEtaOpeningAngle = new TH2F
441  ("hPrimEtaOpeningAngle","Angle between all primary #gamma pair vs E_{#eta}",nptbins,ptmin,ptmax,100,0,0.5);
442  fhPrimEtaOpeningAngle->SetYTitle("#theta(rad)");
443  fhPrimEtaOpeningAngle->SetXTitle("E_{#eta} (GeV)");
444  outputContainer->Add(fhPrimEtaOpeningAngle) ;
445 
446  fhPrimEtaOpeningAngleAsym = new TH2F
447  ("hPrimEtaOpeningAngleAsym","Angle between all primary #gamma pair vs #it{Asymmetry}, #it{p}_{T}>5 GeV/#it{c}",100,0,1,100,0,0.5);
448  fhPrimEtaOpeningAngleAsym->SetXTitle("|#it{A}|=| (#it{E}_{1}-#it{E}_{2}) / (#it{E}_{1}+#it{E}_{2}) |");
449  fhPrimEtaOpeningAngleAsym->SetYTitle("#theta(rad)");
450  outputContainer->Add(fhPrimEtaOpeningAngleAsym) ;
451 
452 
453  fhPrimEtaCosOpeningAngle = new TH2F
454  ("hPrimEtaCosOpeningAngle","Cosinus of angle between all primary #gamma pair vs E_{#eta}",nptbins,ptmin,ptmax,100,-1,1);
455  fhPrimEtaCosOpeningAngle->SetYTitle("cos (#theta) ");
456  fhPrimEtaCosOpeningAngle->SetXTitle("#it{E}_{ #eta} (GeV)");
457  outputContainer->Add(fhPrimEtaCosOpeningAngle) ;
458  }
459 
460  // Primary origin
461  if(fFillOriginHisto)
462  {
463  // Pi0
464  fhPrimPi0PtOrigin = new TH2F("hPrimPi0PtOrigin","Primary #pi^{0} #it{p}_{T} vs origin",nptbins,ptmin,ptmax,11,0,11) ;
465  fhPrimPi0PtOrigin->SetXTitle("#it{p}_{T} (GeV/#it{c})");
466  fhPrimPi0PtOrigin->SetYTitle("Origin");
467  fhPrimPi0PtOrigin->GetYaxis()->SetBinLabel(1 ,"Status 21");
468  fhPrimPi0PtOrigin->GetYaxis()->SetBinLabel(2 ,"Quark");
469  fhPrimPi0PtOrigin->GetYaxis()->SetBinLabel(3 ,"qq Resonances ");
470  fhPrimPi0PtOrigin->GetYaxis()->SetBinLabel(4 ,"Resonances");
471  fhPrimPi0PtOrigin->GetYaxis()->SetBinLabel(5 ,"#rho");
472  fhPrimPi0PtOrigin->GetYaxis()->SetBinLabel(6 ,"#omega");
473  fhPrimPi0PtOrigin->GetYaxis()->SetBinLabel(7 ,"K");
474  fhPrimPi0PtOrigin->GetYaxis()->SetBinLabel(8 ,"Other");
475  fhPrimPi0PtOrigin->GetYaxis()->SetBinLabel(9 ,"#eta");
476  fhPrimPi0PtOrigin->GetYaxis()->SetBinLabel(10 ,"#eta prime");
477  outputContainer->Add(fhPrimPi0PtOrigin) ;
478 
479  fhPrimNotResonancePi0PtOrigin = new TH2F("hPrimNotResonancePi0PtOrigin","Primary #pi^{0} #it{p}_{T} vs origin",nptbins,ptmin,ptmax,11,0,11) ;
480  fhPrimNotResonancePi0PtOrigin->SetXTitle("#it{p}_{T} (GeV/#it{c})");
481  fhPrimNotResonancePi0PtOrigin->SetYTitle("Origin");
482  fhPrimNotResonancePi0PtOrigin->GetYaxis()->SetBinLabel(1 ,"Status 21");
483  fhPrimNotResonancePi0PtOrigin->GetYaxis()->SetBinLabel(2 ,"Quark");
484  fhPrimNotResonancePi0PtOrigin->GetYaxis()->SetBinLabel(3 ,"qq Resonances");
485  fhPrimNotResonancePi0PtOrigin->GetYaxis()->SetBinLabel(4 ,"Resonances");
486  fhPrimNotResonancePi0PtOrigin->GetYaxis()->SetBinLabel(5 ,"#rho");
487  fhPrimNotResonancePi0PtOrigin->GetYaxis()->SetBinLabel(6 ,"#omega");
488  fhPrimNotResonancePi0PtOrigin->GetYaxis()->SetBinLabel(7 ,"K");
489  fhPrimNotResonancePi0PtOrigin->GetYaxis()->SetBinLabel(8 ,"Other");
490  fhPrimNotResonancePi0PtOrigin->GetYaxis()->SetBinLabel(9 ,"#eta");
491  fhPrimNotResonancePi0PtOrigin->GetYaxis()->SetBinLabel(10 ,"#eta prime");
492  outputContainer->Add(fhPrimNotResonancePi0PtOrigin) ;
493 
494  fhPrimPi0PtStatus = new TH2F("hPrimPi0PtStatus","Primary #pi^{0} #it{p}_{T} vs status",nptbins,ptmin,ptmax,101,-50,50) ;
495  fhPrimPi0PtStatus->SetXTitle("#it{p}_{T} (GeV/#it{c})");
496  fhPrimPi0PtStatus->SetYTitle("Status");
497  outputContainer->Add(fhPrimPi0PtStatus) ;
498 
499  // Eta
500  fhPrimEtaPtOrigin = new TH2F("hPrimEtaPtOrigin","Primary #pi^{0} #it{p}_{T} vs origin",nptbins,ptmin,ptmax,7,0,7) ;
501  fhPrimEtaPtOrigin->SetXTitle("#it{p}_{T} (GeV/#it{c})");
502  fhPrimEtaPtOrigin->SetYTitle("Origin");
503  fhPrimEtaPtOrigin->GetYaxis()->SetBinLabel(1 ,"Status 21");
504  fhPrimEtaPtOrigin->GetYaxis()->SetBinLabel(2 ,"Quark");
505  fhPrimEtaPtOrigin->GetYaxis()->SetBinLabel(3 ,"qq Resonances");
506  fhPrimEtaPtOrigin->GetYaxis()->SetBinLabel(4 ,"Resonances");
507  fhPrimEtaPtOrigin->GetYaxis()->SetBinLabel(5 ,"Other");
508  fhPrimEtaPtOrigin->GetYaxis()->SetBinLabel(6 ,"#eta prime ");
509  outputContainer->Add(fhPrimEtaPtOrigin) ;
510 
511  // Production vertex
512  fhPrimPi0ProdVertex = new TH2F("hPrimPi0ProdVertex","generated #pi^{0} #it{p}_{T} vs production vertex",
513  200,0.,20.,5000,0,500) ;
514  fhPrimPi0ProdVertex->SetXTitle("#it{p}_{T} (GeV/#it{c})");
515  fhPrimPi0ProdVertex->SetYTitle("#it{R} (cm)");
516  outputContainer->Add(fhPrimPi0ProdVertex) ;
517 
518  fhPrimEtaProdVertex = new TH2F("hPrimEtaProdVertex","generated #eta #it{p}_{T} vs production vertex",
519  200,0.,20.,5000,0,500) ;
520  fhPrimEtaProdVertex->SetXTitle("#it{p}_{T} (GeV/#it{c})");
521  fhPrimEtaProdVertex->SetYTitle("#it{R} (cm)");
522  outputContainer->Add(fhPrimEtaProdVertex) ;
523  }
524 
525  if(fFillArmenterosThetaStar)
526  {
527  TString ebin[] = {"8 < E < 12 GeV","12 < E < 16 GeV", "16 < E < 20 GeV", "E > 20 GeV" };
528  Int_t narmbins = 400;
529  Float_t armmin = 0;
530  Float_t armmax = 0.4;
531 
532  for(Int_t i = 0; i < 4; i++)
533  {
534  fhArmPrimPi0[i] = new TH2F(Form("hArmenterosPrimPi0EBin%d",i),
535  Form("Armenteros of primary #pi^{0}, %s",ebin[i].Data()),
536  200, -1, 1, narmbins,armmin,armmax);
537  fhArmPrimPi0[i]->SetYTitle("#it{p}_{T}^{Arm}");
538  fhArmPrimPi0[i]->SetXTitle("#alpha^{Arm}");
539  outputContainer->Add(fhArmPrimPi0[i]) ;
540 
541  fhArmPrimEta[i] = new TH2F(Form("hArmenterosPrimEtaEBin%d",i),
542  Form("Armenteros of primary #eta, %s",ebin[i].Data()),
543  200, -1, 1, narmbins,armmin,armmax);
544  fhArmPrimEta[i]->SetYTitle("#it{p}_{T}^{Arm}");
545  fhArmPrimEta[i]->SetXTitle("#alpha^{Arm}");
546  outputContainer->Add(fhArmPrimEta[i]) ;
547  }
548 
549  // Same as asymmetry ...
550  fhCosThStarPrimPi0 = new TH2F
551  ("hCosThStarPrimPi0","cos(#theta *) for primary #pi^{0}",nptbins,ptmin,ptmax,200,-1,1);
552  fhCosThStarPrimPi0->SetYTitle("cos(#theta *)");
553  fhCosThStarPrimPi0->SetXTitle("E_{ #pi^{0}} (GeV)");
554  outputContainer->Add(fhCosThStarPrimPi0) ;
555 
556  fhCosThStarPrimEta = new TH2F
557  ("hCosThStarPrimEta","cos(#theta *) for primary #eta",nptbins,ptmin,ptmax,200,-1,1);
558  fhCosThStarPrimEta->SetYTitle("cos(#theta *)");
559  fhCosThStarPrimEta->SetXTitle("E_{ #eta} (GeV)");
560  outputContainer->Add(fhCosThStarPrimEta) ;
561  }
562 
563  if(fFillOnlyMCAcceptanceHisto) return outputContainer;
564  }
565 
566  //
567  // Create mixed event containers
568  //
569  fEventsList = new TList*[GetNCentrBin()*GetNZvertBin()*GetNRPBin()] ;
570 
571  for(Int_t ic=0; ic<GetNCentrBin(); ic++)
572  {
573  for(Int_t iz=0; iz<GetNZvertBin(); iz++)
574  {
575  for(Int_t irp=0; irp<GetNRPBin(); irp++)
576  {
577  Int_t bin = GetEventMixBin(ic,iz,irp);
578  fEventsList[bin] = new TList() ;
579  fEventsList[bin]->SetOwner(kFALSE);
580  }
581  }
582  }
583 
584  // Init the number of modules, set in the class AliCalorimeterUtils
585  fNModules = GetCaloUtils()->GetNumberOfSuperModulesUsed();
586  if(GetCalorimeter()==kPHOS && fNModules > 4) fNModules = 4;
587 
588  fhReMod = new TH2F*[fNModules] ;
589  fhMiMod = new TH2F*[fNModules] ;
590 
591  if(GetCalorimeter() == kPHOS)
592  {
593  fhReDiffPHOSMod = new TH2F*[fNModules] ;
594  fhMiDiffPHOSMod = new TH2F*[fNModules] ;
595  }
596  else
597  {
598  fhReSameSectorEMCALMod = new TH2F*[fNModules/2] ;
599  fhReSameSideEMCALMod = new TH2F*[fNModules-2] ;
600  fhMiSameSectorEMCALMod = new TH2F*[fNModules/2] ;
601  fhMiSameSideEMCALMod = new TH2F*[fNModules-2] ;
602  }
603 
604 
605  fhRe1 = new TH2F*[GetNCentrBin()*fNPIDBits*fNAsymCuts] ;
606  fhMi1 = new TH2F*[GetNCentrBin()*fNPIDBits*fNAsymCuts] ;
607  if(fFillBadDistHisto)
608  {
609  fhRe2 = new TH2F*[GetNCentrBin()*fNPIDBits*fNAsymCuts] ;
610  fhRe3 = new TH2F*[GetNCentrBin()*fNPIDBits*fNAsymCuts] ;
611  fhMi2 = new TH2F*[GetNCentrBin()*fNPIDBits*fNAsymCuts] ;
612  fhMi3 = new TH2F*[GetNCentrBin()*fNPIDBits*fNAsymCuts] ;
613  }
614  if(fMakeInvPtPlots)
615  {
616  fhReInvPt1 = new TH2F*[GetNCentrBin()*fNPIDBits*fNAsymCuts] ;
617  fhMiInvPt1 = new TH2F*[GetNCentrBin()*fNPIDBits*fNAsymCuts] ;
618  if(fFillBadDistHisto){
619  fhReInvPt2 = new TH2F*[GetNCentrBin()*fNPIDBits*fNAsymCuts] ;
620  fhReInvPt3 = new TH2F*[GetNCentrBin()*fNPIDBits*fNAsymCuts] ;
621  fhMiInvPt2 = new TH2F*[GetNCentrBin()*fNPIDBits*fNAsymCuts] ;
622  fhMiInvPt3 = new TH2F*[GetNCentrBin()*fNPIDBits*fNAsymCuts] ;
623  }
624  }
625 
626  if(fCheckConversion)
627  {
628  fhReConv = new TH2F("hReConv","Real Pair with one recombined conversion ",nptbins,ptmin,ptmax,nmassbins,massmin,massmax) ;
629  fhReConv->SetXTitle("#it{p}_{T} (GeV/#it{c})");
630  fhReConv->SetYTitle("#it{M}_{#gamma,#gamma} (GeV/#it{c}^{2})");
631  outputContainer->Add(fhReConv) ;
632 
633  fhReConv2 = new TH2F("hReConv2","Real Pair with 2 recombined conversion ",nptbins,ptmin,ptmax,nmassbins,massmin,massmax) ;
634  fhReConv2->SetXTitle("#it{p}_{T} (GeV/#it{c})");
635  fhReConv2->SetYTitle("#it{M}_{#gamma,#gamma} (GeV/#it{c}^{2})");
636  outputContainer->Add(fhReConv2) ;
637 
638  if(DoOwnMix())
639  {
640  fhMiConv = new TH2F("hMiConv","Mixed Pair with one recombined conversion ",nptbins,ptmin,ptmax,nmassbins,massmin,massmax) ;
641  fhMiConv->SetXTitle("#it{p}_{T} (GeV/#it{c})");
642  fhMiConv->SetYTitle("#it{M}_{#gamma,#gamma} (GeV/#it{c}^{2})");
643  outputContainer->Add(fhMiConv) ;
644 
645  fhMiConv2 = new TH2F("hMiConv2","Mixed Pair with 2 recombined conversion ",nptbins,ptmin,ptmax,nmassbins,massmin,massmax) ;
646  fhMiConv2->SetXTitle("#it{p}_{T} (GeV/#it{c})");
647  fhMiConv2->SetYTitle("#it{M}_{#gamma,#gamma} (GeV/#it{c}^{2})");
648  outputContainer->Add(fhMiConv2) ;
649  }
650  }
651 
652  for(Int_t ic=0; ic<GetNCentrBin(); ic++)
653  {
654  for(Int_t ipid=0; ipid<fNPIDBits; ipid++)
655  {
656  for(Int_t iasym=0; iasym<fNAsymCuts; iasym++)
657  {
658  Int_t index = ((ic*fNPIDBits)+ipid)*fNAsymCuts + iasym;
659  //printf("cen %d, pid %d, asy %d, Index %d\n",ic,ipid,iasym,index);
660  //Distance to bad module 1
661  snprintf(key, buffersize,"hRe_cen%d_pidbit%d_asy%d_dist1",ic,ipid,iasym) ;
662  snprintf(title, buffersize,"Real #it{M}_{#gamma#gamma} distr. for centrality=%d, PID bit=%d and asymmetry %1.2f, dist bad 1",
663  ic,fPIDBits[ipid], fAsymCuts[iasym]) ;
664  fhRe1[index] = new TH2F(key,title,nptbins,ptmin,ptmax,nmassbins,massmin,massmax) ;
665  fhRe1[index]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
666  fhRe1[index]->SetYTitle("#it{M}_{#gamma,#gamma} (GeV/#it{c}^{2})");
667  //printf("name: %s\n ",fhRe1[index]->GetName());
668  outputContainer->Add(fhRe1[index]) ;
669 
670  if(fFillBadDistHisto)
671  {
672  //Distance to bad module 2
673  snprintf(key, buffersize,"hRe_cen%d_pidbit%d_asy%d_dist2",ic,ipid,iasym) ;
674  snprintf(title, buffersize,"Real #it{M}_{#gamma#gamma} distr. for centrality=%d, PID bit=%d and asymmetry %1.2f, dist bad 2",
675  ic,fPIDBits[ipid], fAsymCuts[iasym]) ;
676  fhRe2[index] = new TH2F(key,title,nptbins,ptmin,ptmax,nmassbins,massmin,massmax) ;
677  fhRe2[index]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
678  fhRe2[index]->SetYTitle("#it{M}_{#gamma,#gamma} (GeV/#it{c}^{2})");
679  outputContainer->Add(fhRe2[index]) ;
680 
681  //Distance to bad module 3
682  snprintf(key, buffersize,"hRe_cen%d_pidbit%d_asy%d_dist3",ic,ipid,iasym) ;
683  snprintf(title, buffersize,"Real #it{M}_{#gamma#gamma} distr. for centrality=%d, PID bit=%d and asymmetry %1.2f, dist bad 3",
684  ic,fPIDBits[ipid], fAsymCuts[iasym]) ;
685  fhRe3[index] = new TH2F(key,title,nptbins,ptmin,ptmax,nmassbins,massmin,massmax) ;
686  fhRe3[index]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
687  fhRe3[index]->SetYTitle("#it{M}_{#gamma,#gamma} (GeV/#it{c}^{2})");
688  outputContainer->Add(fhRe3[index]) ;
689  }
690 
691  //Inverse pT
692  if(fMakeInvPtPlots)
693  {
694  //Distance to bad module 1
695  snprintf(key, buffersize,"hReInvPt_cen%d_pidbit%d_asy%d_dist1",ic,ipid,iasym) ;
696  snprintf(title, buffersize,"Real #it{M}_{#gamma#gamma} distr. for centrality=%d, PID bit=%d and asymmetry %1.2f, dist bad 1",
697  ic,fPIDBits[ipid], fAsymCuts[iasym]) ;
698  fhReInvPt1[index] = new TH2F(key,title,nptbins,ptmin,ptmax,nmassbins,massmin,massmax) ;
699  fhReInvPt1[index]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
700  fhReInvPt1[index]->SetYTitle("#it{M}_{#gamma,#gamma} (GeV/#it{c}^{2})");
701  outputContainer->Add(fhReInvPt1[index]) ;
702 
703  if(fFillBadDistHisto){
704  //Distance to bad module 2
705  snprintf(key, buffersize,"hReInvPt_cen%d_pidbit%d_asy%d_dist2",ic,ipid,iasym) ;
706  snprintf(title, buffersize,"Real #it{M}_{#gamma#gamma} distr. for centrality=%d, PID bit=%d and asymmetry %1.2f, dist bad 2",
707  ic,fPIDBits[ipid], fAsymCuts[iasym]) ;
708  fhReInvPt2[index] = new TH2F(key,title,nptbins,ptmin,ptmax,nmassbins,massmin,massmax) ;
709  fhReInvPt2[index]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
710  fhReInvPt2[index]->SetYTitle("#it{M}_{#gamma,#gamma} (GeV/#it{c}^{2})");
711  outputContainer->Add(fhReInvPt2[index]) ;
712 
713  //Distance to bad module 3
714  snprintf(key, buffersize,"hReInvPt_cen%d_pidbit%d_asy%d_dist3",ic,ipid,iasym) ;
715  snprintf(title, buffersize,"Real #it{M}_{#gamma#gamma} distr. for centrality=%d, PID bit=%d and asymmetry %1.2f, dist bad 3",
716  ic,fPIDBits[ipid], fAsymCuts[iasym]) ;
717  fhReInvPt3[index] = new TH2F(key,title,nptbins,ptmin,ptmax,nmassbins,massmin,massmax) ;
718  fhReInvPt3[index]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
719  fhReInvPt3[index]->SetYTitle("#it{M}_{#gamma,#gamma} (GeV/#it{c}^{2})");
720  outputContainer->Add(fhReInvPt3[index]) ;
721  }
722  }
723 
724  if(DoOwnMix())
725  {
726  //Distance to bad module 1
727  snprintf(key, buffersize,"hMi_cen%d_pidbit%d_asy%d_dist1",ic,ipid,iasym) ;
728  snprintf(title, buffersize,"Mixed #it{M}_{#gamma#gamma} distr. for centrality=%d, PID bit=%d and asymmetry %1.2f, dist bad 1",
729  ic,fPIDBits[ipid], fAsymCuts[iasym]) ;
730  fhMi1[index] = new TH2F(key,title,nptbins,ptmin,ptmax,nmassbins,massmin,massmax) ;
731  fhMi1[index]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
732  fhMi1[index]->SetYTitle("#it{M}_{#gamma,#gamma} (GeV/#it{c}^{2})");
733  outputContainer->Add(fhMi1[index]) ;
734  if(fFillBadDistHisto){
735  //Distance to bad module 2
736  snprintf(key, buffersize,"hMi_cen%d_pidbit%d_asy%d_dist2",ic,ipid,iasym) ;
737  snprintf(title, buffersize,"Mixed #it{M}_{#gamma#gamma} distr. for centrality=%d, PID bit=%d and asymmetry %1.2f, dist bad 2",
738  ic,fPIDBits[ipid], fAsymCuts[iasym]) ;
739  fhMi2[index] = new TH2F(key,title,nptbins,ptmin,ptmax,nmassbins,massmin,massmax) ;
740  fhMi2[index]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
741  fhMi2[index]->SetYTitle("#it{M}_{#gamma,#gamma} (GeV/#it{c}^{2})");
742  outputContainer->Add(fhMi2[index]) ;
743 
744  //Distance to bad module 3
745  snprintf(key, buffersize,"hMi_cen%d_pidbit%d_asy%d_dist3",ic,ipid,iasym) ;
746  snprintf(title, buffersize,"Mixed #it{M}_{#gamma#gamma} distr. for centrality=%d, PID bit=%d and asymmetry %1.2f, dist bad 3",
747  ic,fPIDBits[ipid], fAsymCuts[iasym]) ;
748  fhMi3[index] = new TH2F(key,title,nptbins,ptmin,ptmax,nmassbins,massmin,massmax) ;
749  fhMi3[index]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
750  fhMi3[index]->SetYTitle("#it{M}_{#gamma,#gamma} (GeV/#it{c}^{2})");
751  outputContainer->Add(fhMi3[index]) ;
752  }
753 
754  // Inverse pT
755  if(fMakeInvPtPlots)
756  {
757  // Distance to bad module 1
758  snprintf(key, buffersize,"hMiInvPt_cen%d_pidbit%d_asy%d_dist1",ic,ipid,iasym) ;
759  snprintf(title, buffersize,"Mixed #it{M}_{#gamma#gamma} distr. for centrality=%d, PID bit=%d and asymmetry %1.2f, dist bad 1",
760  ic,fPIDBits[ipid], fAsymCuts[iasym]) ;
761  fhMiInvPt1[index] = new TH2F(key,title,nptbins,ptmin,ptmax,nmassbins,massmin,massmax) ;
762  fhMiInvPt1[index]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
763  fhMiInvPt1[index]->SetYTitle("#it{M}_{#gamma,#gamma} (GeV/#it{c}^{2})");
764  outputContainer->Add(fhMiInvPt1[index]) ;
765  if(fFillBadDistHisto){
766  // Distance to bad module 2
767  snprintf(key, buffersize,"hMiInvPt_cen%d_pidbit%d_asy%d_dist2",ic,ipid,iasym) ;
768  snprintf(title, buffersize,"Mixed #it{M}_{#gamma#gamma} distr. for centrality=%d, PID bit=%d and asymmetry %1.2f, dist bad 2",
769  ic,fPIDBits[ipid], fAsymCuts[iasym]) ;
770  fhMiInvPt2[index] = new TH2F(key,title,nptbins,ptmin,ptmax,nmassbins,massmin,massmax) ;
771  fhMiInvPt2[index]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
772  fhMiInvPt2[index]->SetYTitle("#it{M}_{#gamma,#gamma} (GeV/#it{c}^{2})");
773  outputContainer->Add(fhMiInvPt2[index]) ;
774 
775  // Distance to bad module 3
776  snprintf(key, buffersize,"hMiInvPt_cen%d_pidbit%d_asy%d_dist3",ic,ipid,iasym) ;
777  snprintf(title, buffersize,"Mixed #it{M}_{#gamma#gamma} distr. for centrality=%d, PID bit=%d and asymmetry %1.2f,dist bad 3",
778  ic,fPIDBits[ipid], fAsymCuts[iasym]) ;
779  fhMiInvPt3[index] = new TH2F(key,title,nptbins,ptmin,ptmax,nmassbins,massmin,massmax) ;
780  fhMiInvPt3[index]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
781  fhMiInvPt3[index]->SetYTitle("#it{M}_{#gamma,#gamma} (GeV/#it{c}^{2})");
782  outputContainer->Add(fhMiInvPt3[index]) ;
783  }
784  }
785  }
786  }
787  }
788  }
789 
790  fhEPairDiffTime = new TH2F("hEPairDiffTime","cluster pair time difference vs #it{p}_{T}",nptbins,ptmin,ptmax, tdbins,tdmin,tdmax);
791  fhEPairDiffTime->SetXTitle("#it{p}_{T,pair} (GeV/#it{c})");
792  fhEPairDiffTime->SetYTitle("#Delta t (ns)");
793  outputContainer->Add(fhEPairDiffTime);
794 
795  if(fFillAsymmetryHisto)
796  {
797  fhRePtAsym = new TH2F("hRePtAsym","#it{Asymmetry} vs #it{p}_{T} , for pairs",nptbins,ptmin,ptmax,nasymbins,asymmin,asymmax) ;
798  fhRePtAsym->SetXTitle("#it{p}_{T} (GeV/#it{c})");
799  fhRePtAsym->SetYTitle("#it{Asymmetry}");
800  outputContainer->Add(fhRePtAsym);
801 
802  fhRePtAsymPi0 = new TH2F("hRePtAsymPi0","#it{Asymmetry} vs #it{p}_{T} , for pairs close to #pi^{0} mass",nptbins,ptmin,ptmax,nasymbins,asymmin,asymmax) ;
803  fhRePtAsymPi0->SetXTitle("#it{p}_{T} (GeV/#it{c})");
804  fhRePtAsymPi0->SetYTitle("Asymmetry");
805  outputContainer->Add(fhRePtAsymPi0);
806 
807  fhRePtAsymEta = new TH2F("hRePtAsymEta","#it{Asymmetry} vs #it{p}_{T} , for pairs close to #eta mass",nptbins,ptmin,ptmax,nasymbins,asymmin,asymmax) ;
808  fhRePtAsymEta->SetXTitle("#it{p}_{T} (GeV/#it{c})");
809  fhRePtAsymEta->SetYTitle("#it{Asymmetry}");
810  outputContainer->Add(fhRePtAsymEta);
811  }
812 
813  if(fMultiCutAna)
814  {
815  fhRePIDBits = new TH2F*[fNPIDBits];
816  for(Int_t ipid=0; ipid<fNPIDBits; ipid++){
817  snprintf(key, buffersize,"hRe_pidbit%d",ipid) ;
818  snprintf(title, buffersize,"Real #it{M}_{#gamma#gamma} distr. for PIDBit=%d",fPIDBits[ipid]) ;
819  fhRePIDBits[ipid] = new TH2F(key,title,nptbins,ptmin,ptmax,nmassbins,massmin,massmax) ;
820  fhRePIDBits[ipid]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
821  fhRePIDBits[ipid]->SetYTitle("#it{M}_{#gamma,#gamma} (GeV/#it{c}^{2})");
822  outputContainer->Add(fhRePIDBits[ipid]) ;
823  }// pid bit loop
824 
825  fhRePtNCellAsymCuts = new TH2F*[fNPtCuts*fNAsymCuts*fNCellNCuts];
826  fhMiPtNCellAsymCuts = new TH2F*[fNPtCuts*fNAsymCuts*fNCellNCuts];
827 
828  if(fFillSMCombinations)
829  for(Int_t iSM = 0; iSM < fNModules; iSM++) fhRePtNCellAsymCutsSM[iSM] = new TH2F*[fNPtCuts*fNAsymCuts*fNCellNCuts];
830 
831 
832  for(Int_t ipt=0; ipt<fNPtCuts; ipt++)
833  {
834  for(Int_t icell=0; icell<fNCellNCuts; icell++)
835  {
836  for(Int_t iasym=0; iasym<fNAsymCuts; iasym++)
837  {
838  snprintf(key, buffersize,"hRe_pt%d_cell%d_asym%d",ipt,icell,iasym) ;
839  snprintf(title, buffersize,"Real #it{M}_{#gamma#gamma} distr. for pt >%2.2f, ncell>%d and asym >%1.2f ",fPtCuts[ipt],fCellNCuts[icell], fAsymCuts[iasym]) ;
840  Int_t index = ((ipt*fNCellNCuts)+icell)*fNAsymCuts + iasym;
841  //printf("ipt %d, icell %d, iassym %d, index %d\n",ipt, icell, iasym, index);
842  fhRePtNCellAsymCuts[index] = new TH2F(key,title,nptbins,ptmin,ptmax,nmassbins,massmin,massmax) ;
843  fhRePtNCellAsymCuts[index]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
844  fhRePtNCellAsymCuts[index]->SetYTitle("#it{M}_{#gamma,#gamma} (GeV/#it{c}^{2})");
845  outputContainer->Add(fhRePtNCellAsymCuts[index]) ;
846 
847  snprintf(key, buffersize,"hMi_pt%d_cell%d_asym%d",ipt,icell,iasym) ;
848  snprintf(title, buffersize,"Mixed #it{M}_{#gamma#gamma} distr. for pt >%2.2f, ncell>%d and asym >%1.2f",fPtCuts[ipt],fCellNCuts[icell], fAsymCuts[iasym]) ;
849  fhMiPtNCellAsymCuts[index] = new TH2F(key,title,nptbins,ptmin,ptmax,nmassbins,massmin,massmax) ;
850  fhMiPtNCellAsymCuts[index]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
851  fhMiPtNCellAsymCuts[index]->SetYTitle("#it{M}_{#gamma,#gamma} (GeV/#it{c}^{2})");
852  outputContainer->Add(fhMiPtNCellAsymCuts[index]) ;
853 
854  if(fFillSMCombinations)
855  {
856  for(Int_t iSM = 0; iSM < fNModules; iSM++)
857  {
858  snprintf(key, buffersize,"hRe_pt%d_cell%d_asym%d_SM%d",ipt,icell,iasym,iSM) ;
859  snprintf(title, buffersize,"Real #it{M}_{#gamma#gamma} distr. for pt >%2.2f, ncell>%d and asym >%1.2f, SM %d ",
860  fPtCuts[ipt],fCellNCuts[icell], fAsymCuts[iasym],iSM) ;
861  fhRePtNCellAsymCutsSM[iSM][index] = new TH2F(key,title,nptbins,ptmin,ptmax,nmassbins,massmin,massmax) ;
862  fhRePtNCellAsymCutsSM[iSM][index]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
863  fhRePtNCellAsymCutsSM[iSM][index]->SetYTitle("#it{M}_{#gamma,#gamma} (GeV/#it{c}^{2})");
864  outputContainer->Add(fhRePtNCellAsymCutsSM[iSM][index]) ;
865  }
866  }
867  }
868  }
869  }
870 
871  if(ntrmbins!=0)
872  {
873  fhRePtMult = new TH3F*[fNAsymCuts] ;
874  for(Int_t iasym = 0; iasym<fNAsymCuts; iasym++)
875  {
876  fhRePtMult[iasym] = new TH3F(Form("hRePtMult_asym%d",iasym),Form("(#it{p}_{T},C,M)_{#gamma#gamma}, A<%1.2f",fAsymCuts[iasym]),
877  nptbins,ptmin,ptmax,ntrmbins,ntrmmin,ntrmmax,nmassbins,massmin,massmax);
878  fhRePtMult[iasym]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
879  fhRePtMult[iasym]->SetYTitle("Track multiplicity");
880  fhRePtMult[iasym]->SetZTitle("#it{M}_{#gamma,#gamma} (GeV/#it{c}^{2})");
881  outputContainer->Add(fhRePtMult[iasym]) ;
882  }
883  }
884  } // multi cuts analysis
885 
886  if(fFillSSCombinations)
887  {
888  fhReSS[0] = new TH2F("hRe_SS_Tight"," 0.01 < #lambda_{0}^{2} < 0.4",
889  nptbins,ptmin,ptmax,nmassbins,massmin,massmax);
890  fhReSS[0]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
891  fhReSS[0]->SetYTitle("#it{M}_{#gamma,#gamma} (GeV/#it{c}^{2})");
892  outputContainer->Add(fhReSS[0]) ;
893 
894 
895  fhReSS[1] = new TH2F("hRe_SS_Loose"," #lambda_{0}^{2} > 0.4",
896  nptbins,ptmin,ptmax,nmassbins,massmin,massmax);
897  fhReSS[1]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
898  fhReSS[1]->SetYTitle("#it{M}_{#gamma,#gamma} (GeV/#it{c}^{2})");
899  outputContainer->Add(fhReSS[1]) ;
900 
901 
902  fhReSS[2] = new TH2F("hRe_SS_Both"," cluster_{1} #lambda_{0}^{2} > 0.4; cluster_{2} 0.01 < #lambda_{0}^{2} < 0.4",
903  nptbins,ptmin,ptmax,nmassbins,massmin,massmax);
904  fhReSS[2]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
905  fhReSS[2]->SetYTitle("#it{M}_{#gamma,#gamma} (GeV/#it{c}^{2})");
906  outputContainer->Add(fhReSS[2]) ;
907  }
908 
909  if(DoOwnMix())
910  {
911  fhEventBin=new TH1I("hEventBin","Number of real pairs per bin(cen,vz,rp)",
912  GetNCentrBin()*GetNZvertBin()*GetNRPBin()+1,0,
913  GetNCentrBin()*GetNZvertBin()*GetNRPBin()+1) ;
914  fhEventBin->SetXTitle("bin");
915  outputContainer->Add(fhEventBin) ;
916 
917 
918  fhEventMixBin=new TH1I("hEventMixBin","Number of mixed pairs per bin(cen,vz,rp)",
919  GetNCentrBin()*GetNZvertBin()*GetNRPBin()+1,0,
920  GetNCentrBin()*GetNZvertBin()*GetNRPBin()+1) ;
921  fhEventMixBin->SetXTitle("bin");
922  outputContainer->Add(fhEventMixBin) ;
923  }
924 
925  if( IsHighMultiplicityAnalysisOn() )
926  {
927  fhCentrality=new TH1F("hCentralityBin","Number of events in centrality bin",GetNCentrBin(),0.,1.*GetNCentrBin()) ;
928  fhCentrality->SetXTitle("Centrality bin");
929  outputContainer->Add(fhCentrality) ;
930 
931  fhCentralityNoPair=new TH1F("hCentralityBinNoPair","Number of events in centrality bin, with no cluster pairs",GetNCentrBin(),0.,1.*GetNCentrBin()) ;
932  fhCentralityNoPair->SetXTitle("Centrality bin");
933  outputContainer->Add(fhCentralityNoPair) ;
934 
935  fhEventPlaneResolution=new TH2F("hEventPlaneResolution","Event plane resolution",GetNCentrBin(),0,GetNCentrBin(),100,0.,TMath::TwoPi()) ;
936  fhEventPlaneResolution->SetYTitle("Resolution");
937  fhEventPlaneResolution->SetXTitle("Centrality Bin");
938  outputContainer->Add(fhEventPlaneResolution) ;
939  }
940 
941  if(fFillAngleHisto)
942  {
943  fhRealOpeningAngle = new TH2F
944  ("hRealOpeningAngle","Angle between all #gamma pair vs E_{#pi^{0}}",nptbins,ptmin,ptmax,300,0,TMath::Pi());
945  fhRealOpeningAngle->SetYTitle("#theta(rad)");
946  fhRealOpeningAngle->SetXTitle("E_{ #pi^{0}} (GeV)");
947  outputContainer->Add(fhRealOpeningAngle) ;
948 
949  fhRealCosOpeningAngle = new TH2F
950  ("hRealCosOpeningAngle","Cosinus of angle between all #gamma pair vs E_{#pi^{0}}",nptbins,ptmin,ptmax,100,0,1);
951  fhRealCosOpeningAngle->SetYTitle("cos (#theta) ");
952  fhRealCosOpeningAngle->SetXTitle("E_{ #pi^{0}} (GeV)");
953  outputContainer->Add(fhRealCosOpeningAngle) ;
954 
955  if(DoOwnMix())
956  {
957  fhMixedOpeningAngle = new TH2F
958  ("hMixedOpeningAngle","Angle between all #gamma pair vs E_{#pi^{0}}, Mixed pairs",nptbins,ptmin,ptmax,300,0,TMath::Pi());
959  fhMixedOpeningAngle->SetYTitle("#theta(rad)");
960  fhMixedOpeningAngle->SetXTitle("E_{ #pi^{0}} (GeV)");
961  outputContainer->Add(fhMixedOpeningAngle) ;
962 
963  fhMixedCosOpeningAngle = new TH2F
964  ("hMixedCosOpeningAngle","Cosinus of angle between all #gamma pair vs E_{#pi^{0}}, Mixed pairs",nptbins,ptmin,ptmax,100,0,1);
965  fhMixedCosOpeningAngle->SetYTitle("cos (#theta) ");
966  fhMixedCosOpeningAngle->SetXTitle("E_{ #pi^{0}} (GeV)");
967  outputContainer->Add(fhMixedCosOpeningAngle) ;
968  }
969  }
970 
971  // Histograms filled only if MC data is requested
972  if(IsDataMC())
973  {
974  fhReMCFromConversion = new TH2F("hReMCFromConversion","Invariant mass of 2 clusters originated in conversions",
975  nptbins,ptmin,ptmax,nmassbins,massmin,massmax);
976  fhReMCFromConversion->SetXTitle("#it{p}_{T} (GeV/#it{c})");
977  fhReMCFromConversion->SetYTitle("#it{M}_{#gamma,#gamma} (GeV/#it{c}^{2})");
978  outputContainer->Add(fhReMCFromConversion) ;
979 
980  fhReMCFromNotConversion = new TH2F("hReMCNotFromConversion","Invariant mass of 2 clusters not originated in conversions",
981  nptbins,ptmin,ptmax,nmassbins,massmin,massmax);
982  fhReMCFromNotConversion->SetXTitle("#it{p}_{T} (GeV/#it{c})");
983  fhReMCFromNotConversion->SetYTitle("#it{M}_{#gamma,#gamma} (GeV/#it{c}^{2})");
984  outputContainer->Add(fhReMCFromNotConversion) ;
985 
986  fhReMCFromMixConversion = new TH2F("hReMCFromMixConversion","Invariant mass of 2 clusters one from conversion and the other not",
987  nptbins,ptmin,ptmax,nmassbins,massmin,massmax);
988  fhReMCFromMixConversion->SetXTitle("#it{p}_{T} (GeV/#it{c})");
989  fhReMCFromMixConversion->SetYTitle("#it{M}_{#gamma,#gamma} (GeV/#it{c}^{2})");
990  outputContainer->Add(fhReMCFromMixConversion) ;
991 
992  if(fFillOriginHisto)
993  {
994 
995  fhMCPi0PtOrigin = new TH2F("hMCPi0PtOrigin","Reconstructed pair from generated #pi^{0} #it{p}_{T} vs origin",nptbins,ptmin,ptmax,11,0,11) ;
996  fhMCPi0PtOrigin->SetXTitle("#it{p}_{T} (GeV/#it{c})");
997  fhMCPi0PtOrigin->SetYTitle("Origin");
998  fhMCPi0PtOrigin->GetYaxis()->SetBinLabel(1 ,"Status 21");
999  fhMCPi0PtOrigin->GetYaxis()->SetBinLabel(2 ,"Quark");
1000  fhMCPi0PtOrigin->GetYaxis()->SetBinLabel(3 ,"qq Resonances");
1001  fhMCPi0PtOrigin->GetYaxis()->SetBinLabel(4 ,"Resonances");
1002  fhMCPi0PtOrigin->GetYaxis()->SetBinLabel(5 ,"#rho");
1003  fhMCPi0PtOrigin->GetYaxis()->SetBinLabel(6 ,"#omega");
1004  fhMCPi0PtOrigin->GetYaxis()->SetBinLabel(7 ,"K");
1005  fhMCPi0PtOrigin->GetYaxis()->SetBinLabel(8 ,"Other");
1006  fhMCPi0PtOrigin->GetYaxis()->SetBinLabel(9 ,"#eta");
1007  fhMCPi0PtOrigin->GetYaxis()->SetBinLabel(10 ,"#eta prime");
1008  outputContainer->Add(fhMCPi0PtOrigin) ;
1009 
1010  fhMCNotResonancePi0PtOrigin = new TH2F("hMCNotResonancePi0PtOrigin","Reconstructed pair from generated #pi^{0} #it{p}_{T} vs origin",nptbins,ptmin,ptmax,11,0,11) ;
1011  fhMCNotResonancePi0PtOrigin->SetXTitle("#it{p}_{T} (GeV/#it{c})");
1012  fhMCNotResonancePi0PtOrigin->SetYTitle("Origin");
1013  fhMCNotResonancePi0PtOrigin->GetYaxis()->SetBinLabel(1 ,"Status 21");
1014  fhMCNotResonancePi0PtOrigin->GetYaxis()->SetBinLabel(2 ,"Quark");
1015  fhMCNotResonancePi0PtOrigin->GetYaxis()->SetBinLabel(3 ,"qq Resonances");
1016  fhMCNotResonancePi0PtOrigin->GetYaxis()->SetBinLabel(4 ,"Resonances");
1017  fhMCNotResonancePi0PtOrigin->GetYaxis()->SetBinLabel(5 ,"#rho");
1018  fhMCNotResonancePi0PtOrigin->GetYaxis()->SetBinLabel(6 ,"#omega");
1019  fhMCNotResonancePi0PtOrigin->GetYaxis()->SetBinLabel(7 ,"K");
1020  fhMCNotResonancePi0PtOrigin->GetYaxis()->SetBinLabel(8 ,"Other");
1021  fhMCNotResonancePi0PtOrigin->GetYaxis()->SetBinLabel(9 ,"#eta");
1022  fhMCNotResonancePi0PtOrigin->GetYaxis()->SetBinLabel(10 ,"#eta prime");
1023  outputContainer->Add(fhMCNotResonancePi0PtOrigin) ;
1024 
1025  fhMCPi0PtStatus = new TH2F("hMCPi0PtStatus","Reconstructed pair from generated #pi^{0} #it{p}_{T} vs status",nptbins,ptmin,ptmax,101,-50,50) ;
1026  fhMCPi0PtStatus->SetXTitle("#it{p}_{T} (GeV/#it{c})");
1027  fhMCPi0PtStatus->SetYTitle("Status");
1028  outputContainer->Add(fhMCPi0PtStatus) ;
1029 
1030  // Eta
1031 
1032  fhMCEtaPtOrigin = new TH2F("hMCEtaPtOrigin","Reconstructed pair from generated #pi^{0} #it{p}_{T} vs origin",nptbins,ptmin,ptmax,7,0,7) ;
1033  fhMCEtaPtOrigin->SetXTitle("#it{p}_{T} (GeV/#it{c})");
1034  fhMCEtaPtOrigin->SetYTitle("Origin");
1035  fhMCEtaPtOrigin->GetYaxis()->SetBinLabel(1 ,"Status 21");
1036  fhMCEtaPtOrigin->GetYaxis()->SetBinLabel(2 ,"Quark");
1037  fhMCEtaPtOrigin->GetYaxis()->SetBinLabel(3 ,"qq Resonances");
1038  fhMCEtaPtOrigin->GetYaxis()->SetBinLabel(4 ,"Resonances");
1039  fhMCEtaPtOrigin->GetYaxis()->SetBinLabel(5 ,"Other");
1040  fhMCEtaPtOrigin->GetYaxis()->SetBinLabel(6 ,"#eta prime");
1041  outputContainer->Add(fhMCEtaPtOrigin) ;
1042 
1043  fhMCPi0ProdVertex = new TH2F("hMCPi0ProdVertex","Selected reco pair from generated #pi^{0} #it{p}_{T} vs production vertex",
1044  200,0.,20.,5000,0,500) ;
1045  fhMCPi0ProdVertex->SetXTitle("#it{p}_{T} (GeV/#it{c})");
1046  fhMCPi0ProdVertex->SetYTitle("#it{R} (cm)");
1047  outputContainer->Add(fhMCPi0ProdVertex) ;
1048 
1049  fhMCEtaProdVertex = new TH2F("hMCEtaProdVertex","Selected reco pair from generated #eta #it{p}_{T} vs production vertex",
1050  200,0.,20.,5000,0,500) ;
1051  fhMCEtaProdVertex->SetXTitle("#it{p}_{T} (GeV/#it{c})");
1052  fhMCEtaProdVertex->SetYTitle("#it{R} (cm)");
1053  outputContainer->Add(fhMCEtaProdVertex) ;
1054 
1055  for(Int_t i = 0; i<13; i++)
1056  {
1057  fhMCOrgMass[i] = new TH2F(Form("hMCOrgMass_%d",i),Form("#it{M} vs #it{p}_{T}, origin %d",i),nptbins,ptmin,ptmax,nmassbins,massmin,massmax) ;
1058  fhMCOrgMass[i]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
1059  fhMCOrgMass[i]->SetYTitle("#it{M}_{#gamma,#gamma} (GeV/#it{c}^{2})");
1060  outputContainer->Add(fhMCOrgMass[i]) ;
1061 
1062  fhMCOrgAsym[i]= new TH2F(Form("hMCOrgAsym_%d",i),Form("#it{Asymmetry} vs #it{p}_{T}, origin %d",i),nptbins,ptmin,ptmax,nasymbins,asymmin,asymmax) ;
1063  fhMCOrgAsym[i]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
1064  fhMCOrgAsym[i]->SetYTitle("A");
1065  outputContainer->Add(fhMCOrgAsym[i]) ;
1066 
1067  fhMCOrgDeltaEta[i] = new TH2F(Form("hMCOrgDeltaEta_%d",i),Form("#Delta #eta of pair vs #it{p}_{T}, origin %d",i),nptbins,ptmin,ptmax,netabins,-1.4,1.4) ;
1068  fhMCOrgDeltaEta[i]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
1069  fhMCOrgDeltaEta[i]->SetYTitle("#Delta #eta");
1070  outputContainer->Add(fhMCOrgDeltaEta[i]) ;
1071 
1072  fhMCOrgDeltaPhi[i]= new TH2F(Form("hMCOrgDeltaPhi_%d",i),Form("#Delta #phi of pair vs #it{p}_{T}, origin %d",i),nptbins,ptmin,ptmax,nphibins,-0.7,0.7) ;
1073  fhMCOrgDeltaPhi[i]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
1074  fhMCOrgDeltaPhi[i]->SetYTitle("#Delta #phi (rad)");
1075  outputContainer->Add(fhMCOrgDeltaPhi[i]) ;
1076  }
1077 
1078  if(fMultiCutAnaSim)
1079  {
1080  fhMCPi0MassPtTrue = new TH2F*[fNPtCuts*fNAsymCuts*fNCellNCuts];
1081  fhMCPi0MassPtRec = new TH2F*[fNPtCuts*fNAsymCuts*fNCellNCuts];
1082  fhMCPi0PtTruePtRec = new TH2F*[fNPtCuts*fNAsymCuts*fNCellNCuts];
1083  fhMCEtaMassPtRec = new TH2F*[fNPtCuts*fNAsymCuts*fNCellNCuts];
1084  fhMCEtaMassPtTrue = new TH2F*[fNPtCuts*fNAsymCuts*fNCellNCuts];
1085  fhMCEtaPtTruePtRec = new TH2F*[fNPtCuts*fNAsymCuts*fNCellNCuts];
1086 
1087  for(Int_t ipt=0; ipt<fNPtCuts; ipt++)
1088  {
1089  for(Int_t icell=0; icell<fNCellNCuts; icell++)
1090  {
1091  for(Int_t iasym=0; iasym<fNAsymCuts; iasym++)
1092  {
1093  Int_t index = ((ipt*fNCellNCuts)+icell)*fNAsymCuts + iasym;
1094 
1095  fhMCPi0MassPtRec[index] = new TH2F(Form("hMCPi0MassPtRec_pt%d_cell%d_asym%d",ipt,icell,iasym),
1096  Form("Reconstructed #it{M} vs reconstructed #it{p}_{T} of true #pi^{0} cluster pairs for #it{p}_{T} >%2.2f, #it{N}^{cluster}_{cell}>%d and |#it{A}|>%1.2f",fPtCuts[ipt],fCellNCuts[icell], fAsymCuts[iasym]),
1097  nptbins,ptmin,ptmax,nmassbins,massmin,massmax) ;
1098  fhMCPi0MassPtRec[index]->SetXTitle("#it{p}_{T, reconstructed} (GeV/#it{c})");
1099  fhMCPi0MassPtRec[index]->SetYTitle("#it{M}_{#gamma,#gamma} (GeV/#it{c}^{2})");
1100  outputContainer->Add(fhMCPi0MassPtRec[index]) ;
1101 
1102  fhMCPi0MassPtTrue[index] = new TH2F(Form("hMCPi0MassPtTrue_pt%d_cell%d_asym%d",ipt,icell,iasym),
1103  Form("Reconstructed #it{M} vs generated #it{p}_{T} of true #pi^{0} cluster pairs for #it{p}_{T} >%2.2f, #it{N}^{cluster}_{cell}>%d and |#it{A}|>%1.2f",fPtCuts[ipt],fCellNCuts[icell], fAsymCuts[iasym]),
1104  nptbins,ptmin,ptmax,nmassbins,massmin,massmax) ;
1105  fhMCPi0MassPtTrue[index]->SetXTitle("#it{p}_{T, generated} (GeV/#it{c})");
1106  fhMCPi0MassPtTrue[index]->SetYTitle("#it{M}_{#gamma,#gamma} (GeV/#it{c}^{2})");
1107  outputContainer->Add(fhMCPi0MassPtTrue[index]) ;
1108 
1109  fhMCPi0PtTruePtRec[index] = new TH2F(Form("hMCPi0PtTruePtRec_pt%d_cell%d_asym%d",ipt,icell,iasym),
1110  Form("Generated vs reconstructed #it{p}_{T} of true #pi^{0} cluster pairs, 0.01 < rec. mass < 0.17 MeV/#it{c}^{2} for #it{p}_{T} >%2.2f, #it{N}^{cluster}_{cell}>%d and |#it{A}|>%1.2f",fPtCuts[ipt],fCellNCuts[icell], fAsymCuts[iasym]),
1111  nptbins,ptmin,ptmax,nptbins,ptmin,ptmax) ;
1112  fhMCPi0PtTruePtRec[index]->SetXTitle("#it{p}_{T, generated} (GeV/#it{c})");
1113  fhMCPi0PtTruePtRec[index]->SetYTitle("#it{p}_{T, reconstructed} (GeV/#it{c})");
1114  outputContainer->Add(fhMCPi0PtTruePtRec[index]) ;
1115 
1116  fhMCEtaMassPtRec[index] = new TH2F(Form("hMCEtaMassPtRec_pt%d_cell%d_asym%d",ipt,icell,iasym),
1117  Form("Reconstructed #it{M} vs reconstructed #it{p}_{T} of true #eta cluster pairs for #it{p}_{T} >%2.2f, #it{N}^{cluster}_{cell}>%d and |#it{A}|>%1.2f",fPtCuts[ipt],fCellNCuts[icell], fAsymCuts[iasym]),
1118  nptbins,ptmin,ptmax,nmassbins,massmin,massmax) ;
1119  fhMCEtaMassPtRec[index]->SetXTitle("#it{p}_{T, generated} (GeV/#it{c})");
1120  fhMCEtaMassPtRec[index]->SetYTitle("#it{M}_{#gamma,#gamma} (GeV/#it{c}^{2})");
1121  outputContainer->Add(fhMCEtaMassPtRec[index]) ;
1122 
1123  fhMCEtaMassPtTrue[index] = new TH2F(Form("hMCEtaMassPtTrue_pt%d_cell%d_asym%d",ipt,icell,iasym),
1124  Form("Reconstructed #it{M} vs generated #it{p}_{T} of true #eta cluster pairs for #it{p}_{T} >%2.2f, #it{N}^{cluster}_{cell}>%d and |#it{A}|>%1.2f",fPtCuts[ipt],fCellNCuts[icell], fAsymCuts[iasym]),
1125  nptbins,ptmin,ptmax,nmassbins,massmin,massmax) ;
1126  fhMCEtaMassPtTrue[index]->SetXTitle("#it{p}_{T, generated} (GeV/#it{c})");
1127  fhMCEtaMassPtTrue[index]->SetYTitle("#it{M}_{#gamma,#gamma} (GeV/#it{c}^{2})");
1128  outputContainer->Add(fhMCEtaMassPtTrue[index]) ;
1129 
1130  fhMCEtaPtTruePtRec[index] = new TH2F(Form("hMCEtaPtTruePtRec_pt%d_cell%d_asym%d",ipt,icell,iasym),
1131  Form("Generated vs reconstructed #it{p}_{T} of true #eta cluster pairs, 0.01 < rec. mass < 0.17 MeV/#it{c}^{2} for #it{p}_{T} >%2.2f, ncell>%d and asym >%1.2f",fPtCuts[ipt],fCellNCuts[icell], fAsymCuts[iasym]),
1132  nptbins,ptmin,ptmax,nptbins,ptmin,ptmax) ;
1133  fhMCEtaPtTruePtRec[index]->SetXTitle("#it{p}_{T, generated} (GeV/#it{c})");
1134  fhMCEtaPtTruePtRec[index]->SetYTitle("#it{p}_{T, reconstructed} (GeV/#it{c})");
1135  outputContainer->Add(fhMCEtaPtTruePtRec[index]) ;
1136  }
1137  }
1138  }
1139  }//multi cut ana
1140  else
1141  {
1142  fhMCPi0MassPtTrue = new TH2F*[1];
1143  fhMCPi0PtTruePtRec = new TH2F*[1];
1144  fhMCEtaMassPtTrue = new TH2F*[1];
1145  fhMCEtaPtTruePtRec = new TH2F*[1];
1146 
1147  fhMCPi0MassPtTrue[0] = new TH2F("hMCPi0MassPtTrue","Reconstructed Mass vs generated p_T of true #pi^{0} cluster pairs",nptbins,ptmin,ptmax,nmassbins,massmin,massmax) ;
1148  fhMCPi0MassPtTrue[0]->SetXTitle("#it{p}_{T, generated} (GeV/#it{c})");
1149  fhMCPi0MassPtTrue[0]->SetYTitle("#it{M}_{#gamma,#gamma} (GeV/#it{c}^{2})");
1150  outputContainer->Add(fhMCPi0MassPtTrue[0]) ;
1151 
1152  fhMCPi0PtTruePtRec[0]= new TH2F("hMCPi0PtTruePtRec","Generated vs reconstructed p_T of true #pi^{0} cluster pairs, 0.01 < rec. mass < 0.17 MeV/#it{c}^{2}",nptbins,ptmin,ptmax,nptbins,ptmin,ptmax) ;
1153  fhMCPi0PtTruePtRec[0]->SetXTitle("#it{p}_{T, generated} (GeV/#it{c})");
1154  fhMCPi0PtTruePtRec[0]->SetYTitle("#it{p}_{T, reconstructed} (GeV/#it{c})");
1155  outputContainer->Add(fhMCPi0PtTruePtRec[0]) ;
1156 
1157  fhMCEtaMassPtTrue[0] = new TH2F("hMCEtaMassPtTrue","Reconstructed Mass vs generated p_T of true #eta cluster pairs",nptbins,ptmin,ptmax,nmassbins,massmin,massmax) ;
1158  fhMCEtaMassPtTrue[0]->SetXTitle("#it{p}_{T, generated} (GeV/#it{c})");
1159  fhMCEtaMassPtTrue[0]->SetYTitle("#it{M}_{#gamma,#gamma} (GeV/#it{c}^{2})");
1160  outputContainer->Add(fhMCEtaMassPtTrue[0]) ;
1161 
1162  fhMCEtaPtTruePtRec[0]= new TH2F("hMCEtaPtTruePtRec","Generated vs reconstructed p_T of true #eta cluster pairs, 0.01 < rec. mass < 0.17 MeV/#it{c}^{2}",nptbins,ptmin,ptmax,nptbins,ptmin,ptmax) ;
1163  fhMCEtaPtTruePtRec[0]->SetXTitle("#it{p}_{T, generated} (GeV/#it{c})");
1164  fhMCEtaPtTruePtRec[0]->SetYTitle("#it{p}_{T, reconstructed} (GeV/#it{c})");
1165  outputContainer->Add(fhMCEtaPtTruePtRec[0]) ;
1166  }
1167  }
1168  }
1169 
1170  if(fFillSMCombinations)
1171  {
1172  TString pairnamePHOS[] = {"(0-1)","(0-2)","(1-2)","(0-3)","(0-4)","(1-3)","(1-4)","(2-3)","(2-4)","(3-4)"};
1173  for(Int_t imod=0; imod<fNModules; imod++)
1174  {
1175  //Module dependent invariant mass
1176  snprintf(key, buffersize,"hReMod_%d",imod) ;
1177  snprintf(title, buffersize,"Real #it{M}_{#gamma#gamma} distr. for Module %d",imod) ;
1178  fhReMod[imod] = new TH2F(key,title,nptbins,ptmin,ptmax,nmassbins,massmin,massmax) ;
1179  fhReMod[imod]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
1180  fhReMod[imod]->SetYTitle("#it{M}_{#gamma,#gamma} (GeV/#it{c}^{2})");
1181  outputContainer->Add(fhReMod[imod]) ;
1182  if(GetCalorimeter()==kPHOS)
1183  {
1184  snprintf(key, buffersize,"hReDiffPHOSMod_%d",imod) ;
1185  snprintf(title, buffersize,"Real pairs PHOS, clusters in different Modules: %s",(pairnamePHOS[imod]).Data()) ;
1186  fhReDiffPHOSMod[imod] = new TH2F(key,title,nptbins,ptmin,ptmax,nmassbins,massmin,massmax) ;
1187  fhReDiffPHOSMod[imod]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
1188  fhReDiffPHOSMod[imod]->SetYTitle("#it{M}_{#gamma,#gamma} (GeV/#it{c}^{2})");
1189  outputContainer->Add(fhReDiffPHOSMod[imod]) ;
1190  }
1191  else
1192  {//EMCAL
1193  if(imod<fNModules/2)
1194  {
1195  snprintf(key, buffersize,"hReSameSectorEMCAL_%d",imod) ;
1196  snprintf(title, buffersize,"Real pairs EMCAL, clusters in same sector, SM(%d,%d)",imod*2,imod*2+1) ;
1197  fhReSameSectorEMCALMod[imod] = new TH2F(key,title,nptbins,ptmin,ptmax,nmassbins,massmin,massmax) ;
1198  fhReSameSectorEMCALMod[imod]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
1199  fhReSameSectorEMCALMod[imod]->SetYTitle("#it{M}_{#gamma,#gamma} (GeV/#it{c}^{2})");
1200  outputContainer->Add(fhReSameSectorEMCALMod[imod]) ;
1201  }
1202  if(imod<fNModules-2)
1203  {
1204  snprintf(key, buffersize,"hReSameSideEMCAL_%d",imod) ;
1205  snprintf(title, buffersize,"Real pairs EMCAL, clusters in same side SM(%d,%d)",imod, imod+2) ;
1206  fhReSameSideEMCALMod[imod] = new TH2F(key,title,nptbins,ptmin,ptmax,nmassbins,massmin,massmax) ;
1207  fhReSameSideEMCALMod[imod]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
1208  fhReSameSideEMCALMod[imod]->SetYTitle("#it{M}_{#gamma,#gamma} (GeV/#it{c}^{2})");
1209  outputContainer->Add(fhReSameSideEMCALMod[imod]) ;
1210  }
1211  }//EMCAL
1212 
1213  if(DoOwnMix())
1214  {
1215  snprintf(key, buffersize,"hMiMod_%d",imod) ;
1216  snprintf(title, buffersize,"Mixed #it{M}_{#gamma#gamma} distr. for Module %d",imod) ;
1217  fhMiMod[imod] = new TH2F(key,title,nptbins,ptmin,ptmax,nmassbins,massmin,massmax) ;
1218  fhMiMod[imod]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
1219  fhMiMod[imod]->SetYTitle("#it{M}_{#gamma,#gamma} (GeV/#it{c}^{2})");
1220  outputContainer->Add(fhMiMod[imod]) ;
1221 
1222  if(GetCalorimeter()==kPHOS)
1223  {
1224  snprintf(key, buffersize,"hMiDiffPHOSMod_%d",imod) ;
1225  snprintf(title, buffersize,"Mixed pairs PHOS, clusters in different Modules: %s",(pairnamePHOS[imod]).Data()) ;
1226  fhMiDiffPHOSMod[imod] = new TH2F(key,title,nptbins,ptmin,ptmax,nmassbins,massmin,massmax) ;
1227  fhMiDiffPHOSMod[imod]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
1228  fhMiDiffPHOSMod[imod]->SetYTitle("#it{M}_{#gamma,#gamma} (GeV/#it{c}^{2})");
1229  outputContainer->Add(fhMiDiffPHOSMod[imod]) ;
1230  }//PHOS
1231  else
1232  {//EMCAL
1233  if(imod<fNModules/2)
1234  {
1235  snprintf(key, buffersize,"hMiSameSectorEMCALMod_%d",imod) ;
1236  snprintf(title, buffersize,"Mixed pairs EMCAL, clusters in same sector, SM(%d,%d)",imod*2,imod*2+1) ;
1237  fhMiSameSectorEMCALMod[imod] = new TH2F(key,title,nptbins,ptmin,ptmax,nmassbins,massmin,massmax) ;
1238  fhMiSameSectorEMCALMod[imod]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
1239  fhMiSameSectorEMCALMod[imod]->SetYTitle("#it{M}_{#gamma,#gamma} (GeV/#it{c}^{2})");
1240  outputContainer->Add(fhMiSameSectorEMCALMod[imod]) ;
1241  }
1242  if(imod<fNModules-2){
1243 
1244  snprintf(key, buffersize,"hMiSameSideEMCALMod_%d",imod) ;
1245  snprintf(title, buffersize,"Mixed pairs EMCAL, clusters in same side SM(%d,%d)",imod, imod+2) ;
1246  fhMiSameSideEMCALMod[imod] = new TH2F(key,title,nptbins,ptmin,ptmax,nmassbins,massmin,massmax) ;
1247  fhMiSameSideEMCALMod[imod]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
1248  fhMiSameSideEMCALMod[imod]->SetYTitle("#it{M}_{#gamma,#gamma} (GeV/#it{c}^{2})");
1249  outputContainer->Add(fhMiSameSideEMCALMod[imod]) ;
1250  }
1251  } // EMCAL
1252  } // own mix
1253  } // loop combinations
1254  } // SM combinations
1255 
1256 // for(Int_t i = 0; i < outputContainer->GetEntries() ; i++){
1257 //
1258 // printf("Histogram %d, name: %s\n ",i, outputContainer->At(i)->GetName());
1259 //
1260 // }
1261 
1262  return outputContainer;
1263 }
1264 
1265 //___________________________________________________
1267 //___________________________________________________
1268 void AliAnaPi0::Print(const Option_t * /*opt*/) const
1269 {
1270  printf("**** Print %s %s ****\n", GetName(), GetTitle() ) ;
1271  AliAnaCaloTrackCorrBaseClass::Print(" ");
1272 
1273  printf("Number of bins in Centrality: %d \n",GetNCentrBin()) ;
1274  printf("Number of bins in Z vert. pos: %d \n",GetNZvertBin()) ;
1275  printf("Number of bins in Reac. Plain: %d \n",GetNRPBin()) ;
1276  printf("Depth of event buffer: %d \n",GetNMaxEvMix()) ;
1277  printf("Pair in same Module: %d \n",fSameSM) ;
1278  printf("Cuts: \n") ;
1279  // printf("Z vertex position: -%2.3f < z < %2.3f \n",GetZvertexCut(),GetZvertexCut()) ; //It crashes here, why?
1280  printf("Number of modules: %d \n",fNModules) ;
1281  printf("Select pairs with their angle: %d, edep %d, min angle %2.3f, max angle %2.3f \n",fUseAngleCut, fUseAngleEDepCut, fAngleCut, fAngleMaxCut) ;
1282  printf("Asymmetry cuts: n = %d, \n",fNAsymCuts) ;
1283  printf("\tasymmetry < ");
1284  for(Int_t i = 0; i < fNAsymCuts; i++) printf("%2.2f ",fAsymCuts[i]);
1285  printf("\n");
1286 
1287  printf("PID selection bits: n = %d, \n",fNPIDBits) ;
1288  printf("\tPID bit = ");
1289  for(Int_t i = 0; i < fNPIDBits; i++) printf("%d ",fPIDBits[i]);
1290  printf("\n");
1291 
1292  if(fMultiCutAna)
1293  {
1294  printf("pT cuts: n = %d, \n",fNPtCuts) ;
1295  printf("\tpT > ");
1296  for(Int_t i = 0; i < fNPtCuts; i++) printf("%2.2f ",fPtCuts[i]);
1297  printf("GeV/c\n");
1298 
1299  printf("N cell in cluster cuts: n = %d, \n",fNCellNCuts) ;
1300  printf("\tnCell > ");
1301  for(Int_t i = 0; i < fNCellNCuts; i++) printf("%d ",fCellNCuts[i]);
1302  printf("\n");
1303 
1304  }
1305  printf("------------------------------------------------------\n") ;
1306 }
1307 
1308 //________________________________________
1310 //________________________________________
1311 void AliAnaPi0::FillAcceptanceHistograms()
1312 {
1313  Double_t mesonY = -100 ;
1314  Double_t mesonE = -1 ;
1315  Double_t mesonPt = -1 ;
1316  Double_t mesonPhi = 100 ;
1317  Double_t mesonYeta= -1 ;
1318 
1319  Int_t pdg = 0 ;
1320  Int_t nprim = 0 ;
1321  Int_t nDaught = 0 ;
1322  Int_t iphot1 = 0 ;
1323  Int_t iphot2 = 0 ;
1324 
1325  Float_t cen = GetEventCentrality();
1326  Float_t ep = GetEventPlaneAngle();
1327 
1328  TParticle * primStack = 0;
1329  AliAODMCParticle * primAOD = 0;
1330 
1331  // Get the ESD MC particles container
1332  AliStack * stack = 0;
1333  if( GetReader()->ReadStack() )
1334  {
1335  stack = GetMCStack();
1336  if(!stack ) return;
1337  nprim = stack->GetNtrack();
1338  }
1339 
1340  // Get the AOD MC particles container
1341  TClonesArray * mcparticles = 0;
1342  if( GetReader()->ReadAODMCParticles() )
1343  {
1344  mcparticles = GetReader()->GetAODMCParticles();
1345  if( !mcparticles ) return;
1346  nprim = mcparticles->GetEntriesFast();
1347  }
1348 
1349  for(Int_t i=0 ; i < nprim; i++)
1350  {
1351  if(GetReader()->AcceptOnlyHIJINGLabels() && !GetReader()->IsHIJINGLabel(i)) continue ;
1352 
1353  if(GetReader()->ReadStack())
1354  {
1355  primStack = stack->Particle(i) ;
1356  if(!primStack)
1357  {
1358  AliWarning("ESD primaries pointer not available!!");
1359  continue;
1360  }
1361 
1362  // If too small skip
1363  if( primStack->Energy() < 0.4 ) continue;
1364 
1365  pdg = primStack->GetPdgCode();
1366  nDaught = primStack->GetNDaughters();
1367  iphot1 = primStack->GetDaughter(0) ;
1368  iphot2 = primStack->GetDaughter(1) ;
1369  if(primStack->Energy() == TMath::Abs(primStack->Pz())) continue ; //Protection against floating point exception
1370 
1371  //printf("i %d, %s %d %s %d \n",i, stack->Particle(i)->GetName(), stack->Particle(i)->GetPdgCode(),
1372  // prim->GetName(), prim->GetPdgCode());
1373 
1374  //Photon kinematics
1375  primStack->Momentum(fPi0Mom);
1376 
1377  mesonY = 0.5*TMath::Log((primStack->Energy()+primStack->Pz())/(primStack->Energy()-primStack->Pz())) ;
1378  }
1379  else
1380  {
1381  primAOD = (AliAODMCParticle *) mcparticles->At(i);
1382  if(!primAOD)
1383  {
1384  AliWarning("AOD primaries pointer not available!!");
1385  continue;
1386  }
1387 
1388  // If too small skip
1389  if( primAOD->E() < 0.4 ) continue;
1390 
1391  pdg = primAOD->GetPdgCode();
1392  nDaught = primAOD->GetNDaughters();
1393  iphot1 = primAOD->GetFirstDaughter() ;
1394  iphot2 = primAOD->GetLastDaughter() ;
1395 
1396  if(primAOD->E() == TMath::Abs(primAOD->Pz())) continue ; //Protection against floating point exception
1397 
1398  //Photon kinematics
1399  fPi0Mom.SetPxPyPzE(primAOD->Px(),primAOD->Py(),primAOD->Pz(),primAOD->E());
1400 
1401  mesonY = 0.5*TMath::Log((primAOD->E()+primAOD->Pz())/(primAOD->E()-primAOD->Pz())) ;
1402  }
1403 
1404  // Select only pi0 or eta
1405  if( pdg != 111 && pdg != 221) continue ;
1406 
1407  mesonPt = fPi0Mom.Pt () ;
1408  mesonE = fPi0Mom.E () ;
1409  mesonYeta= fPi0Mom.Eta() ;
1410  mesonPhi = fPi0Mom.Phi() ;
1411  if( mesonPhi < 0 ) mesonPhi+=TMath::TwoPi();
1412  mesonPhi *= TMath::RadToDeg();
1413 
1414  if(pdg == 111)
1415  {
1416  if(TMath::Abs(mesonY) < 1.0)
1417  {
1418  fhPrimPi0E ->Fill(mesonE , GetEventWeight()) ;
1419  fhPrimPi0Pt ->Fill(mesonPt, GetEventWeight()) ;
1420  fhPrimPi0Phi->Fill(mesonPt, mesonPhi, GetEventWeight()) ;
1421 
1422  fhPrimPi0YetaYcut->Fill(mesonPt, mesonYeta, GetEventWeight()) ;
1423 
1424  if( IsHighMultiplicityAnalysisOn() )
1425  {
1426  fhPrimPi0PtCentrality->Fill(mesonPt, cen, GetEventWeight()) ;
1427  fhPrimPi0PtEventPlane->Fill(mesonPt, ep , GetEventWeight()) ;
1428  }
1429  }
1430 
1431  fhPrimPi0Y ->Fill(mesonPt, mesonY , GetEventWeight()) ;
1432  fhPrimPi0Yeta->Fill(mesonPt, mesonYeta, GetEventWeight()) ;
1433  }
1434  else if(pdg == 221)
1435  {
1436  if(TMath::Abs(mesonY) < 1.0)
1437  {
1438  fhPrimEtaE ->Fill(mesonE , GetEventWeight()) ;
1439  fhPrimEtaPt ->Fill(mesonPt, GetEventWeight()) ;
1440  fhPrimEtaPhi->Fill(mesonPt, mesonPhi, GetEventWeight()) ;
1441 
1442  fhPrimEtaYetaYcut->Fill(mesonPt, mesonYeta, GetEventWeight()) ;
1443 
1444  if( IsHighMultiplicityAnalysisOn() )
1445  {
1446  fhPrimEtaPtCentrality->Fill(mesonPt, cen, GetEventWeight()) ;
1447  fhPrimEtaPtEventPlane->Fill(mesonPt, ep , GetEventWeight()) ;
1448  }
1449  }
1450 
1451  fhPrimEtaY ->Fill(mesonPt, mesonY , GetEventWeight()) ;
1452  fhPrimEtaYeta->Fill(mesonPt, mesonYeta, GetEventWeight()) ;
1453  }
1454 
1455  //Origin of meson
1456  if(fFillOriginHisto && TMath::Abs(mesonY) < 0.7)
1457  {
1458  Int_t momindex = -1;
1459  Int_t mompdg = -1;
1460  Int_t momstatus = -1;
1461  Int_t status = -1;
1462  Float_t momR = 0;
1463  if(GetReader()->ReadStack()) momindex = primStack->GetFirstMother();
1464  if(GetReader()->ReadAODMCParticles()) momindex = primAOD ->GetMother();
1465 
1466  if(momindex >= 0 && momindex < nprim)
1467  {
1468  if(GetReader()->ReadStack())
1469  {
1470  status = primStack->GetStatusCode();
1471  TParticle* mother = stack->Particle(momindex);
1472  mompdg = TMath::Abs(mother->GetPdgCode());
1473  momstatus = mother->GetStatusCode();
1474  momR = mother->R();
1475  }
1476 
1477  if(GetReader()->ReadAODMCParticles())
1478  {
1479  status = primAOD->GetStatus();
1480  AliAODMCParticle* mother = (AliAODMCParticle*) mcparticles->At(momindex);
1481  mompdg = TMath::Abs(mother->GetPdgCode());
1482  momstatus = mother->GetStatus();
1483  momR = TMath::Sqrt(mother->Xv()*mother->Xv()+mother->Yv()*mother->Yv());
1484  }
1485 
1486  if(pdg == 111)
1487  {
1488  fhPrimPi0ProdVertex->Fill(mesonPt, momR , GetEventWeight());
1489  fhPrimPi0PtStatus ->Fill(mesonPt, status, GetEventWeight());
1490 
1491 
1492  if (momstatus == 21) fhPrimPi0PtOrigin->Fill(mesonPt, 0.5, GetEventWeight());//parton
1493  else if(mompdg < 22 ) fhPrimPi0PtOrigin->Fill(mesonPt, 1.5, GetEventWeight());//quark
1494  else if(mompdg > 2100 && mompdg < 2210)
1495  fhPrimPi0PtOrigin->Fill(mesonPt, 2.5, GetEventWeight());// resonances
1496  else if(mompdg == 221) fhPrimPi0PtOrigin->Fill(mesonPt, 8.5, GetEventWeight());//eta
1497  else if(mompdg == 331) fhPrimPi0PtOrigin->Fill(mesonPt, 9.5, GetEventWeight());//eta prime
1498  else if(mompdg == 213) fhPrimPi0PtOrigin->Fill(mesonPt, 4.5, GetEventWeight());//rho
1499  else if(mompdg == 223) fhPrimPi0PtOrigin->Fill(mesonPt, 5.5, GetEventWeight());//omega
1500  else if(mompdg >= 310 && mompdg <= 323)
1501  fhPrimPi0PtOrigin->Fill(mesonPt, 6.5, GetEventWeight());//k0S, k+-,k*
1502  else if(mompdg == 130) fhPrimPi0PtOrigin->Fill(mesonPt, 6.5, GetEventWeight());//k0L
1503  else if(momstatus == 11 || momstatus == 12 )
1504  fhPrimPi0PtOrigin->Fill(mesonPt, 3.5, GetEventWeight());//resonances
1505  else fhPrimPi0PtOrigin->Fill(mesonPt, 7.5, GetEventWeight());//other?
1506 
1507  //printf("Prim Pi0: index %d, pt %2.2f Prod vertex %3.3f, origin pdg %d, origin status %d, origin UI %d\n",
1508  // momindex, mesonPt,mother->R(),mompdg,momstatus,mother->GetUniqueID());
1509 
1510  if(status!=11)
1511  {
1512  if (momstatus == 21) fhPrimNotResonancePi0PtOrigin->Fill(mesonPt, 0.5, GetEventWeight());//parton
1513  else if(mompdg < 22 ) fhPrimNotResonancePi0PtOrigin->Fill(mesonPt, 1.5, GetEventWeight());//quark
1514  else if(mompdg > 2100 && mompdg < 2210)
1515  fhPrimNotResonancePi0PtOrigin->Fill(mesonPt, 2.5, GetEventWeight());// resonances
1516  else if(mompdg == 221) fhPrimNotResonancePi0PtOrigin->Fill(mesonPt, 8.5, GetEventWeight());//eta
1517  else if(mompdg == 331) fhPrimNotResonancePi0PtOrigin->Fill(mesonPt, 9.5, GetEventWeight());//eta prime
1518  else if(mompdg == 213) fhPrimNotResonancePi0PtOrigin->Fill(mesonPt, 4.5, GetEventWeight());//rho
1519  else if(mompdg == 223) fhPrimNotResonancePi0PtOrigin->Fill(mesonPt, 5.5, GetEventWeight());//omega
1520  else if(mompdg >= 310 && mompdg <= 323)
1521  fhPrimNotResonancePi0PtOrigin->Fill(mesonPt, 6.5, GetEventWeight());//k0S, k+-,k*
1522  else if(mompdg == 130) fhPrimNotResonancePi0PtOrigin->Fill(mesonPt, 6.5, GetEventWeight());//k0L
1523  else if(momstatus == 11 || momstatus == 12 )
1524  fhPrimNotResonancePi0PtOrigin->Fill(mesonPt, 3.5, GetEventWeight());//resonances
1525  else fhPrimNotResonancePi0PtOrigin->Fill(mesonPt, 7.5, GetEventWeight());//other?
1526  }
1527 
1528  }//pi0
1529  else
1530  {
1531  if (momstatus == 21 ) fhPrimEtaPtOrigin->Fill(mesonPt, 0.5, GetEventWeight());//parton
1532  else if(mompdg < 22 ) fhPrimEtaPtOrigin->Fill(mesonPt, 1.5, GetEventWeight());//quark
1533  else if(mompdg > 2100 && mompdg < 2210)
1534  fhPrimEtaPtOrigin->Fill(mesonPt, 2.5, GetEventWeight());//qq resonances
1535  else if(mompdg == 331) fhPrimEtaPtOrigin->Fill(mesonPt, 5.5, GetEventWeight());//eta prime
1536  else if(momstatus == 11 || momstatus == 12 )
1537  fhPrimEtaPtOrigin->Fill(mesonPt, 3.5, GetEventWeight());//resonances
1538  else fhPrimEtaPtOrigin->Fill(mesonPt, 4.5, GetEventWeight());//stable, conversions?
1539  //printf("Other Meson pdg %d, Mother %s, pdg %d, status %d\n",pdg, TDatabasePDG::Instance()->GetParticle(mompdg)->GetName(),mompdg, momstatus );
1540 
1541  fhPrimEtaProdVertex->Fill(mesonPt, momR, GetEventWeight());
1542  }
1543  } // pi0 has mother
1544  }
1545 
1546  //Check if both photons hit Calorimeter
1547  if(nDaught != 2 ) continue; //Only interested in 2 gamma decay
1548 
1549  if(iphot1 < 0 || iphot1 >= nprim || iphot2 < 0 || iphot2 >= nprim) continue ;
1550 
1551  Int_t pdg1 = 0;
1552  Int_t pdg2 = 0;
1553  Bool_t inacceptance1 = kTRUE;
1554  Bool_t inacceptance2 = kTRUE;
1555 
1556  if(GetReader()->ReadStack())
1557  {
1558  TParticle * phot1 = stack->Particle(iphot1) ;
1559  TParticle * phot2 = stack->Particle(iphot2) ;
1560 
1561  if(!phot1 || !phot2) continue ;
1562 
1563  pdg1 = phot1->GetPdgCode();
1564  pdg2 = phot2->GetPdgCode();
1565 
1566  phot1->Momentum(fPhotonMom1);
1567  phot2->Momentum(fPhotonMom2);
1568 
1569  // Check if photons hit the Calorimeter acceptance
1570  if(IsRealCaloAcceptanceOn())
1571  {
1572  if( !GetCaloUtils()->IsMCParticleInCalorimeterAcceptance( GetCalorimeter(), phot1 )) inacceptance1 = kFALSE ;
1573  if( !GetCaloUtils()->IsMCParticleInCalorimeterAcceptance( GetCalorimeter(), phot2 )) inacceptance2 = kFALSE ;
1574  }
1575  }
1576 
1577  if(GetReader()->ReadAODMCParticles())
1578  {
1579  AliAODMCParticle * phot1 = (AliAODMCParticle *) mcparticles->At(iphot1) ;
1580  AliAODMCParticle * phot2 = (AliAODMCParticle *) mcparticles->At(iphot2) ;
1581 
1582  if(!phot1 || !phot2) continue ;
1583 
1584  pdg1 = phot1->GetPdgCode();
1585  pdg2 = phot2->GetPdgCode();
1586 
1587  fPhotonMom1.SetPxPyPzE(phot1->Px(),phot1->Py(),phot1->Pz(),phot1->E());
1588  fPhotonMom2.SetPxPyPzE(phot2->Px(),phot2->Py(),phot2->Pz(),phot2->E());
1589 
1590  // Check if photons hit the Calorimeter acceptance
1591  if(IsRealCaloAcceptanceOn())
1592  {
1593  if( !GetCaloUtils()->IsMCParticleInCalorimeterAcceptance( GetCalorimeter(), phot1 )) inacceptance1 = kFALSE ;
1594  if( !GetCaloUtils()->IsMCParticleInCalorimeterAcceptance( GetCalorimeter(), phot2 )) inacceptance2 = kFALSE ;
1595  }
1596  }
1597 
1598  if( pdg1 != 22 || pdg2 !=22) continue ;
1599 
1600  // Check if photons hit desired acceptance in the fidutial borders fixed in the analysis
1601  if(IsFiducialCutOn())
1602  {
1603  if( inacceptance1 && !GetFiducialCut()->IsInFiducialCut(fPhotonMom1.Eta(), fPhotonMom1.Phi(), GetCalorimeter()) ) inacceptance1 = kFALSE ;
1604  if( inacceptance2 && !GetFiducialCut()->IsInFiducialCut(fPhotonMom2.Eta(), fPhotonMom2.Phi(), GetCalorimeter()) ) inacceptance2 = kFALSE ;
1605  }
1606 
1607  if(fFillArmenterosThetaStar) FillArmenterosThetaStar(pdg);
1608 
1609  if(GetCalorimeter()==kEMCAL && inacceptance1 && inacceptance2 && fCheckAccInSector)
1610  {
1611  Int_t absID1=0;
1612  Int_t absID2=0;
1613 
1614  Float_t photonPhi1 = fPhotonMom1.Phi();
1615  Float_t photonPhi2 = fPhotonMom2.Phi();
1616 
1617  if(photonPhi1 < 0) photonPhi1+=TMath::TwoPi();
1618  if(photonPhi2 < 0) photonPhi2+=TMath::TwoPi();
1619 
1620  GetEMCALGeometry()->GetAbsCellIdFromEtaPhi(fPhotonMom1.Eta(),photonPhi1,absID1);
1621  GetEMCALGeometry()->GetAbsCellIdFromEtaPhi(fPhotonMom2.Eta(),photonPhi2,absID2);
1622 
1623  Int_t sm1 = GetEMCALGeometry()->GetSuperModuleNumber(absID1);
1624  Int_t sm2 = GetEMCALGeometry()->GetSuperModuleNumber(absID2);
1625 
1626  Int_t j=0;
1627  Bool_t sameSector = kFALSE;
1628  for(Int_t isector = 0; isector < fNModules/2; isector++)
1629  {
1630  j=2*isector;
1631  if((sm1==j && sm2==j+1) || (sm1==j+1 && sm2==j)) sameSector = kTRUE;
1632  }
1633 
1634  if(sm1!=sm2 && !sameSector)
1635  {
1636  inacceptance1 = kFALSE;
1637  inacceptance2 = kFALSE;
1638  }
1639  //if(sm1!=sm2)printf("sm1 %d, sm2 %d, same sector %d, in acceptance %d\n",sm1,sm2,sameSector,inacceptance);
1640  // if(GetEMCALGeometry()->Impact(phot1) && GetEMCALGeometry()->Impact(phot2))
1641  // inacceptance = kTRUE;
1642  }
1643 
1644  AliDebug(2,Form("Accepted in %s?: m (%2.2f,%2.2f,%2.2f), p1 (%2.2f,%2.2f,%2.2f), p2 (%2.2f,%2.2f,%2.2f) : in1 %d, in2 %d",
1645  GetCalorimeterString().Data(),
1646  mesonPt,mesonYeta,mesonPhi,
1647  fPhotonMom1.Pt(),fPhotonMom1.Eta(),fPhotonMom1.Phi()*TMath::RadToDeg(),
1648  fPhotonMom2.Pt(),fPhotonMom2.Eta(),fPhotonMom2.Phi()*TMath::RadToDeg(),
1649  inacceptance1, inacceptance2));
1650 
1651  if(inacceptance1 && inacceptance2)
1652  {
1653  Float_t asym = TMath::Abs((fPhotonMom1.E()-fPhotonMom2.E()) / (fPhotonMom1.E()+fPhotonMom2.E()));
1654  Double_t angle = fPhotonMom1.Angle(fPhotonMom2.Vect());
1655 
1656  AliDebug(2,Form("\t ACCEPTED pdg %d: pt %2.2f, phi %2.2f, eta %2.2f",pdg,mesonPt,mesonPhi,mesonYeta));
1657 
1658  if(pdg==111)
1659  {
1660  fhPrimPi0AccE ->Fill(mesonE , GetEventWeight()) ;
1661  fhPrimPi0AccPt ->Fill(mesonPt, GetEventWeight()) ;
1662  fhPrimPi0AccPhi ->Fill(mesonPt, mesonPhi , GetEventWeight()) ;
1663  fhPrimPi0AccY ->Fill(mesonPt, mesonY , GetEventWeight()) ;
1664  fhPrimPi0AccYeta->Fill(mesonPt, mesonYeta, GetEventWeight()) ;
1665 
1666  if( IsHighMultiplicityAnalysisOn() )
1667  {
1668  fhPrimPi0AccPtCentrality->Fill(mesonPt, cen, GetEventWeight()) ;
1669  fhPrimPi0AccPtEventPlane->Fill(mesonPt, ep , GetEventWeight()) ;
1670  }
1671 
1672  if(fFillAngleHisto)
1673  {
1674  fhPrimPi0OpeningAngle ->Fill(mesonPt, angle, GetEventWeight());
1675  if(mesonPt > 5)fhPrimPi0OpeningAngleAsym->Fill(asym, angle, GetEventWeight());
1676  fhPrimPi0CosOpeningAngle ->Fill(mesonPt, TMath::Cos(angle), GetEventWeight());
1677  }
1678  }
1679  else if(pdg==221)
1680  {
1681  fhPrimEtaAccE ->Fill(mesonE , GetEventWeight()) ;
1682  fhPrimEtaAccPt ->Fill(mesonPt, GetEventWeight()) ;
1683  fhPrimEtaAccPhi ->Fill(mesonPt, mesonPhi , GetEventWeight()) ;
1684  fhPrimEtaAccY ->Fill(mesonPt, mesonY , GetEventWeight()) ;
1685  fhPrimEtaAccYeta->Fill(mesonPt, mesonYeta, GetEventWeight()) ;
1686 
1687  if( IsHighMultiplicityAnalysisOn() )
1688  {
1689  fhPrimEtaAccPtCentrality->Fill(mesonPt, cen, GetEventWeight()) ;
1690  fhPrimEtaAccPtEventPlane->Fill(mesonPt, ep , GetEventWeight()) ;
1691  }
1692 
1693  if(fFillAngleHisto)
1694  {
1695  fhPrimEtaOpeningAngle ->Fill(mesonPt, angle, GetEventWeight());
1696  if(mesonPt > 5)fhPrimEtaOpeningAngleAsym->Fill(asym, angle, GetEventWeight());
1697  fhPrimEtaCosOpeningAngle ->Fill(mesonPt, TMath::Cos(angle), GetEventWeight());
1698  }
1699  }
1700  } // Accepted
1701  } // loop on primaries
1702 }
1703 
1704 //________________________________________________
1706 //________________________________________________
1707 void AliAnaPi0::FillArmenterosThetaStar(Int_t pdg)
1708 {
1709  // Get pTArm and AlphaArm
1710  Float_t momentumSquaredMother = fPi0Mom.P()*fPi0Mom.P();
1711  Float_t momentumDaughter1AlongMother = 0.;
1712  Float_t momentumDaughter2AlongMother = 0.;
1713 
1714  if (momentumSquaredMother > 0.)
1715  {
1716  momentumDaughter1AlongMother = (fPhotonMom1.Px()*fPi0Mom.Px() + fPhotonMom1.Py()*fPi0Mom.Py()+ fPhotonMom1.Pz()*fPi0Mom.Pz()) / sqrt(momentumSquaredMother);
1717  momentumDaughter2AlongMother = (fPhotonMom2.Px()*fPi0Mom.Px() + fPhotonMom2.Py()*fPi0Mom.Py()+ fPhotonMom2.Pz()*fPi0Mom.Pz()) / sqrt(momentumSquaredMother);
1718  }
1719 
1720  Float_t momentumSquaredDaughter1 = fPhotonMom1.P()*fPhotonMom1.P();
1721  Float_t ptArmSquared = momentumSquaredDaughter1 - momentumDaughter1AlongMother*momentumDaughter1AlongMother;
1722 
1723  Float_t pTArm = 0.;
1724  if (ptArmSquared > 0.)
1725  pTArm = sqrt(ptArmSquared);
1726 
1727  Float_t alphaArm = 0.;
1728  if(momentumDaughter1AlongMother + momentumDaughter2AlongMother > 0)
1729  alphaArm = (momentumDaughter1AlongMother -momentumDaughter2AlongMother) / (momentumDaughter1AlongMother + momentumDaughter2AlongMother);
1730 
1731  fPhotonMom1Boost = fPhotonMom1;
1732  fPhotonMom1Boost.Boost(-fPi0Mom.BoostVector());
1733  Float_t cosThStar=TMath::Cos(fPhotonMom1Boost.Vect().Angle(fPi0Mom.Vect()));
1734 
1735  Float_t en = fPi0Mom.Energy();
1736  Int_t ebin = -1;
1737  if(en > 8 && en <= 12) ebin = 0;
1738  if(en > 12 && en <= 16) ebin = 1;
1739  if(en > 16 && en <= 20) ebin = 2;
1740  if(en > 20) ebin = 3;
1741  if(ebin < 0 || ebin > 3) return ;
1742 
1743  if(pdg==111)
1744  {
1745  fhCosThStarPrimPi0->Fill(en, cosThStar, GetEventWeight());
1746  fhArmPrimPi0[ebin]->Fill(alphaArm, pTArm , GetEventWeight());
1747  }
1748  else
1749  {
1750  fhCosThStarPrimEta->Fill(en, cosThStar, GetEventWeight());
1751  fhArmPrimEta[ebin]->Fill(alphaArm,pTArm , GetEventWeight());
1752  }
1753 
1754  // if(GetDebug() > 2 )
1755  // {
1756  // Float_t asym = 0;
1757  // if(fPhotonMom1.E()+fPhotonMom2.E() > 0) asym = TMath::Abs(fPhotonMom1.E()-fPhotonMom2.E())/(fPhotonMom1.E()+fPhotonMom2.E());
1758  //
1759  // printf("AliAnaPi0::FillArmenterosThetaStar() - E %f, alphaArm %f, pTArm %f, cos(theta*) %f, asymmetry %1.3f\n",
1760  // en,alphaArm,pTArm,cosThStar,asym);
1761  // }
1762 }
1763 
1764 //_______________________________________________________________________________________
1769 //_______________________________________________________________________________________
1770 void AliAnaPi0::FillMCVersusRecDataHistograms(Int_t index1, Int_t index2,
1771  Float_t pt1, Float_t pt2,
1772  Int_t ncell1, Int_t ncell2,
1773  Double_t mass, Double_t pt, Double_t asym,
1774  Double_t deta, Double_t dphi)
1775 {
1776  Int_t ancPDG = 0;
1777  Int_t ancStatus = 0;
1778  Int_t ancLabel = GetMCAnalysisUtils()->CheckCommonAncestor(index1, index2,
1779  GetReader(), ancPDG, ancStatus,fPi0Mom, fProdVertex);
1780 
1781  Int_t momindex = -1;
1782  Int_t mompdg = -1;
1783  Int_t momstatus = -1;
1784  Int_t status = -1;
1785  Float_t prodR = -1;
1786  Int_t mcIndex = -1;
1787 
1788  if(ancLabel > -1)
1789  {
1790  AliDebug(1,Form("Common ancestor label %d, pdg %d, name %s, status %d",
1791  ancLabel,ancPDG,TDatabasePDG::Instance()->GetParticle(ancPDG)->GetName(),ancStatus));
1792 
1793  if(ancPDG==22)
1794  {//gamma
1795  mcIndex = 0;
1796  }
1797  else if(TMath::Abs(ancPDG)==11)
1798  {//e
1799  mcIndex = 1;
1800  }
1801  else if(ancPDG==111)
1802  {//Pi0
1803  mcIndex = 2;
1804  if(fMultiCutAnaSim)
1805  {
1806  for(Int_t ipt=0; ipt<fNPtCuts; ipt++)
1807  {
1808  for(Int_t icell=0; icell<fNCellNCuts; icell++)
1809  {
1810  for(Int_t iasym=0; iasym<fNAsymCuts; iasym++)
1811  {
1812  Int_t index = ((ipt*fNCellNCuts)+icell)*fNAsymCuts + iasym;
1813  if(pt1 > fPtCuts[ipt] && pt2 > fPtCuts[ipt] &&
1814  asym < fAsymCuts[iasym] &&
1815  ncell1 >= fCellNCuts[icell] && ncell2 >= fCellNCuts[icell])
1816  {
1817  fhMCPi0MassPtRec [index]->Fill(pt, mass, GetEventWeight());
1818  fhMCPi0MassPtTrue[index]->Fill(fPi0Mom.Pt(),mass, GetEventWeight());
1819  if(mass < 0.17 && mass > 0.1)
1820  fhMCPi0PtTruePtRec[index]->Fill(fPi0Mom.Pt(), pt, GetEventWeight());
1821  }//pass the different cuts
1822  }// pid bit cut loop
1823  }// icell loop
1824  }// pt cut loop
1825  }// Multi cut ana sim
1826  else
1827  {
1828  fhMCPi0MassPtTrue[0]->Fill(fPi0Mom.Pt(), mass, GetEventWeight());
1829 
1830  if(mass < 0.17 && mass > 0.1)
1831  {
1832  fhMCPi0PtTruePtRec[0]->Fill(fPi0Mom.Pt(), pt, GetEventWeight());
1833 
1834  //Int_t uniqueId = -1;
1835  if(GetReader()->ReadStack())
1836  {
1837  TParticle* ancestor = GetMCStack()->Particle(ancLabel);
1838  status = ancestor->GetStatusCode();
1839  momindex = ancestor->GetFirstMother();
1840  if(momindex < 0) return;
1841  TParticle* mother = GetMCStack()->Particle(momindex);
1842  mompdg = TMath::Abs(mother->GetPdgCode());
1843  momstatus = mother->GetStatusCode();
1844  prodR = mother->R();
1845  //uniqueId = mother->GetUniqueID();
1846  }
1847  else
1848  {
1849  TClonesArray * mcparticles = GetReader()->GetAODMCParticles();
1850  AliAODMCParticle* ancestor = (AliAODMCParticle *) mcparticles->At(ancLabel);
1851  status = ancestor->GetStatus();
1852  momindex = ancestor->GetMother();
1853  if(momindex < 0) return;
1854  AliAODMCParticle* mother = (AliAODMCParticle *) mcparticles->At(momindex);
1855  mompdg = TMath::Abs(mother->GetPdgCode());
1856  momstatus = mother->GetStatus();
1857  prodR = TMath::Sqrt(mother->Xv()*mother->Xv()+mother->Yv()*mother->Yv());
1858  //uniqueId = mother->GetUniqueID();
1859  }
1860 
1861  // printf("Reco Pi0: pt %2.2f Prod vertex %3.3f, origin pdg %d, origin status %d, origin UI %d\n",
1862  // pt,prodR,mompdg,momstatus,uniqueId);
1863 
1864  fhMCPi0ProdVertex->Fill(pt, prodR , GetEventWeight());
1865  fhMCPi0PtStatus ->Fill(pt, status, GetEventWeight());
1866 
1867  if (momstatus == 21) fhMCPi0PtOrigin->Fill(pt, 0.5, GetEventWeight());//parton
1868  else if(mompdg < 22 ) fhMCPi0PtOrigin->Fill(pt, 1.5, GetEventWeight());//quark
1869  else if(mompdg > 2100 && mompdg < 2210)
1870  fhMCPi0PtOrigin->Fill(pt, 2.5, GetEventWeight());// resonances
1871  else if(mompdg == 221) fhMCPi0PtOrigin->Fill(pt, 8.5, GetEventWeight());//eta
1872  else if(mompdg == 331) fhMCPi0PtOrigin->Fill(pt, 9.5, GetEventWeight());//eta prime
1873  else if(mompdg == 213) fhMCPi0PtOrigin->Fill(pt, 4.5, GetEventWeight());//rho
1874  else if(mompdg == 223) fhMCPi0PtOrigin->Fill(pt, 5.5, GetEventWeight());//omega
1875  else if(mompdg >= 310 && mompdg <= 323)
1876  fhMCPi0PtOrigin->Fill(pt, 6.5, GetEventWeight());//k0S, k+-,k*
1877  else if(mompdg == 130) fhMCPi0PtOrigin->Fill(pt, 6.5, GetEventWeight());//k0L
1878  else if(momstatus == 11 || momstatus == 12 )
1879  fhMCPi0PtOrigin->Fill(pt, 3.5, GetEventWeight());//resonances
1880  else fhMCPi0PtOrigin->Fill(pt, 7.5, GetEventWeight());//other?
1881 
1882  if(status!=11)
1883  {
1884  if (momstatus == 21) fhMCNotResonancePi0PtOrigin->Fill(pt, 0.5, GetEventWeight());//parton
1885  else if(mompdg < 22 ) fhMCNotResonancePi0PtOrigin->Fill(pt, 1.5, GetEventWeight());//quark
1886  else if(mompdg > 2100 && mompdg < 2210)
1887  fhMCNotResonancePi0PtOrigin->Fill(pt, 2.5, GetEventWeight());// resonances
1888  else if(mompdg == 221) fhMCNotResonancePi0PtOrigin->Fill(pt, 8.5, GetEventWeight());//eta
1889  else if(mompdg == 331) fhMCNotResonancePi0PtOrigin->Fill(pt, 9.5, GetEventWeight());//eta prime
1890  else if(mompdg == 213) fhMCNotResonancePi0PtOrigin->Fill(pt, 4.5, GetEventWeight());//rho
1891  else if(mompdg == 223) fhMCNotResonancePi0PtOrigin->Fill(pt, 5.5, GetEventWeight());//omega
1892  else if(mompdg >= 310 && mompdg <= 323)
1893  fhMCNotResonancePi0PtOrigin->Fill(pt, 6.5, GetEventWeight());//k0S, k+-,k*
1894  else if(mompdg == 130) fhMCNotResonancePi0PtOrigin->Fill(pt, 6.5, GetEventWeight());//k0L
1895  else if(momstatus == 11 || momstatus == 12 )
1896  fhMCNotResonancePi0PtOrigin->Fill(pt, 3.5, GetEventWeight());//resonances
1897  else fhMCNotResonancePi0PtOrigin->Fill(pt, 7.5, GetEventWeight());//other?
1898  }
1899  }//pi0 mass region
1900  }
1901  }
1902  else if(ancPDG==221)
1903  {//Eta
1904  mcIndex = 3;
1905  if(fMultiCutAnaSim)
1906  {
1907  for(Int_t ipt=0; ipt<fNPtCuts; ipt++)
1908  {
1909  for(Int_t icell=0; icell<fNCellNCuts; icell++)
1910  {
1911  for(Int_t iasym=0; iasym<fNAsymCuts; iasym++)
1912  {
1913  Int_t index = ((ipt*fNCellNCuts)+icell)*fNAsymCuts + iasym;
1914  if(pt1 > fPtCuts[ipt] && pt2 > fPtCuts[ipt] &&
1915  asym < fAsymCuts[iasym] &&
1916  ncell1 >= fCellNCuts[icell] && ncell2 >= fCellNCuts[icell])
1917  {
1918  fhMCEtaMassPtRec [index]->Fill(pt , mass, GetEventWeight());
1919  fhMCEtaMassPtTrue[index]->Fill(fPi0Mom.Pt(), mass, GetEventWeight());
1920  if(mass < 0.65 && mass > 0.45)
1921  fhMCEtaPtTruePtRec[index]->Fill(fPi0Mom.Pt(), pt, GetEventWeight());
1922  }//pass the different cuts
1923  }// pid bit cut loop
1924  }// icell loop
1925  }// pt cut loop
1926  } //Multi cut ana sim
1927  else
1928  {
1929  fhMCEtaMassPtTrue[0]->Fill(fPi0Mom.Pt(), mass, GetEventWeight());
1930  if(mass < 0.65 && mass > 0.45) fhMCEtaPtTruePtRec[0]->Fill(fPi0Mom.Pt(), pt, GetEventWeight());
1931 
1932  if(GetReader()->ReadStack())
1933  {
1934  TParticle* ancestor = GetMCStack()->Particle(ancLabel);
1935  momindex = ancestor->GetFirstMother();
1936  if(momindex < 0) return;
1937  TParticle* mother = GetMCStack()->Particle(momindex);
1938  mompdg = TMath::Abs(mother->GetPdgCode());
1939  momstatus = mother->GetStatusCode();
1940  prodR = mother->R();
1941  }
1942  else
1943  {
1944  TClonesArray * mcparticles = GetReader()->GetAODMCParticles();
1945  AliAODMCParticle* ancestor = (AliAODMCParticle *) mcparticles->At(ancLabel);
1946  momindex = ancestor->GetMother();
1947  if(momindex < 0) return;
1948  AliAODMCParticle* mother = (AliAODMCParticle *) mcparticles->At(momindex);
1949  mompdg = TMath::Abs(mother->GetPdgCode());
1950  momstatus = mother->GetStatus();
1951  prodR = TMath::Sqrt(mother->Xv()*mother->Xv()+mother->Yv()*mother->Yv());
1952  }
1953 
1954  fhMCEtaProdVertex->Fill(pt, prodR, GetEventWeight());
1955 
1956  if (momstatus == 21 ) fhMCEtaPtOrigin->Fill(pt, 0.5, GetEventWeight());//parton
1957  else if(mompdg < 22 ) fhMCEtaPtOrigin->Fill(pt, 1.5, GetEventWeight());//quark
1958  else if(mompdg > 2100 && mompdg < 2210)
1959  fhMCEtaPtOrigin->Fill(pt, 2.5, GetEventWeight());//qq resonances
1960  else if(mompdg == 331) fhMCEtaPtOrigin->Fill(pt, 5.5, GetEventWeight());//eta prime
1961  else if(momstatus == 11 || momstatus == 12 )
1962  fhMCEtaPtOrigin->Fill(pt, 3.5, GetEventWeight());//resonances
1963  else fhMCEtaPtOrigin->Fill(pt, 4.5, GetEventWeight());//stable, conversions?
1964  //printf("Other Meson pdg %d, Mother %s, pdg %d, status %d\n",pdg, TDatabasePDG::Instance()->GetParticle(mompdg)->GetName(),mompdg, momstatus );
1965 
1966  }// eta mass region
1967  }
1968  else if(ancPDG==-2212)
1969  {//AProton
1970  mcIndex = 4;
1971  }
1972  else if(ancPDG==-2112)
1973  {//ANeutron
1974  mcIndex = 5;
1975  }
1976  else if(TMath::Abs(ancPDG)==13)
1977  {//muons
1978  mcIndex = 6;
1979  }
1980  else if (TMath::Abs(ancPDG) > 100 && ancLabel > 7)
1981  {
1982  if(ancStatus==1)
1983  {//Stable particles, converted? not decayed resonances
1984  mcIndex = 6;
1985  }
1986  else
1987  {//resonances and other decays, more hadron conversions?
1988  mcIndex = 7;
1989  }
1990  }
1991  else
1992  {//Partons, colliding protons, strings, intermediate corrections
1993  if(ancStatus==11 || ancStatus==12)
1994  {//String fragmentation
1995  mcIndex = 8;
1996  }
1997  else if (ancStatus==21){
1998  if(ancLabel < 2)
1999  {//Colliding protons
2000  mcIndex = 11;
2001  }//colliding protons
2002  else if(ancLabel < 6)
2003  {//partonic initial states interactions
2004  mcIndex = 9;
2005  }
2006  else if(ancLabel < 8)
2007  {//Final state partons radiations?
2008  mcIndex = 10;
2009  }
2010  // else {
2011  // printf("AliAnaPi0::FillMCVersusRecDataHistograms() - Check ** Common ancestor label %d, pdg %d, name %s, status %d; \n",
2012  // ancLabel,ancPDG,TDatabasePDG::Instance()->GetParticle(ancPDG)->GetName(),ancStatus);
2013  // }
2014  }//status 21
2015  //else {
2016  // printf("AliAnaPi0::FillMCVersusRecDataHistograms() - Check *** Common ancestor label %d, pdg %d, name %s, status %d; \n",
2017  // ancLabel,ancPDG,TDatabasePDG::Instance()->GetParticle(ancPDG)->GetName(),ancStatus);
2018  // }
2019  }
2020  }//ancLabel > -1
2021  else
2022  { //ancLabel <= -1
2023  //printf("Not related at all label = %d\n",ancLabel);
2024  AliDebug(1,"Common ancestor not found");
2025 
2026  mcIndex = 12;
2027  }
2028 
2029  if(mcIndex >= 0 && mcIndex < 13)
2030  {
2031  fhMCOrgMass [mcIndex]->Fill(pt, mass, GetEventWeight());
2032  fhMCOrgAsym [mcIndex]->Fill(pt, asym, GetEventWeight());
2033  fhMCOrgDeltaEta[mcIndex]->Fill(pt, deta, GetEventWeight());
2034  fhMCOrgDeltaPhi[mcIndex]->Fill(pt, dphi, GetEventWeight());
2035  }
2036 }
2037 
2038 //__________________________________________
2042 //__________________________________________
2043 void AliAnaPi0::MakeAnalysisFillHistograms()
2044 {
2045  // In case of simulated data, fill acceptance histograms
2046  if(IsDataMC())
2047  {
2048  FillAcceptanceHistograms();
2049 
2050  if(fFillOnlyMCAcceptanceHisto) return;
2051  }
2052 
2053 //if (GetReader()->GetEventNumber()%10000 == 0)
2054 // printf("--- Event %d ---\n",GetReader()->GetEventNumber());
2055 
2056  if(!GetInputAODBranch())
2057  {
2058  AliFatal(Form("No input aod photons in AOD with name branch < %s >, STOP",GetInputAODName().Data()));
2059  return;
2060  }
2061 
2062  // Init some variables
2063  Int_t nPhot = GetInputAODBranch()->GetEntriesFast() ;
2064 
2065  AliDebug(1,Form("Photon entries %d", nPhot));
2066 
2067  // If less than photon 2 entries in the list, skip this event
2068  if ( nPhot < 2 )
2069  {
2070  AliDebug(1,Form("nPhotons %d, cent bin %d continue to next event",nPhot, GetEventCentrality()));
2071 
2072  if(GetNCentrBin() > 1)
2073  fhCentralityNoPair->Fill(GetEventCentrality() * GetNCentrBin() / GetReader()->GetCentralityOpt(), GetEventWeight());
2074 
2075  return ;
2076  }
2077 
2078  Int_t ncentr = GetNCentrBin();
2079  if(GetNCentrBin()==0) ncentr = 1;
2080 
2081  //Init variables
2082  Int_t module1 = -1;
2083  Int_t module2 = -1;
2084  Double_t vert[] = {0.0, 0.0, 0.0} ; //vertex
2085  Int_t evtIndex1 = 0 ;
2086  Int_t currentEvtIndex = -1;
2087  Int_t curCentrBin = GetEventCentralityBin();
2088  //Int_t curVzBin = GetEventVzBin();
2089  //Int_t curRPBin = GetEventRPBin();
2090  Int_t eventbin = GetEventMixBin();
2091 
2092  if(eventbin > GetNCentrBin()*GetNZvertBin()*GetNRPBin())
2093  {
2094  AliWarning(Form("Mix Bin not expected: cen bin %d, z bin %d, rp bin %d, total bin %d, Event Centrality %d, z vertex %2.3f, Reaction Plane %2.3f",
2095  GetEventCentralityBin(),GetEventVzBin(), GetEventRPBin(),eventbin,GetEventCentrality(),vert[2],GetEventPlaneAngle()));
2096  return;
2097  }
2098 
2099  // Get shower shape information of clusters
2100 // TObjArray *clusters = 0;
2101 // if (GetCalorimeter()==kEMCAL) clusters = GetEMCALClusters();
2102 // else if(GetCalorimeter()==kPHOS ) clusters = GetPHOSClusters() ;
2103 
2104  //---------------------------------
2105  // First loop on photons/clusters
2106  //---------------------------------
2107  for(Int_t i1=0; i1<nPhot-1; i1++)
2108  {
2109  AliAODPWG4Particle * p1 = (AliAODPWG4Particle*) (GetInputAODBranch()->At(i1)) ;
2110 
2111  // Select photons within a pT range
2112  if ( p1->Pt() < GetMinPt() || p1->Pt() > GetMaxPt() ) continue ;
2113 
2114  //printf("AliAnaPi0::MakeAnalysisFillHistograms() : cluster1 id %d/%d\n",i1,nPhot-1);
2115 
2116  // get the event index in the mixed buffer where the photon comes from
2117  // in case of mixing with analysis frame, not own mixing
2118  evtIndex1 = GetEventIndex(p1, vert) ;
2119  if ( evtIndex1 == -1 )
2120  return ;
2121  if ( evtIndex1 == -2 )
2122  continue ;
2123 
2124  // Only effective in case of mixed event frame
2125  if(TMath::Abs(vert[2]) > GetZvertexCut()) continue ; //vertex cut
2126 
2127  if (evtIndex1 != currentEvtIndex)
2128  {
2129  // Fill event bin info
2130  if( DoOwnMix() ) fhEventBin->Fill(eventbin, GetEventWeight()) ;
2131 
2132  if( IsHighMultiplicityAnalysisOn() )
2133  {
2134  fhCentrality->Fill(curCentrBin, GetEventWeight());
2135 
2136  if( GetEventPlane() )
2137  fhEventPlaneResolution->Fill(curCentrBin, TMath::Cos(2.*GetEventPlane()->GetQsubRes()), GetEventWeight());
2138  }
2139 
2140  currentEvtIndex = evtIndex1 ;
2141  }
2142 
2143  //printf("AliAnaPi0::MakeAnalysisFillHistograms(): Photon 1 Evt %d Vertex : %f,%f,%f\n",evtIndex1, GetVertex(evtIndex1)[0] ,GetVertex(evtIndex1)[1],GetVertex(evtIndex1)[2]);
2144 
2145  //Get the momentum of this cluster
2146  fPhotonMom1.SetPxPyPzE(p1->Px(),p1->Py(),p1->Pz(),p1->E());
2147 
2148  //Get (Super)Module number of this cluster
2149  module1 = p1->GetSModNumber();// GetModuleNumber(p1);
2150 
2151  //------------------------------------------
2152  // Recover original cluster
2153 // Int_t iclus1 = -1 ;
2154 // AliVCluster * cluster1 = FindCluster(clusters,p1->GetCaloLabel(0),iclus1);
2155 // if(!cluster1) AliWarning("Cluster1 not found!");
2156 
2157  //---------------------------------
2158  // Second loop on photons/clusters
2159  //---------------------------------
2160  for(Int_t i2=i1+1; i2<nPhot; i2++)
2161  {
2162  AliAODPWG4Particle * p2 = (AliAODPWG4Particle*) (GetInputAODBranch()->At(i2)) ;
2163 
2164  // Select photons within a pT range
2165  if ( p2->Pt() < GetMinPt() || p2->Pt() > GetMaxPt() ) continue ;
2166 
2167  //printf("AliAnaPi0::MakeAnalysisFillHistograms() : cluster2 i %d/%d\n",i2,nPhot);
2168 
2169  //In case of mixing frame, check we are not in the same event as the first cluster
2170  Int_t evtIndex2 = GetEventIndex(p2, vert) ;
2171  if ( evtIndex2 == -1 )
2172  return ;
2173  if ( evtIndex2 == -2 )
2174  continue ;
2175  if (GetMixedEvent() && (evtIndex1 == evtIndex2))
2176  continue ;
2177 
2178 // //------------------------------------------
2179 // // Recover original cluster
2180 // Int_t iclus2 = -1;
2181 // AliVCluster * cluster2 = FindCluster(clusters,p2->GetCaloLabel(0),iclus2,iclus1+1);
2182 // // start new loop from iclus1+1 to gain some time
2183 // if(!cluster2) AliWarning("Cluster2 not found!");
2184 //
2185 // // Get the TOF,l0 and ncells from the clusters
2186 // Float_t tof1 = -1;
2187 // Float_t l01 = -1;
2188 // Int_t ncell1 = 0;
2189 // if(cluster1)
2190 // {
2191 // tof1 = cluster1->GetTOF()*1e9;
2192 // l01 = cluster1->GetM02();
2193 // ncell1 = cluster1->GetNCells();
2194 // //printf("cluster1: E %2.2f (%2.2f), l0 %2.2f, tof %2.2f\n",cluster1->E(),p1->E(),l01,tof1);
2195 // }
2196 // //else printf("cluster1 not available: calo label %d / %d, cluster ID %d\n",
2197 // // p1->GetCaloLabel(0),(GetReader()->GetInputEvent())->GetNumberOfCaloClusters()-1,cluster1->GetID());
2198 //
2199 // Float_t tof2 = -1;
2200 // Float_t l02 = -1;
2201 // Int_t ncell2 = 0;
2202 // if(cluster2)
2203 // {
2204 // tof2 = cluster2->GetTOF()*1e9;
2205 // l02 = cluster2->GetM02();
2206 // ncell2 = cluster2->GetNCells();
2207 // //printf("cluster2: E %2.2f (%2.2f), l0 %2.2f, tof %2.2f\n",cluster2->E(),p2->E(),l02,tof2);
2208 // }
2209 // //else printf("cluster2 not available: calo label %d / %d, cluster ID %d\n",
2210 // // p2->GetCaloLabel(0),(GetReader()->GetInputEvent())->GetNumberOfCaloClusters()-1,cluster2->GetID());
2211 //
2212 // if(cluster1 && cluster2)
2213 // {
2214 // Double_t t12diff = tof1-tof2;
2215 // if(TMath::Abs(t12diff) > GetPairTimeCut()) continue;
2216 // }
2217 
2218  Float_t tof1 = p1->GetTime();
2219  Float_t l01 = p1->GetM02();
2220  Int_t ncell1 = p1->GetNCells();
2221  //printf("cluster1: E %2.2f, l0 %2.2f, tof %2.2f\n",p1->E(),l01,tof1);
2222 
2223  Float_t tof2 = p2->GetTime();
2224  Float_t l02 = p2->GetM02();
2225  Int_t ncell2 = p2->GetNCells();
2226  //printf("cluster2: E %2.2f, l0 %2.2f, tof %2.2f\n",p2->E(),l02,tof2);
2227 
2228  Double_t t12diff = tof1-tof2;
2229  fhEPairDiffTime->Fill((fPhotonMom1 + fPhotonMom2).Pt(), t12diff, GetEventWeight());
2230  if(TMath::Abs(t12diff) > GetPairTimeCut()) continue;
2231 
2232  //------------------------------------------
2233 
2234  //printf("AliAnaPi0::MakeAnalysisFillHistograms(): Photon 2 Evt %d Vertex : %f,%f,%f\n",evtIndex2, GetVertex(evtIndex2)[0] ,GetVertex(evtIndex2)[1],GetVertex(evtIndex2)[2]);
2235 
2236  // Get the momentum of this cluster
2237  fPhotonMom2.SetPxPyPzE(p2->Px(),p2->Py(),p2->Pz(),p2->E());
2238 
2239  // Get module number
2240  module2 = p2->GetSModNumber(); //GetModuleNumber(p2);
2241 
2242  //---------------------------------
2243  // Get pair kinematics
2244  //---------------------------------
2245  Double_t m = (fPhotonMom1 + fPhotonMom2).M() ;
2246  Double_t pt = (fPhotonMom1 + fPhotonMom2).Pt();
2247  Double_t deta = fPhotonMom1.Eta() - fPhotonMom2.Eta();
2248  Double_t dphi = fPhotonMom1.Phi() - fPhotonMom2.Phi();
2249  Double_t a = TMath::Abs(p1->E()-p2->E())/(p1->E()+p2->E()) ;
2250 
2251  AliDebug(2,Form("E: fPhotonMom1 %f, fPhotonMom2 %f; Pair: pT %f, mass %f, a %f", p1->E(), p2->E(), (fPhotonMom1 + fPhotonMom2).E(),m,a));
2252 
2253  //--------------------------------
2254  // Opening angle selection
2255  //--------------------------------
2256  // Check if opening angle is too large or too small compared to what is expected
2257  Double_t angle = fPhotonMom1.Angle(fPhotonMom2.Vect());
2258  if(fUseAngleEDepCut && !GetNeutralMesonSelection()->IsAngleInWindow((fPhotonMom1+fPhotonMom2).E(),angle+0.05))
2259  {
2260  AliDebug(2,Form("Real pair angle %f not in E %f window",angle, (fPhotonMom1+fPhotonMom2).E()));
2261  continue;
2262  }
2263 
2264  if(fUseAngleCut && (angle < fAngleCut || angle > fAngleMaxCut))
2265  {
2266  AliDebug(2,Form("Real pair cut %f < angle %f < cut %f",fAngleCut, angle, fAngleMaxCut));
2267  continue;
2268  }
2269 
2270  //-------------------------------------------------------------------------------------------------
2271  // Fill module dependent histograms, put a cut on assymmetry on the first available cut in the array
2272  //-------------------------------------------------------------------------------------------------
2273  if(a < fAsymCuts[0] && fFillSMCombinations)
2274  {
2275  if(module1==module2 && module1 >=0 && module1<fNModules)
2276  fhReMod[module1]->Fill(pt, m, GetEventWeight()) ;
2277 
2278  if (GetCalorimeter() == kEMCAL )
2279  {
2280  // Same sector
2281  Int_t j=0;
2282  for(Int_t i = 0; i < fNModules/2; i++)
2283  {
2284  j=2*i;
2285  if((module1==j && module2==j+1) || (module1==j+1 && module2==j)) fhReSameSectorEMCALMod[i]->Fill(pt, m, GetEventWeight()) ;
2286  }
2287 
2288  // Same side
2289  for(Int_t i = 0; i < fNModules-2; i++){
2290  if((module1==i && module2==i+2) || (module1==i+2 && module2==i)) fhReSameSideEMCALMod[i]->Fill(pt, m, GetEventWeight());
2291  }
2292  } // EMCAL
2293  else
2294  { // PHOS
2295  if((module1==0 && module2==1) || (module1==1 && module2==0)) fhReDiffPHOSMod[0]->Fill(pt, m, GetEventWeight()) ;
2296  if((module1==0 && module2==2) || (module1==2 && module2==0)) fhReDiffPHOSMod[1]->Fill(pt, m, GetEventWeight()) ;
2297  if((module1==1 && module2==2) || (module1==2 && module2==1)) fhReDiffPHOSMod[2]->Fill(pt, m, GetEventWeight()) ;
2298  } // PHOS
2299  }
2300 
2301  // In case we want only pairs in same (super) module, check their origin.
2302  Bool_t ok = kTRUE;
2303 
2304  if(fSameSM && module1!=module2) ok=kFALSE;
2305 
2306  if(ok)
2307  {
2308  // Check if one of the clusters comes from a conversion
2309  if(fCheckConversion)
2310  {
2311  if (p1->IsTagged() && p2->IsTagged()) fhReConv2->Fill(pt, m, GetEventWeight());
2312  else if(p1->IsTagged() || p2->IsTagged()) fhReConv ->Fill(pt, m, GetEventWeight());
2313  }
2314 
2315  // Fill shower shape cut histograms
2316  if(fFillSSCombinations)
2317  {
2318  if ( l01 > 0.01 && l01 < 0.4 &&
2319  l02 > 0.01 && l02 < 0.4 ) fhReSS[0]->Fill(pt, m, GetEventWeight()); // Tight
2320  else if( l01 > 0.4 && l02 > 0.4 ) fhReSS[1]->Fill(pt, m, GetEventWeight()); // Loose
2321  else if( l01 > 0.01 && l01 < 0.4 && l02 > 0.4 ) fhReSS[2]->Fill(pt, m, GetEventWeight()); // Both
2322  else if( l02 > 0.01 && l02 < 0.4 && l01 > 0.4 ) fhReSS[2]->Fill(pt, m, GetEventWeight()); // Both
2323  }
2324 
2325  // Fill histograms for different bad channel distance, centrality, assymmetry cut and pid bit
2326  for(Int_t ipid=0; ipid<fNPIDBits; ipid++)
2327  {
2328  if((p1->IsPIDOK(fPIDBits[ipid],AliCaloPID::kPhoton)) && (p2->IsPIDOK(fPIDBits[ipid],AliCaloPID::kPhoton)))
2329  {
2330  for(Int_t iasym=0; iasym < fNAsymCuts; iasym++)
2331  {
2332  if(a < fAsymCuts[iasym])
2333  {
2334  Int_t index = ((curCentrBin*fNPIDBits)+ipid)*fNAsymCuts + iasym;
2335  //printf("index %d :(cen %d * nPID %d + ipid %d)*nasym %d + iasym %d - max index %d\n",index,curCentrBin,fNPIDBits,ipid,fNAsymCuts,iasym, curCentrBin*fNPIDBits*fNAsymCuts);
2336 
2337  if(index < 0 || index >= ncentr*fNPIDBits*fNAsymCuts) continue ;
2338 
2339  fhRe1 [index]->Fill(pt, m, GetEventWeight());
2340  if(fMakeInvPtPlots)fhReInvPt1[index]->Fill(pt, m, 1./pt * GetEventWeight()) ;
2341 
2342  if(fFillBadDistHisto)
2343  {
2344  if(p1->DistToBad()>0 && p2->DistToBad()>0)
2345  {
2346  fhRe2 [index]->Fill(pt, m, GetEventWeight()) ;
2347  if(fMakeInvPtPlots)fhReInvPt2[index]->Fill(pt, m, 1./pt * GetEventWeight()) ;
2348 
2349  if(p1->DistToBad()>1 && p2->DistToBad()>1)
2350  {
2351  fhRe3 [index]->Fill(pt, m, GetEventWeight()) ;
2352  if(fMakeInvPtPlots)fhReInvPt3[index]->Fill(pt, m, 1./pt * GetEventWeight()) ;
2353  }// bad 3
2354  }// bad2
2355  }// Fill bad dist histos
2356  }//assymetry cut
2357  }// asymmetry cut loop
2358  }// bad 1
2359  }// pid bit loop
2360 
2361  // Fill histograms with opening angle
2362  if(fFillAngleHisto)
2363  {
2364  fhRealOpeningAngle ->Fill(pt, angle, GetEventWeight());
2365  fhRealCosOpeningAngle->Fill(pt, TMath::Cos(angle), GetEventWeight());
2366  }
2367 
2368  // Fill histograms with pair assymmetry
2369  if(fFillAsymmetryHisto)
2370  {
2371  fhRePtAsym->Fill(pt, a, GetEventWeight());
2372  if(m > 0.10 && m < 0.17) fhRePtAsymPi0->Fill(pt, a, GetEventWeight());
2373  if(m > 0.45 && m < 0.65) fhRePtAsymEta->Fill(pt, a, GetEventWeight());
2374  }
2375 
2376  //---------
2377  // MC data
2378  //---------
2379  // Do some MC checks on the origin of the pair, is there any common ancestor and if there is one, who?
2380  if(IsDataMC())
2381  {
2382  if(GetMCAnalysisUtils()->CheckTagBit(p1->GetTag(),AliMCAnalysisUtils::kMCConversion) &&
2383  GetMCAnalysisUtils()->CheckTagBit(p2->GetTag(),AliMCAnalysisUtils::kMCConversion))
2384  {
2385  fhReMCFromConversion->Fill(pt, m, GetEventWeight());
2386  }
2387  else if(!GetMCAnalysisUtils()->CheckTagBit(p1->GetTag(),AliMCAnalysisUtils::kMCConversion) &&
2388  !GetMCAnalysisUtils()->CheckTagBit(p2->GetTag(),AliMCAnalysisUtils::kMCConversion))
2389  {
2390  fhReMCFromNotConversion->Fill(pt, m, GetEventWeight());
2391  }
2392  else
2393  {
2394  fhReMCFromMixConversion->Fill(pt, m, GetEventWeight());
2395  }
2396 
2397  if(fFillOriginHisto)
2398  FillMCVersusRecDataHistograms(p1->GetLabel(), p2->GetLabel(),p1->Pt(), p2->Pt(),ncell1, ncell2, m, pt, a,deta, dphi);
2399  }
2400 
2401  //-----------------------
2402  // Multi cuts analysis
2403  //-----------------------
2404  if(fMultiCutAna)
2405  {
2406  // Histograms for different PID bits selection
2407  for(Int_t ipid=0; ipid<fNPIDBits; ipid++)
2408  {
2409  if(p1->IsPIDOK(fPIDBits[ipid],AliCaloPID::kPhoton) &&
2410  p2->IsPIDOK(fPIDBits[ipid],AliCaloPID::kPhoton)) fhRePIDBits[ipid]->Fill(pt, m, GetEventWeight()) ;
2411 
2412  //printf("ipt %d, ipid%d, name %s\n",ipt, ipid, fhRePtPIDCuts[ipt*fNPIDBitsBits+ipid]->GetName());
2413  } // pid bit cut loop
2414 
2415  // Several pt,ncell and asymmetry cuts
2416  for(Int_t ipt=0; ipt<fNPtCuts; ipt++)
2417  {
2418  for(Int_t icell=0; icell<fNCellNCuts; icell++)
2419  {
2420  for(Int_t iasym=0; iasym<fNAsymCuts; iasym++)
2421  {
2422  Int_t index = ((ipt*fNCellNCuts)+icell)*fNAsymCuts + iasym;
2423  if(p1->E() > fPtCuts[ipt] && p2->E() > fPtCuts[ipt] &&
2424  a < fAsymCuts[iasym] &&
2425  ncell1 >= fCellNCuts[icell] && ncell2 >= fCellNCuts[icell])
2426  {
2427  fhRePtNCellAsymCuts[index]->Fill(pt, m, GetEventWeight()) ;
2428  //printf("ipt %d, icell%d, iasym %d, name %s\n",ipt, icell, iasym, fhRePtNCellAsymCuts[((ipt*fNCellNCuts)+icell)*fNAsymCuts + iasym]->GetName());
2429  if(fFillSMCombinations && module1==module2)
2430  {
2431  fhRePtNCellAsymCutsSM[module1][index]->Fill(pt, m, GetEventWeight()) ;
2432  }
2433  }
2434  }// pid bit cut loop
2435  }// icell loop
2436  }// pt cut loop
2437 
2438  if(GetHistogramRanges()->GetHistoTrackMultiplicityBins())
2439  {
2440  for(Int_t iasym = 0; iasym < fNAsymCuts; iasym++)
2441  {
2442  if(a < fAsymCuts[iasym]) fhRePtMult[iasym]->Fill(pt, GetTrackMultiplicity(), m, GetEventWeight()) ;
2443  }
2444  }
2445  }// multiple cuts analysis
2446  }// ok if same sm
2447  }// second same event particle
2448  }// first cluster
2449 
2450  //-------------------------------------------------------------
2451  // Mixing
2452  //-------------------------------------------------------------
2453  if(DoOwnMix())
2454  {
2455  // Recover events in with same characteristics as the current event
2456 
2457  // Check that the bin exists, if not (bad determination of RP, centrality or vz bin) do nothing
2458  if(eventbin < 0) return ;
2459 
2460  TList * evMixList=fEventsList[eventbin] ;
2461 
2462  if(!evMixList)
2463  {
2464  AliWarning(Form("Mix event list not available, bin %d",eventbin));
2465  return;
2466  }
2467 
2468  Int_t nMixed = evMixList->GetSize() ;
2469  for(Int_t ii=0; ii<nMixed; ii++)
2470  {
2471  TClonesArray* ev2= (TClonesArray*) (evMixList->At(ii));
2472  Int_t nPhot2=ev2->GetEntriesFast() ;
2473  Double_t m = -999;
2474  AliDebug(1,Form("Mixed event %d photon entries %d, centrality bin %d",ii, nPhot2, GetEventCentralityBin()));
2475 
2476  fhEventMixBin->Fill(eventbin, GetEventWeight()) ;
2477 
2478  //---------------------------------
2479  // First loop on photons/clusters
2480  //---------------------------------
2481  for(Int_t i1=0; i1<nPhot; i1++)
2482  {
2483  AliAODPWG4Particle * p1 = (AliAODPWG4Particle*) (GetInputAODBranch()->At(i1)) ;
2484 
2485  if(fSameSM && GetModuleNumber(p1)!=module1) continue;
2486 
2487  //Get kinematics of cluster and (super) module of this cluster
2488  fPhotonMom1.SetPxPyPzE(p1->Px(),p1->Py(),p1->Pz(),p1->E());
2489  module1 = GetModuleNumber(p1);
2490 
2491  //---------------------------------
2492  // First loop on photons/clusters
2493  //---------------------------------
2494  for(Int_t i2=0; i2<nPhot2; i2++)
2495  {
2496  AliAODPWG4Particle * p2 = (AliAODPWG4Particle*) (ev2->At(i2)) ;
2497 
2498  // Get kinematics of second cluster and calculate those of the pair
2499  fPhotonMom2.SetPxPyPzE(p2->Px(),p2->Py(),p2->Pz(),p2->E());
2500  m = (fPhotonMom1+fPhotonMom2).M() ;
2501  Double_t pt = (fPhotonMom1 + fPhotonMom2).Pt();
2502  Double_t a = TMath::Abs(p1->E()-p2->E())/(p1->E()+p2->E()) ;
2503 
2504  // Check if opening angle is too large or too small compared to what is expected
2505  Double_t angle = fPhotonMom1.Angle(fPhotonMom2.Vect());
2506  if(fUseAngleEDepCut && !GetNeutralMesonSelection()->IsAngleInWindow((fPhotonMom1+fPhotonMom2).E(),angle+0.05))
2507  {
2508  AliDebug(2,Form("Mix pair angle %f not in E %f window",angle, (fPhotonMom1+fPhotonMom2).E()));
2509  continue;
2510  }
2511 
2512  if(fUseAngleCut && (angle < fAngleCut || angle > fAngleMaxCut))
2513  {
2514  AliDebug(2,Form("Mix pair angle %f < cut %f",angle,fAngleCut));
2515  continue;
2516  }
2517 
2518  AliDebug(2,Form("Mixed Event: pT: fPhotonMom1 %2.2f, fPhotonMom2 %2.2f; Pair: pT %2.2f, mass %2.3f, a %2.3f",p1->Pt(), p2->Pt(), pt,m,a));
2519 
2520  // In case we want only pairs in same (super) module, check their origin.
2521  module2 = GetModuleNumber(p2);
2522 
2523  //-------------------------------------------------------------------------------------------------
2524  //Fill module dependent histograms, put a cut on assymmetry on the first available cut in the array
2525  //-------------------------------------------------------------------------------------------------
2526  if(a < fAsymCuts[0] && fFillSMCombinations)
2527  {
2528  if(module1==module2 && module1 >=0 && module1<fNModules)
2529  fhMiMod[module1]->Fill(pt, m, GetEventWeight()) ;
2530 
2531  if(GetCalorimeter()==kEMCAL)
2532  {
2533  // Same sector
2534  Int_t j=0;
2535  for(Int_t i = 0; i < fNModules/2; i++)
2536  {
2537  j=2*i;
2538  if((module1==j && module2==j+1) || (module1==j+1 && module2==j)) fhMiSameSectorEMCALMod[i]->Fill(pt, m, GetEventWeight()) ;
2539  }
2540 
2541  // Same side
2542  for(Int_t i = 0; i < fNModules-2; i++)
2543  {
2544  if((module1==i && module2==i+2) || (module1==i+2 && module2==i)) fhMiSameSideEMCALMod[i]->Fill(pt, m, GetEventWeight());
2545  }
2546  } // EMCAL
2547  else
2548  { // PHOS
2549  if((module1==0 && module2==1) || (module1==1 && module2==0)) fhMiDiffPHOSMod[0]->Fill(pt, m, GetEventWeight()) ;
2550  if((module1==0 && module2==2) || (module1==2 && module2==0)) fhMiDiffPHOSMod[1]->Fill(pt, m, GetEventWeight()) ;
2551  if((module1==1 && module2==2) || (module1==2 && module2==1)) fhMiDiffPHOSMod[2]->Fill(pt, m, GetEventWeight()) ;
2552  } // PHOS
2553  }
2554 
2555  Bool_t ok = kTRUE;
2556  if(fSameSM && module1!=module2) ok=kFALSE;
2557  if(ok){
2558 
2559  // Check if one of the clusters comes from a conversion
2560  if(fCheckConversion)
2561  {
2562  if (p1->IsTagged() && p2->IsTagged()) fhMiConv2->Fill(pt, m, GetEventWeight());
2563  else if(p1->IsTagged() || p2->IsTagged()) fhMiConv ->Fill(pt, m, GetEventWeight());
2564  }
2565 
2566  // Fill histograms for different bad channel distance, centrality, assymmetry cut and pid bit
2567  for(Int_t ipid=0; ipid<fNPIDBits; ipid++)
2568  {
2569  if((p1->IsPIDOK(ipid,AliCaloPID::kPhoton)) && (p2->IsPIDOK(ipid,AliCaloPID::kPhoton)))
2570  {
2571  for(Int_t iasym=0; iasym < fNAsymCuts; iasym++)
2572  {
2573  if(a < fAsymCuts[iasym])
2574  {
2575  Int_t index = ((curCentrBin*fNPIDBits)+ipid)*fNAsymCuts + iasym;
2576 
2577  if(index < 0 || index >= ncentr*fNPIDBits*fNAsymCuts) continue ;
2578 
2579  fhMi1[index]->Fill(pt, m, GetEventWeight()) ;
2580 
2581  if(fMakeInvPtPlots)fhMiInvPt1[index]->Fill(pt, m, 1./pt * GetEventWeight()) ;
2582 
2583  if(fFillBadDistHisto)
2584  {
2585  if(p1->DistToBad()>0 && p2->DistToBad()>0)
2586  {
2587  fhMi2[index]->Fill(pt, m, GetEventWeight()) ;
2588  if(fMakeInvPtPlots)fhMiInvPt2[index]->Fill(pt, m, 1./pt * GetEventWeight()) ;
2589 
2590  if(p1->DistToBad()>1 && p2->DistToBad()>1)
2591  {
2592  fhMi3[index]->Fill(pt, m, GetEventWeight()) ;
2593  if(fMakeInvPtPlots)fhMiInvPt3[index]->Fill(pt, m, 1./pt * GetEventWeight()) ;
2594  }
2595  }
2596  }// Fill bad dist histo
2597  }//Asymmetry cut
2598  }// Asymmetry loop
2599  }//PID cut
2600  }//loop for histograms
2601 
2602  //-----------------------
2603  // Multi cuts analysis
2604  //-----------------------
2605  if(fMultiCutAna)
2606  {
2607  // Several pt,ncell and asymmetry cuts
2608 
2609  for(Int_t ipt=0; ipt<fNPtCuts; ipt++)
2610  {
2611  for(Int_t icell=0; icell<fNCellNCuts; icell++)
2612  {
2613  for(Int_t iasym=0; iasym<fNAsymCuts; iasym++)
2614  {
2615  Int_t index = ((ipt*fNCellNCuts)+icell)*fNAsymCuts + iasym;
2616  if(p1->Pt() > fPtCuts[ipt] && p2->Pt() > fPtCuts[ipt] &&
2617  a < fAsymCuts[iasym] //&&
2618  //p1->GetBtag() >= fCellNCuts[icell] && p2->GetBtag() >= fCellNCuts[icell] // trick, correct it.
2619  )
2620  {
2621  fhMiPtNCellAsymCuts[index]->Fill(pt, m, GetEventWeight()) ;
2622  //printf("ipt %d, icell%d, iasym %d, name %s\n",ipt, icell, iasym, fhRePtNCellAsymCuts[((ipt*fNCellNCuts)+icell)*fNAsymCuts + iasym]->GetName());
2623  }
2624  }// pid bit cut loop
2625  }// icell loop
2626  }// pt cut loop
2627  } // Multi cut ana
2628 
2629  // Fill histograms with opening angle
2630  if(fFillAngleHisto)
2631  {
2632  fhMixedOpeningAngle ->Fill(pt, angle, GetEventWeight());
2633  fhMixedCosOpeningAngle->Fill(pt, TMath::Cos(angle), GetEventWeight());
2634  }
2635  }//ok
2636  }// second cluster loop
2637  }//first cluster loop
2638  }//loop on mixed events
2639 
2640  //--------------------------------------------------------
2641  // Add the current event to the list of events for mixing
2642  //--------------------------------------------------------
2643 
2644  TClonesArray *currentEvent = new TClonesArray(*GetInputAODBranch());
2645  // Add current event to buffer and Remove redundant events
2646  if( currentEvent->GetEntriesFast() > 0 )
2647  {
2648  evMixList->AddFirst(currentEvent) ;
2649  currentEvent=0 ; //Now list of particles belongs to buffer and it will be deleted with buffer
2650  if( evMixList->GetSize() >= GetNMaxEvMix() )
2651  {
2652  TClonesArray * tmp = (TClonesArray*) (evMixList->Last()) ;
2653  evMixList->RemoveLast() ;
2654  delete tmp ;
2655  }
2656  }
2657  else
2658  { // empty event
2659  delete currentEvent ;
2660  currentEvent=0 ;
2661  }
2662  }// DoOwnMix
2663 
2664  AliDebug(1,"End fill histograms");
2665 }
2666 
2667 //________________________________________________________________________
2671 //________________________________________________________________________
2672 Int_t AliAnaPi0::GetEventIndex(AliAODPWG4Particle * part, Double_t * vert)
2673 {
2674  Int_t evtIndex = -1 ;
2675 
2676  if(GetReader()->GetDataType()!=AliCaloTrackReader::kMC)
2677  {
2678  if (GetMixedEvent())
2679  {
2680  evtIndex = GetMixedEvent()->EventIndexForCaloCluster(part->GetCaloLabel(0)) ;
2681  GetVertex(vert,evtIndex);
2682 
2683  if(TMath::Abs(vert[2])> GetZvertexCut())
2684  evtIndex = -2 ; //Event can not be used (vertex, centrality,... cuts not fulfilled)
2685  }
2686  else
2687  {
2688  // Single event
2689  GetVertex(vert);
2690 
2691  if(TMath::Abs(vert[2])> GetZvertexCut())
2692  evtIndex = -1 ; //Event can not be used (vertex, centrality,... cuts not fulfilled)
2693  else
2694  evtIndex = 0 ;
2695  }
2696  } // No MC reader
2697  else
2698  {
2699  evtIndex = 0;
2700  vert[0] = 0. ;
2701  vert[1] = 0. ;
2702  vert[2] = 0. ;
2703  }
2704 
2705  return evtIndex ;
2706 }
2707 
AliAnaPi0::fhPrimEtaY
TH2F * fhPrimEtaY
! Rapidity distribution of primary particles vs pT
Definition: AliAnaPi0.h:359
AliHistogramRanges::GetHistoPtMax
Float_t GetHistoPtMax() const
Definition: AliHistogramRanges.h:39
AliAnaPi0::fhRe3
TH2F ** fhRe3
[GetNCentrBin()*fNPIDBits*fNAsymCuts]
Definition: AliAnaPi0.h:260
AliAnaPi0::fhPrimEtaAccPtCentrality
TH2F * fhPrimEtaAccPtCentrality
! primary eta with accepted daughters reconstructed centrality vs pT
Definition: AliAnaPi0.h:371
ClassImp
ClassImp(AliAnalysisTaskTriggerRates) AliAnalysisTaskTriggerRates
Definition: AliAnalysisTaskTriggerRates.cxx:34
AliAnaPi0::fhPrimPi0Y
TH2F * fhPrimPi0Y
! Rapidity distribution of primary particles vs pT
Definition: AliAnaPi0.h:338
AliAnaPi0::fPhotonMom1Boost
TLorentzVector fPhotonMom1Boost
! Photon cluster momentum, temporary array
Definition: AliAnaPi0.h:195
AliAnaPi0::fNModules
Int_t fNModules
[GetNCentrBin()*GetNZvertBin()*GetNRPBin()]
Definition: AliAnaPi0.h:160
AliAnaPi0::fhPrimEtaYetaYcut
TH2F * fhPrimEtaYetaYcut
! PseudoRapidity distribution of primary particles vs pT, Y<1
Definition: AliAnaPi0.h:362
AliAnaPi0::fhPrimPi0AccPt
TH1F * fhPrimPi0AccPt
! Spectrum of primary with accepted daughters
Definition: AliAnaPi0.h:337
AliHistogramRanges::GetHistoPtMin
Float_t GetHistoPtMin() const
Definition: AliHistogramRanges.h:38
AliAnaPi0::fhPrimEtaAccPhi
TH2F * fhPrimEtaAccPhi
! Azimutal distribution of primary with accepted daughters vs pT
Definition: AliAnaPi0.h:365
AliAnaPi0::fPi0Mom
TLorentzVector fPi0Mom
! Pi0 cluster momentum, temporary array
Definition: AliAnaPi0.h:197
AliAnaPi0::fhReMCFromMixConversion
TH2F * fhReMCFromMixConversion
! Invariant mass of 2 clusters one from conversion and the other not
Definition: AliAnaPi0.h:422
AliCaloTrackReader.h
AliAnaPi0::fPIDBits
Int_t fPIDBits[10]
Array with different PID bits.
Definition: AliAnaPi0.h:177
AliAnaCaloTrackCorrBaseClass::AddToHistogramsName
virtual void AddToHistogramsName(TString add)
Definition: AliAnaCaloTrackCorrBaseClass.h:85
AliAnaPi0::fhReConv
TH2F * fhReConv
[fNModules-1]
Definition: AliAnaPi0.h:239
AliAnaPi0::fFillArmenterosThetaStar
Bool_t fFillArmenterosThetaStar
Fill armenteros histograms.
Definition: AliAnaPi0.h:189
AliAnaPi0::fNPIDBits
Int_t fNPIDBits
Number of possible PID bit combinations.
Definition: AliAnaPi0.h:176
AliAnaPi0::fPhotonMom2
TLorentzVector fPhotonMom2
! Photon cluster momentum, temporary array
Definition: AliAnaPi0.h:196
title
const char * title
Definition: MakeQAPdf.C:26
AliAnaPi0::fCellNCuts
Int_t fCellNCuts[10]
Array with different cell number cluster cuts.
Definition: AliAnaPi0.h:175
AliAnaPi0::fhMCOrgDeltaEta
TH2F * fhMCOrgDeltaEta[13]
! Delta Eta vs pt of real pairs, check common origin of pair
Definition: AliAnaPi0.h:387
AliCaloPID.h
AliAnaPi0::fhReSS
TH2F * fhReSS[3]
[fNAsymCuts]
Definition: AliAnaPi0.h:307
AliAnaPi0::fhPrimEtaCosOpeningAngle
TH2F * fhPrimEtaCosOpeningAngle
! Cosinus of opening angle of pair version pair energy, eta primaries
Definition: AliAnaPi0.h:368
AliAnaPi0::fhPrimPi0PtStatus
TH2F * fhPrimPi0PtStatus
! Spectrum of generated pi0 vs pi0 status
Definition: AliAnaPi0.h:379
AliAnaPi0::fCheckConversion
Bool_t fCheckConversion
Fill histograms with tagged photons as conversion.
Definition: AliAnaPi0.h:183
AliAnaCaloTrackCorrBaseClass::GetPairTimeCut
virtual Float_t GetPairTimeCut() const
Time cut in ns.
Definition: AliAnaCaloTrackCorrBaseClass.h:221
AliAnaCaloTrackCorrBaseClass::GetVertex
virtual void GetVertex(Double_t vertex[3]) const
Definition: AliAnaCaloTrackCorrBaseClass.h:271
AliHistogramRanges::GetHistoTrackMultiplicityMax
Int_t GetHistoTrackMultiplicityMax() const
Definition: AliHistogramRanges.h:133
AliAnaPi0::fhPrimPi0E
TH1F * fhPrimPi0E
! Spectrum of Primary
Definition: AliAnaPi0.h:334
AliAnaPi0::fUseAngleEDepCut
Bool_t fUseAngleEDepCut
Select pairs depending on their opening angle.
Definition: AliAnaPi0.h:163
AliAnaPi0::fhMiSameSectorEMCALMod
TH2F ** fhMiSameSectorEMCALMod
[fNModules-2]
Definition: AliAnaPi0.h:232
AliAnaPi0
Selected photon clusters invariant mass analysis.
Definition: AliAnaPi0.h:38
AliAnaPi0::fhCentrality
TH1F * fhCentrality
! Histogram with centrality bins with at least one pare
Definition: AliAnaPi0.h:319
AliAnaCaloTrackCorrBaseClass::SetInputAODName
virtual void SetInputAODName(TString name)
Definition: AliAnaCaloTrackCorrBaseClass.h:124
AliAnaPi0::fhPrimPi0PtCentrality
TH2F * fhPrimPi0PtCentrality
! primary pi0 reconstructed centrality vs pT
Definition: AliAnaPi0.h:348
AliAnaPi0::fhPrimEtaPtCentrality
TH2F * fhPrimEtaPtCentrality
! primary eta reconstructed centrality vs pT
Definition: AliAnaPi0.h:369
AliAnaCaloTrackCorrBaseClass::GetZvertexCut
virtual Float_t GetZvertexCut() const
Maximal number of events for mixin.
Definition: AliAnaCaloTrackCorrBaseClass.h:235
AliAnaCaloTrackCorrBaseClass::GetEventPlaneAngle
virtual Double_t GetEventPlaneAngle() const
Definition: AliAnaCaloTrackCorrBaseClass.h:112
AliAnaPi0::fhReMCFromConversion
TH2F * fhReMCFromConversion
! Invariant mass of 2 clusters originated in conversions
Definition: AliAnaPi0.h:420
AliAnaPi0::fhReSameSideEMCALMod
TH2F ** fhReSameSideEMCALMod
[fNModules]
Definition: AliAnaPi0.h:217
AliAnaPi0::fhPrimPi0OpeningAngle
TH2F * fhPrimPi0OpeningAngle
! Opening angle of pair versus pair energy, primaries
Definition: AliAnaPi0.h:345
AliAnaCaloTrackCorrBaseClass::GetNeutralMesonSelection
virtual AliNeutralMesonSelection * GetNeutralMesonSelection()
Definition: AliAnaCaloTrackCorrBaseClass.h:332
AliAnaPi0::fhRealOpeningAngle
TH2F * fhRealOpeningAngle
! Opening angle of pair versus pair energy
Definition: AliAnaPi0.h:326
AliAnaPi0::fhMCOrgMass
TH2F * fhMCOrgMass[13]
! Mass vs pt of real pairs, check common origin of pair
Definition: AliAnaPi0.h:385
AliHistogramRanges::GetHistoMassBins
Int_t GetHistoMassBins() const
Definition: AliHistogramRanges.h:97
AliHistogramRanges::GetHistoPhiBins
Int_t GetHistoPhiBins() const
Definition: AliHistogramRanges.h:60
AliAnaPi0::fhRePtAsymPi0
TH2F * fhRePtAsymPi0
! REAL two-photon pt vs asymmetry, close to pi0 mass
Definition: AliAnaPi0.h:312
AliFiducialCut.h
AliAnaPi0::fPhotonMom1
TLorentzVector fPhotonMom1
! Photon cluster momentum, temporary array
Definition: AliAnaPi0.h:194
AliAnaPi0::fhRealCosOpeningAngle
TH2F * fhRealCosOpeningAngle
! Cosinus of opening angle of pair version pair energy
Definition: AliAnaPi0.h:327
AliAnaPi0::fhRe1
TH2F ** fhRe1
REAL two-photon invariant mass distribution for different centralities and Asymmetry.
Definition: AliAnaPi0.h:245
AliHistogramRanges::GetHistoMassMin
Float_t GetHistoMassMin() const
Definition: AliHistogramRanges.h:98
AliAnaPi0::fhPrimPi0ProdVertex
TH2F * fhPrimPi0ProdVertex
! Spectrum of primary pi0 vs production vertex
Definition: AliAnaPi0.h:417
AliAnaPi0::fhReDiffPHOSMod
TH2F ** fhReDiffPHOSMod
[fNModules/2]
Definition: AliAnaPi0.h:223
AliAnaPi0::fhMiConv
TH2F * fhMiConv
! MIXED two-photon invariant mass distribution one of the pair was 2 clusters with small mass ...
Definition: AliAnaPi0.h:240
AliHistogramRanges::GetHistoTrackMultiplicityMin
Int_t GetHistoTrackMultiplicityMin() const
Definition: AliHistogramRanges.h:132
AliAnaPi0::fAsymCuts
Float_t fAsymCuts[10]
Array with different assymetry cuts.
Definition: AliAnaPi0.h:173
AliAnaPi0::fPtCuts
Float_t fPtCuts[10]
Array with different pt cuts.
Definition: AliAnaPi0.h:171
AliAnaPi0::fhMCOrgAsym
TH2F * fhMCOrgAsym[13]
! Asymmetry vs pt of real pairs, check common origin of pair
Definition: AliAnaPi0.h:386
AliAnaPi0::fhRePtNCellAsymCuts
TH2F ** fhRePtNCellAsymCuts
[GetNCentrBin()*fNPIDBits*fNAsymCuts]
Definition: AliAnaPi0.h:293
AliAnaPi0::fhPrimPi0PtOrigin
TH2F * fhPrimPi0PtOrigin
! Spectrum of generated pi0 vs mother
Definition: AliAnaPi0.h:376
AliHistogramRanges::GetHistoDiffTimeMin
Float_t GetHistoDiffTimeMin() const
Definition: AliHistogramRanges.h:278
AliCalorimeterUtils::IsMCParticleInCalorimeterAcceptance
Bool_t IsMCParticleInCalorimeterAcceptance(Int_t calo, TParticle *particle)
Definition: AliCalorimeterUtils.cxx:1391
AliAnaPi0::fhMCPi0MassPtRec
TH2F ** fhMCPi0MassPtRec
Real pi0 pairs, reconstructed mass vs reconstructed pt of original pair.
Definition: AliAnaPi0.h:393
AliAnaPi0::fhPrimPi0AccPtCentrality
TH2F * fhPrimPi0AccPtCentrality
! primary pi0 with accepted daughters reconstructed centrality vs pT
Definition: AliAnaPi0.h:350
AliHistogramRanges::GetHistoPhiMin
Float_t GetHistoPhiMin() const
Definition: AliHistogramRanges.h:61
AliAnaPi0::fNAsymCuts
Int_t fNAsymCuts
Number of assymmetry cuts.
Definition: AliAnaPi0.h:172
AliHistogramRanges::GetHistoDiffTimeMax
Float_t GetHistoDiffTimeMax() const
Definition: AliHistogramRanges.h:279
AliAnaPi0::fhMixedCosOpeningAngle
TH2F * fhMixedCosOpeningAngle
! Cosinus of opening angle of pair version pair energy
Definition: AliAnaPi0.h:329
AliAnaPi0::fProdVertex
TVector3 fProdVertex
! Production vertex, temporary array
Definition: AliAnaPi0.h:198
AliAnaPi0::fhReInvPt2
TH2F ** fhReInvPt2
[GetNCentrBin()*fNPIDBits*fNAsymCuts]
Definition: AliAnaPi0.h:276
AliAnaPi0::fhMiDiffPHOSMod
TH2F ** fhMiDiffPHOSMod
[fNModules/2]
Definition: AliAnaPi0.h:235
AliAnaPi0::fhMi2
TH2F ** fhMi2
[GetNCentrBin()*fNPIDBits*fNAsymCuts]
Definition: AliAnaPi0.h:256
AliAnaPi0::GetCreateOutputObjects
TList * GetCreateOutputObjects()
Definition: AliAnaPi0.cxx:234
AliAnaPi0::fhReInvPt3
TH2F ** fhReInvPt3
[GetNCentrBin()*fNPIDBits*fNAsymCuts]
Definition: AliAnaPi0.h:284
AliAnaPi0::fFillOriginHisto
Bool_t fFillOriginHisto
Fill histograms depending on their origin.
Definition: AliAnaPi0.h:188
AliAnaPi0::fhPrimPi0Phi
TH2F * fhPrimPi0Phi
! Azimutal distribution of primary particles vs pT
Definition: AliAnaPi0.h:343
AliNeutralMesonSelection::IsAngleInWindow
Bool_t IsAngleInWindow(Float_t angle, Float_t e) const
Definition: AliNeutralMesonSelection.cxx:194
AliHistogramRanges::GetHistoAsymmetryMax
Float_t GetHistoAsymmetryMax() const
Definition: AliHistogramRanges.h:108
AliHistogramRanges::GetHistoMassMax
Float_t GetHistoMassMax() const
Definition: AliHistogramRanges.h:99
AliAnaCaloTrackCorrBaseClass
Base class for CaloTrackCorr analysis algorithms.
Definition: AliAnaCaloTrackCorrBaseClass.h:58
AliNeutralMesonSelection.h
AliAnaPi0::fFillAngleHisto
Bool_t fFillAngleHisto
Fill histograms with pair opening angle.
Definition: AliAnaPi0.h:186
AliAnaCaloTrackCorrBaseClass::GetCalorimeterString
virtual TString GetCalorimeterString() const
Definition: AliAnaCaloTrackCorrBaseClass.h:169
AliAnaPi0::fFillBadDistHisto
Bool_t fFillBadDistHisto
Do plots for different distances to bad channels.
Definition: AliAnaPi0.h:184
AliAnaCaloTrackCorrBaseClass::IsRealCaloAcceptanceOn
virtual Bool_t IsRealCaloAcceptanceOn() const
Definition: AliAnaCaloTrackCorrBaseClass.h:183
AliAnaCaloTrackCorrBaseClass::GetFiducialCut
virtual AliFiducialCut * GetFiducialCut()
Definition: AliAnaCaloTrackCorrBaseClass.h:324
AliAnaCaloTrackCorrBaseClass::GetInputAODBranch
virtual TClonesArray * GetInputAODBranch() const
Definition: AliAnaCaloTrackCorrBaseClass.h:137
AliAnaPi0::fhMCEtaPtTruePtRec
TH2F ** fhMCEtaPtTruePtRec
[fNPtCuts*fNAsymCuts*fNCellNCuts]
Definition: AliAnaPi0.h:408
AliAnaPi0::fMultiCutAna
Bool_t fMultiCutAna
Do analysis with several or fixed cut.
Definition: AliAnaPi0.h:168
AliCaloTrackReader::GetAODMCParticles
virtual TClonesArray * GetAODMCParticles() const
Definition: AliCaloTrackReader.cxx:699
AliAnaCaloTrackCorrBaseClass::GetHistogramRanges
virtual AliHistogramRanges * GetHistogramRanges()
Definition: AliAnaCaloTrackCorrBaseClass.h:326
AliAnaPi0::fhReMod
TH2F ** fhReMod
REAL two-photon invariant mass distribution for different calorimeter modules.
Definition: AliAnaPi0.h:214
AliCaloTrackReader::ReadStack
Bool_t ReadStack() const
Definition: AliCaloTrackReader.h:595
AliAnaPi0::fhCentralityNoPair
TH1F * fhCentralityNoPair
! Histogram with centrality bins with no pair
Definition: AliAnaPi0.h:320
AliAnaPi0::fFillAsymmetryHisto
Bool_t fFillAsymmetryHisto
Fill histograms with asymmetry vs pt.
Definition: AliAnaPi0.h:187
AliHistogramRanges::GetHistoDiffTimeBins
Int_t GetHistoDiffTimeBins() const
Definition: AliHistogramRanges.h:277
AliAnaPi0::GetEventIndex
Int_t GetEventIndex(AliAODPWG4Particle *part, Double_t *vert)
Definition: AliAnaPi0.cxx:2672
AliAnaPi0::FillAcceptanceHistograms
void FillAcceptanceHistograms()
Fill acceptance histograms if MC data is available.
Definition: AliAnaPi0.cxx:1311
AliAnaPi0::fhReSameSectorEMCALMod
TH2F ** fhReSameSectorEMCALMod
[fNModules-2]
Definition: AliAnaPi0.h:220
AliAnaPi0::fhMCEtaMassPtTrue
TH2F ** fhMCEtaMassPtTrue
[fNPtCuts*fNAsymCuts*fNCellNCuts]
Definition: AliAnaPi0.h:405
AliAnaPi0::fhReInvPt1
TH2F ** fhReInvPt1
[GetNCentrBin()*fNPIDBits*fNAsymCuts]
Definition: AliAnaPi0.h:269
AliAnaPi0::fNCellNCuts
Int_t fNCellNCuts
Number of cuts with number of cells in cluster.
Definition: AliAnaPi0.h:174
AliAnaPi0::fhMiInvPt1
TH2F ** fhMiInvPt1
[GetNCentrBin()*fNPIDBits*fNAsymCuts]
Definition: AliAnaPi0.h:272
AliHistogramRanges::GetHistoAsymmetryBins
Int_t GetHistoAsymmetryBins() const
Definition: AliHistogramRanges.h:106
AliAnaPi0::fhPrimEtaAccYeta
TH2F * fhPrimEtaAccYeta
! PseudoRapidity distribution of primary with accepted daughters vs pT
Definition: AliAnaPi0.h:363
AliMCAnalysisUtils::CheckCommonAncestor
Int_t CheckCommonAncestor(Int_t index1, Int_t index2, const AliCaloTrackReader *reader, Int_t &ancPDG, Int_t &ancStatus, TLorentzVector &momentum, TVector3 &v)
Definition: AliMCAnalysisUtils.cxx:65
AliAnaPi0::fAngleMaxCut
Float_t fAngleMaxCut
Select pairs with opening angle smaller than a threshold.
Definition: AliAnaPi0.h:165
AliAnaPi0::fhEventMixBin
TH1I * fhEventMixBin
! Number of mixed pairs in a particular bin (cen,vz,rp)
Definition: AliAnaPi0.h:318
AliAnaPi0::FillArmenterosThetaStar
void FillArmenterosThetaStar(Int_t pdg)
Fill armenteros plots.
Definition: AliAnaPi0.cxx:1707
AliAnaCaloTrackCorrBaseClass::GetEMCALGeometry
virtual AliEMCALGeometry * GetEMCALGeometry() const
Definition: AliAnaCaloTrackCorrBaseClass.h:336
AliAnaPi0::fhPrimPi0CosOpeningAngle
TH2F * fhPrimPi0CosOpeningAngle
! Cosinus of opening angle of pair version pair energy, pi0 primaries
Definition: AliAnaPi0.h:347
AliAnaPi0::fhMCPi0PtOrigin
TH2F * fhMCPi0PtOrigin
[fNPtCuts*fNAsymCuts*fNCellNCuts]
Definition: AliAnaPi0.h:410
AliAnaPi0::fhMiSameSideEMCALMod
TH2F ** fhMiSameSideEMCALMod
[fNModules]
Definition: AliAnaPi0.h:229
AliAnaPi0::fFillOnlyMCAcceptanceHisto
Bool_t fFillOnlyMCAcceptanceHisto
Do analysis only of MC kinematics input.
Definition: AliAnaPi0.h:190
AliFiducialCut::IsInFiducialCut
Bool_t IsInFiducialCut(Float_t eta, Float_t phi, Int_t det) const
Definition: AliFiducialCut.cxx:78
AliAnaPi0::fMakeInvPtPlots
Bool_t fMakeInvPtPlots
Do plots with inverse pt weight.
Definition: AliAnaPi0.h:180
AliAnaPi0.h
AliAnaPi0::fhPrimEtaPtEventPlane
TH2F * fhPrimEtaPtEventPlane
! primary eta reconstructed event plane vs pT
Definition: AliAnaPi0.h:370
AliAnaPi0::fhPrimPi0AccE
TH1F * fhPrimPi0AccE
! Spectrum of primary with accepted daughters
Definition: AliAnaPi0.h:336
AliAnaPi0::fNPtCuts
Int_t fNPtCuts
Number of pt cuts.
Definition: AliAnaPi0.h:170
AliAnaPi0::fhPrimEtaE
TH1F * fhPrimEtaE
! Spectrum of Primary
Definition: AliAnaPi0.h:355
AliAnaPi0::fhPrimEtaAccPtEventPlane
TH2F * fhPrimEtaAccPtEventPlane
! primary eta with accepted daughters reconstructed event plane vs pT
Definition: AliAnaPi0.h:372
AliAnaCaloTrackCorrBaseClass::GetCaloUtils
virtual AliCalorimeterUtils * GetCaloUtils() const
Definition: AliAnaCaloTrackCorrBaseClass.h:322
AliAnaPi0::fhPrimNotResonancePi0PtOrigin
TH2F * fhPrimNotResonancePi0PtOrigin
! Spectrum of generated pi0 vs mother
Definition: AliAnaPi0.h:378
AliAnaPi0::fhPrimEtaAccY
TH2F * fhPrimEtaAccY
! Rapidity distribution of primary with accepted daughters vs pT
Definition: AliAnaPi0.h:360
AliAnaPi0::fhPrimEtaPt
TH1F * fhPrimEtaPt
! Spectrum of Primary
Definition: AliAnaPi0.h:356
AliAnaPi0::fMultiCutAnaSim
Bool_t fMultiCutAnaSim
Do analysis with several or fixed cut, in the simulation related part.
Definition: AliAnaPi0.h:169
AliAnaPi0::fhPrimEtaAccPt
TH1F * fhPrimEtaAccPt
! Spectrum of primary with accepted daughters
Definition: AliAnaPi0.h:358
AliAnaPi0::fhArmPrimPi0
TH2F * fhArmPrimPi0[4]
! Armenteros plots for primary pi0 in 6 energy bins
Definition: AliAnaPi0.h:424
AliAnaPi0::Print
void Print(const Option_t *opt) const
Print some relevant parameters set for the analysis.
Definition: AliAnaPi0.cxx:1268
AliAnaPi0::fhPrimPi0AccPtEventPlane
TH2F * fhPrimPi0AccPtEventPlane
! primary pi0 with accepted daughters reconstructed event plane vs pT
Definition: AliAnaPi0.h:351
AliAnaCaloTrackCorrBaseClass::GetEventPlane
virtual AliEventplane * GetEventPlane() const
Definition: AliAnaCaloTrackCorrBaseClass.h:111
AliCalorimeterUtils::GetNumberOfSuperModulesUsed
Int_t GetNumberOfSuperModulesUsed() const
Definition: AliCalorimeterUtils.h:350
AliAnaCaloTrackCorrBaseClass::IsHighMultiplicityAnalysisOn
virtual Bool_t IsHighMultiplicityAnalysisOn() const
Definition: AliAnaCaloTrackCorrBaseClass.h:200
AliAnaPi0::fhPrimEtaOpeningAngle
TH2F * fhPrimEtaOpeningAngle
! Opening angle of pair versus pair energy, eta primaries
Definition: AliAnaPi0.h:366
AliAnaPi0::fhMiPtNCellAsymCuts
TH2F ** fhMiPtNCellAsymCuts
[fNPtCuts*fNAsymCuts*fNCellNCuts]
Definition: AliAnaPi0.h:296
AliAnaPi0::fhPrimPi0OpeningAngleAsym
TH2F * fhPrimPi0OpeningAngleAsym
! Opening angle of pair versus pair E asymmetry, pi0 primaries
Definition: AliAnaPi0.h:346
AliAnaPi0::fhPrimPi0AccY
TH2F * fhPrimPi0AccY
! Rapidity distribution of primary with accepted daughters vs pT
Definition: AliAnaPi0.h:339
AliAnaPi0::fhPrimPi0AccPhi
TH2F * fhPrimPi0AccPhi
! Azimutal distribution of primary with accepted daughters vs pT
Definition: AliAnaPi0.h:344
AliAnaPi0::fAngleCut
Float_t fAngleCut
Select pairs with opening angle larger than a threshold.
Definition: AliAnaPi0.h:164
AliAnaCaloTrackCorrBaseClass::GetEventWeight
virtual Double_t GetEventWeight() const
Definition: AliAnaCaloTrackCorrBaseClass.h:243
AliAnaPi0::fhPrimPi0AccYeta
TH2F * fhPrimPi0AccYeta
! PseudoRapidity distribution of primary with accepted daughters vs pT
Definition: AliAnaPi0.h:342
AliAnaCaloTrackCorrBaseClass::GetNCentrBin
virtual Int_t GetNCentrBin() const
Number of bins in reaction plain.
Definition: AliAnaCaloTrackCorrBaseClass.h:232
AliHistogramRanges::GetHistoTrackMultiplicityBins
Int_t GetHistoTrackMultiplicityBins() const
Definition: AliHistogramRanges.h:131
AliAnaPi0::fhMi1
TH2F ** fhMi1
[GetNCentrBin()*fNPIDBits*fNAsymCuts]
Definition: AliAnaPi0.h:248
AliAnaPi0::fhCosThStarPrimEta
TH2F * fhCosThStarPrimEta
! cos(theta*) plots vs E for primary eta, same as asymmetry ...
Definition: AliAnaPi0.h:427
AliAnaPi0::GetAnalysisCuts
TObjString * GetAnalysisCuts()
Save parameters used for analysis.
Definition: AliAnaPi0.cxx:186
AliAnaPi0::MakeAnalysisFillHistograms
void MakeAnalysisFillHistograms()
Definition: AliAnaPi0.cxx:2043
AliHistogramRanges::GetHistoEtaMin
Float_t GetHistoEtaMin() const
Definition: AliHistogramRanges.h:89
AliAnaPi0::fFillSMCombinations
Bool_t fFillSMCombinations
Fill histograms with different cluster pairs in SM combinations.
Definition: AliAnaPi0.h:182
AliAnaPi0::fhReMCFromNotConversion
TH2F * fhReMCFromNotConversion
! Invariant mass of 2 clusters not originated in conversions
Definition: AliAnaPi0.h:421
AliAnaPi0::fhPrimPi0YetaYcut
TH2F * fhPrimPi0YetaYcut
! PseudoRapidity distribution of primary particles vs pT, Y<1
Definition: AliAnaPi0.h:341
AliAnaCaloTrackCorrBaseClass::GetNMaxEvMix
virtual Int_t GetNMaxEvMix() const
Number of bins in track multiplicity.
Definition: AliAnaCaloTrackCorrBaseClass.h:234
AliAnaPi0::fhPrimPi0Yeta
TH2F * fhPrimPi0Yeta
! PseudoRapidity distribution of primary particles vs pT
Definition: AliAnaPi0.h:340
AliAnaPi0::fhRePtMult
TH3F ** fhRePtMult
[fNPIDBits]
Definition: AliAnaPi0.h:305
AliAnaCaloTrackCorrBaseClass::GetModuleNumber
virtual Int_t GetModuleNumber(AliAODPWG4Particle *part) const
Definition: AliAnaCaloTrackCorrBaseClass.h:301
AliAnaPi0::fhMCNotResonancePi0PtOrigin
TH2F * fhMCNotResonancePi0PtOrigin
! Mass of reconstructed pi0 pairs in calorimeter vs mother origin ID, pi0 status 1.
Definition: AliAnaPi0.h:412
AliAnaPi0::fhRePtAsym
TH2F * fhRePtAsym
! REAL two-photon pt vs asymmetry
Definition: AliAnaPi0.h:311
AliHistogramRanges::GetHistoEtaMax
Float_t GetHistoEtaMax() const
Definition: AliHistogramRanges.h:90
AliAnaPi0::fhEPairDiffTime
TH2F * fhEPairDiffTime
! E pair vs Pair of clusters time difference vs E
Definition: AliAnaPi0.h:429
AliAnaPi0::InitParameters
void InitParameters()
Definition: AliAnaPi0.cxx:153
AliAnaPi0::fhMCEtaPtOrigin
TH2F * fhMCEtaPtOrigin
! Mass of reconstructed eta pairs in calorimeter vs mother origin ID.
Definition: AliAnaPi0.h:411
AliAnaPi0::fhMCPi0PtStatus
TH2F * fhMCPi0PtStatus
! Mass of reconstructed pi0 pairs in calorimeter vs mother status.
Definition: AliAnaPi0.h:413
AliHistogramRanges::GetHistoPtBins
Int_t GetHistoPtBins() const
Definition: AliHistogramRanges.h:37
AliAnaPi0::fhMiConv2
TH2F * fhMiConv2
! MIXED two-photon invariant mass distribution both pair photons recombined from 2 clusters with smal...
Definition: AliAnaPi0.h:242
AliAnaPi0::~AliAnaPi0
virtual ~AliAnaPi0()
Destructor.
Definition: AliAnaPi0.cxx:127
AliAnaPi0::fhPrimEtaProdVertex
TH2F * fhPrimEtaProdVertex
! Spectrum of primary eta vs production vertex
Definition: AliAnaPi0.h:418
AliAnaPi0::fSameSM
Bool_t fSameSM
Select only pairs in same SM;.
Definition: AliAnaPi0.h:181
AliAnaPi0::fCheckAccInSector
Bool_t fCheckAccInSector
Check that the decay pi0 falls in the same SM or sector.
Definition: AliAnaPi0.h:192
AliAnaPi0::fhEventBin
TH1I * fhEventBin
! Number of real pairs in a particular bin (cen,vz,rp)
Definition: AliAnaPi0.h:317
AliAnaPi0::fhPrimPi0Pt
TH1F * fhPrimPi0Pt
! Spectrum of Primary
Definition: AliAnaPi0.h:335
AliAnaPi0::fhRe2
TH2F ** fhRe2
[GetNCentrBin()*fNPIDBits*fNAsymCuts]
Definition: AliAnaPi0.h:252
AliAnaPi0::fhPrimEtaAccE
TH1F * fhPrimEtaAccE
! Spectrum of primary with accepted daughters
Definition: AliAnaPi0.h:357
AliAnaPi0::fhMCEtaMassPtRec
TH2F ** fhMCEtaMassPtRec
[fNPtCuts*fNAsymCuts*fNCellNCuts]
Definition: AliAnaPi0.h:402
AliAnaPi0::fUseAngleCut
Bool_t fUseAngleCut
Select pairs depending on their opening angle.
Definition: AliAnaPi0.h:162
AliAnaPi0::fhRePtAsymEta
TH2F * fhRePtAsymEta
! REAL two-photon pt vs asymmetry, close to eta mass
Definition: AliAnaPi0.h:313
AliAnaCaloTrackCorrBaseClass::GetNRPBin
virtual Int_t GetNRPBin() const
Number of bins in vertex.
Definition: AliAnaCaloTrackCorrBaseClass.h:231
AliAnaPi0::fhEventPlaneResolution
TH2F * fhEventPlaneResolution
! Histogram with Event plane resolution vs centrality
Definition: AliAnaPi0.h:322
AliAnaPi0::fhMixedOpeningAngle
TH2F * fhMixedOpeningAngle
! Opening angle of pair versus pair energy
Definition: AliAnaPi0.h:328
AliAnaPi0::FillMCVersusRecDataHistograms
void FillMCVersusRecDataHistograms(Int_t index1, Int_t index2, Float_t pt1, Float_t pt2, Int_t ncells1, Int_t ncells2, Double_t mass, Double_t pt, Double_t asym, Double_t deta, Double_t dphi)
Definition: AliAnaPi0.cxx:1770
AliAnaPi0::fhMCPi0PtTruePtRec
TH2F ** fhMCPi0PtTruePtRec
[fNPtCuts*fNAsymCuts*fNCellNCuts]
Definition: AliAnaPi0.h:399
AliAnaPi0::fhPrimEtaYeta
TH2F * fhPrimEtaYeta
! PseudoRapidity distribution of primary particles vs pT
Definition: AliAnaPi0.h:361
AliAnaCaloTrackCorrBaseClass::GetMCAnalysisUtils
virtual AliMCAnalysisUtils * GetMCAnalysisUtils()
Definition: AliAnaCaloTrackCorrBaseClass.h:330
AliHistogramRanges::GetHistoAsymmetryMin
Float_t GetHistoAsymmetryMin() const
Definition: AliHistogramRanges.h:107
AliAnaPi0::fhCosThStarPrimPi0
TH2F * fhCosThStarPrimPi0
! cos(theta*) plots vs E for primary pi0, same as asymmetry ...
Definition: AliAnaPi0.h:426
AliAnaCaloTrackCorrBaseClass::Print
virtual void Print(const Option_t *) const
Print some relevant parameters set for the analysis.
Definition: AliAnaCaloTrackCorrBaseClass.cxx:611
AliAnaPi0::fhRePIDBits
TH2F ** fhRePIDBits
[fNPtCuts*fNAsymCuts*fNCellNCutsfNModules]
Definition: AliAnaPi0.h:302
AliAnaPi0::fhMCPi0ProdVertex
TH2F * fhMCPi0ProdVertex
! Spectrum of selected pi0 vs production vertex
Definition: AliAnaPi0.h:415
AliAnaPi0::fhPrimPi0PtEventPlane
TH2F * fhPrimPi0PtEventPlane
! primary pi0 reconstructed event plane vs pT
Definition: AliAnaPi0.h:349
AliAnaPi0::fhRePtNCellAsymCutsSM
TH2F ** fhRePtNCellAsymCutsSM[12]
[fNPtCuts*fNAsymCuts*fNCellNCuts]
Definition: AliAnaPi0.h:299
AliAnaPi0::fEventsList
TList ** fEventsList
Containers for photons in stored events.
Definition: AliAnaPi0.h:158
AliAnaPi0::fhMCEtaProdVertex
TH2F * fhMCEtaProdVertex
! Spectrum of selected eta vs production vertex
Definition: AliAnaPi0.h:416
AliAnaPi0::fhPrimEtaOpeningAngleAsym
TH2F * fhPrimEtaOpeningAngleAsym
! Opening angle of pair versus pair E asymmetry, eta primaries
Definition: AliAnaPi0.h:367
AliHistogramRanges::GetHistoPhiMax
Float_t GetHistoPhiMax() const
Definition: AliHistogramRanges.h:62
AliAnaCaloTrackCorrBaseClass::GetReader
virtual AliCaloTrackReader * GetReader() const
Definition: AliAnaCaloTrackCorrBaseClass.h:334
AliAnaPi0::fhMiInvPt3
TH2F ** fhMiInvPt3
[GetNCentrBin()*fNPIDBits*fNAsymCuts]
Definition: AliAnaPi0.h:288
AliAnaPi0::AliAnaPi0
AliAnaPi0()
Default Constructor. Initialized parameters with default values.
Definition: AliAnaPi0.cxx:53
AliHistogramRanges::GetHistoEtaBins
Int_t GetHistoEtaBins() const
Definition: AliHistogramRanges.h:88
AliAnaPi0::fhMi3
TH2F ** fhMi3
[GetNCentrBin()*fNPIDBits*fNAsymCuts]
Definition: AliAnaPi0.h:264
AliAnaPi0::fhArmPrimEta
TH2F * fhArmPrimEta[4]
! Armenteros plots for primary eta in 6 energy bins
Definition: AliAnaPi0.h:425
AliAnaPi0::fFillSSCombinations
Bool_t fFillSSCombinations
Do invariant mass for different combination of shower shape clusters.
Definition: AliAnaPi0.h:185
AliAnaPi0::fhPrimEtaPtOrigin
TH2F * fhPrimEtaPtOrigin
! Spectrum of generated eta vs mother
Definition: AliAnaPi0.h:377
AliAnaPi0::fhMiInvPt2
TH2F ** fhMiInvPt2
[GetNCentrBin()*fNPIDBits*fNAsymCuts]
Definition: AliAnaPi0.h:280
AliMCAnalysisUtils::CheckTagBit
Bool_t CheckTagBit(Int_t tag, UInt_t test) const
Definition: AliMCAnalysisUtils.h:104
AliAnaPi0::fhMCOrgDeltaPhi
TH2F * fhMCOrgDeltaPhi[13]
! Delta Phi vs pt of real pairs, check common origin of pair
Definition: AliAnaPi0.h:388
AliAnaPi0::fhMiMod
TH2F ** fhMiMod
[fNModules]
Definition: AliAnaPi0.h:226
AliAnaCaloTrackCorrBaseClass::GetMixedEvent
virtual AliMixedEvent * GetMixedEvent() const
Definition: AliAnaCaloTrackCorrBaseClass.h:265
AliAnaPi0::fhReConv2
TH2F * fhReConv2
! REAL two-photon invariant mass distribution both pair photons recombined from 2 clusters with small...
Definition: AliAnaPi0.h:241
AliAnaPi0::fhPrimEtaPhi
TH2F * fhPrimEtaPhi
! Azimutal distribution of primary particles vs pT
Definition: AliAnaPi0.h:364
AliAnaPi0::fhMCPi0MassPtTrue
TH2F ** fhMCPi0MassPtTrue
[fNPtCuts*fNAsymCuts*fNCellNCuts]
Definition: AliAnaPi0.h:396