1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
17 //_________________________________________________________________________
18 // Class to check results from simulations or reconstructed real data.
19 // Fill few histograms and do some checking plots
21 //-- Author: Gustavo Conesa (INFN-LNF)
22 //_________________________________________________________________________
25 // --- ROOT system ---
26 //#include "Riostream.h"
27 #include "TObjArray.h"
28 #include "TParticle.h"
29 #include "TDatabasePDG.h"
37 #include <TObjString.h>
39 //---- AliRoot system ----
40 #include "AliAnaCalorimeterQA.h"
41 #include "AliCaloTrackReader.h"
43 #include "AliVCaloCells.h"
44 #include "AliFiducialCut.h"
45 #include "AliAODTrack.h"
46 #include "AliVCluster.h"
47 #include "AliVEvent.h"
48 #include "AliVEventHandler.h"
49 #include "AliAnalysisManager.h"
50 #include "AliAODMCParticle.h"
51 #include "AliMCAnalysisUtils.h"
52 #include "AliAODPid.h"
53 #include "AliExternalTrackParam.h"
55 ClassImp(AliAnaCalorimeterQA)
57 //____________________________________________________________________________
58 AliAnaCalorimeterQA::AliAnaCalorimeterQA() :
59 AliAnaPartCorrBaseClass(), fCalorimeter(""), fStyleMacro(""),
60 fFillAllPosHisto(kFALSE), fFillAllPosHisto2(kTRUE),
61 fFillAllTH12(kFALSE), fFillAllTH3(kTRUE),
62 fFillAllTMHisto(kTRUE), fFillAllPi0Histo(kTRUE),
63 fCorrelate(kTRUE), fNModules(12), fNRCU(2),
64 fTimeCutMin(-1), fTimeCutMax(9999999),
65 fEMCALCellAmpMin(0), fPHOSCellAmpMin(0),
66 fhE(0), fhPt(0), fhPhi(0), fhEta(0), fhEtaPhiE(0),
67 fhECharged(0), fhPtCharged(0), fhPhiCharged(0), fhEtaCharged(0), fhEtaPhiECharged(0),
70 fhIM(0 ), fhIMCellCut(0), fhAsym(0),
71 fhNCellsPerCluster(0), fhNCellsPerClusterNoCut(0), fhNCellsPerClusterMIP(0), fhNCellsPerClusterMIPCharged(0),
72 fhNCellsvsClusterMaxCellDiffE0(0), fhNCellsvsClusterMaxCellDiffE2(0), fhNCellsvsClusterMaxCellDiffE6(0),
76 fhClusterTimeEnergy(0), fhCellTimeSpreadRespectToCellMax(0),
77 fhCellIdCellLargeTimeSpread(0), fhClusterPairDiffTimeE(0),
79 fhClusterMaxCellCloseCellRatio(0), fhClusterMaxCellCloseCellDiff(0),
80 fhClusterMaxCellDiff(0), fhClusterMaxCellDiffNoCut(0),
81 //fhClusterMaxCellDiffDivLambda0(0),
82 fhLambda0vsClusterMaxCellDiffE0(0), fhLambda0vsClusterMaxCellDiffE2(0), fhLambda0vsClusterMaxCellDiffE6(0),
86 fhBadClusterEnergy(0), fhBadClusterTimeEnergy(0), fhBadClusterPairDiffTimeE(0),
87 fhBadClusterMaxCellCloseCellRatio(0), fhBadClusterMaxCellDiff(0),
90 fhRNCells(0), fhXNCells(0), fhYNCells(0), fhZNCells(0),
91 fhRE(0), fhXE(0), fhYE(0), fhZE(0),
93 fhRCellE(0), fhXCellE(0), fhYCellE(0), fhZCellE(0),
95 fhDeltaCellClusterRNCells(0),fhDeltaCellClusterXNCells(0),fhDeltaCellClusterYNCells(0),fhDeltaCellClusterZNCells(0),
96 fhDeltaCellClusterRE(0), fhDeltaCellClusterXE(0), fhDeltaCellClusterYE(0), fhDeltaCellClusterZE(0),
98 fhNCells(0), fhAmplitude(0), fhAmpId(0), fhEtaPhiAmp(0),
99 fhTime(0), fhTimeId(0), fhTimeAmp(0),
100 //fhT0Time(0), fhT0TimeId(0), fhT0TimeAmp(0),
101 fhCaloCorrNClusters(0), fhCaloCorrEClusters(0), fhCaloCorrNCells(0), fhCaloCorrECells(0),
102 fhCaloV0SCorrNClusters(0), fhCaloV0SCorrEClusters(0), fhCaloV0SCorrNCells(0), fhCaloV0SCorrECells(0),
103 fhCaloV0MCorrNClusters(0), fhCaloV0MCorrEClusters(0), fhCaloV0MCorrNCells(0), fhCaloV0MCorrECells(0),
104 fhCaloTrackMCorrNClusters(0), fhCaloTrackMCorrEClusters(0), fhCaloTrackMCorrNCells(0), fhCaloTrackMCorrECells(0),
105 //Super-Module dependent histgrams
106 fhEMod(0), fhNClustersMod(0), fhNCellsPerClusterMod(0), fhNCellsPerClusterModNoCut(0), fhNCellsMod(0),
107 fhGridCellsMod(0), fhGridCellsEMod(0), fhGridCellsTimeMod(0),
108 fhAmplitudeMod(0), fhAmplitudeModFraction(0), fhTimeAmpPerRCU(0),
109 //fhT0TimeAmpPerRCU(0), fhTimeCorrRCU(0),
110 fhIMMod(0), fhIMCellCutMod(0),
113 fhDeltaE(0), fhDeltaPt(0), fhDeltaPhi(0), fhDeltaEta(0),
114 fhRatioE(0), fhRatioPt(0), fhRatioPhi(0), fhRatioEta(0),
115 fh2E(0), fh2Pt(0), fh2Phi(0), fh2Eta(0),
118 fhGenGamPt(0), fhGenGamEta(0), fhGenGamPhi(0),
119 fhGenPi0Pt(0), fhGenPi0Eta(0), fhGenPi0Phi(0),
120 fhGenEtaPt(0), fhGenEtaEta(0), fhGenEtaPhi(0),
121 fhGenOmegaPt(0), fhGenOmegaEta(0), fhGenOmegaPhi(0),
122 fhGenElePt(0), fhGenEleEta(0), fhGenElePhi(0),
123 fhEMVxyz(0), fhEMR(0), fhHaVxyz(0), fhHaR(0),
124 fhGamE(0), fhGamPt(0), fhGamPhi(0), fhGamEta(0),
125 fhGamDeltaE(0), fhGamDeltaPt(0), fhGamDeltaPhi(0), fhGamDeltaEta(0),
126 fhGamRatioE(0), fhGamRatioPt(0), fhGamRatioPhi(0), fhGamRatioEta(0),
127 fhEleE(0), fhElePt(0), fhElePhi(0), fhEleEta(0),
128 fhPi0E(0), fhPi0Pt(0), fhPi0Phi(0), fhPi0Eta(0),
129 fhNeHadE(0), fhNeHadPt(0), fhNeHadPhi(0), fhNeHadEta(0),
130 fhChHadE(0), fhChHadPt(0), fhChHadPhi(0), fhChHadEta(0),
131 fhGamECharged(0), fhGamPtCharged(0), fhGamPhiCharged(0), fhGamEtaCharged(0),
132 fhEleECharged(0), fhElePtCharged(0), fhElePhiCharged(0), fhEleEtaCharged(0),
133 fhPi0ECharged(0), fhPi0PtCharged(0), fhPi0PhiCharged(0), fhPi0EtaCharged(0),
134 fhNeHadECharged(0), fhNeHadPtCharged(0), fhNeHadPhiCharged(0), fhNeHadEtaCharged(0),
135 fhChHadECharged(0), fhChHadPtCharged(0), fhChHadPhiCharged(0), fhChHadEtaCharged(0),
136 fhGenGamAccE(0), fhGenGamAccPt(0), fhGenGamAccEta(0), fhGenGamAccPhi(0),
137 fhGenPi0AccE(0), fhGenPi0AccPt(0), fhGenPi0AccEta(0), fhGenPi0AccPhi(0),
138 fh1pOverE(0), fh1dR(0), fh2EledEdx(0), fh2MatchdEdx(0),
139 fhMCEle1pOverE(0), fhMCEle1dR(0), fhMCEle2MatchdEdx(0),
140 fhMCChHad1pOverE(0), fhMCChHad1dR(0), fhMCChHad2MatchdEdx(0),
141 fhMCNeutral1pOverE(0), fhMCNeutral1dR(0), fhMCNeutral2MatchdEdx(0),fh1pOverER02(0),
142 fhMCEle1pOverER02(0), fhMCChHad1pOverER02(0), fhMCNeutral1pOverER02(0)
146 //Initialize parameters
150 //________________________________________________________________________
151 TObjString * AliAnaCalorimeterQA::GetAnalysisCuts()
153 //Save parameters used for analysis
154 TString parList ; //this will be list of parameters used for this analysis.
155 const Int_t buffersize = 255;
156 char onePar[buffersize] ;
158 snprintf(onePar,buffersize,"--- AliAnaCalorimeterQA ---\n") ;
160 snprintf(onePar,buffersize,"Calorimeter: %s\n",fCalorimeter.Data()) ;
162 snprintf(onePar,buffersize,"Time Cut : %2.2f < T < %2.2f ns \n",fTimeCutMin, fTimeCutMax) ;
164 snprintf(onePar,buffersize,"PHOS Cell Amplitude > %2.2f GeV, EMCAL Cell Amplitude > %2.2f GeV \n",fPHOSCellAmpMin, fEMCALCellAmpMin) ;
166 //Get parameters set in base class.
167 //parList += GetBaseParametersList() ;
169 //Get parameters set in FiducialCut class (not available yet)
170 //parlist += GetFidCut()->GetFidCutParametersList()
172 return new TObjString(parList) ;
176 //________________________________________________________________________
177 TList * AliAnaCalorimeterQA::GetCreateOutputObjects()
179 // Create histograms to be saved in output file and
180 // store them in outputContainer
182 TList * outputContainer = new TList() ;
183 outputContainer->SetName("QAHistos") ;
186 Int_t nptbins = GetHistoPtBins(); Float_t ptmax = GetHistoPtMax(); Float_t ptmin = GetHistoPtMin();
187 Int_t nfineptbins = GetHistoFinePtBins(); Float_t ptfinemax = GetHistoFinePtMax(); Float_t ptfinemin = GetHistoFinePtMin();
188 Int_t nphibins = GetHistoPhiBins(); Float_t phimax = GetHistoPhiMax(); Float_t phimin = GetHistoPhiMin();
189 Int_t netabins = GetHistoEtaBins(); Float_t etamax = GetHistoEtaMax(); Float_t etamin = GetHistoEtaMin();
190 Int_t nmassbins = GetHistoMassBins(); Float_t massmax = GetHistoMassMax(); Float_t massmin = GetHistoMassMin();
191 Int_t nasymbins = GetHistoAsymmetryBins(); Float_t asymmax = GetHistoAsymmetryMax(); Float_t asymmin = GetHistoAsymmetryMin();
192 Int_t nPoverEbins = GetHistoPOverEBins(); Float_t pOverEmax = GetHistoPOverEMax(); Float_t pOverEmin = GetHistoPOverEMin();
193 Int_t ndedxbins = GetHistodEdxBins(); Float_t dedxmax = GetHistodEdxMax(); Float_t dedxmin = GetHistodEdxMin();
194 Int_t ndRbins = GetHistodRBins(); Float_t dRmax = GetHistodRMax(); Float_t dRmin = GetHistodRMin();
195 Int_t ntimebins = GetHistoTimeBins(); Float_t timemax = GetHistoTimeMax(); Float_t timemin = GetHistoTimeMin();
196 Int_t nbins = GetHistoNClusterCellBins(); Int_t nmax = GetHistoNClusterCellMax(); Int_t nmin = GetHistoNClusterCellMin();
197 Int_t nratiobins = GetHistoRatioBins(); Float_t ratiomax = GetHistoRatioMax(); Float_t ratiomin = GetHistoRatioMin();
198 Int_t nvdistbins = GetHistoVertexDistBins(); Float_t vdistmax = GetHistoVertexDistMax(); Float_t vdistmin = GetHistoVertexDistMin();
199 Int_t rbins = GetHistoRBins(); Float_t rmax = GetHistoRMax(); Float_t rmin = GetHistoRMin();
200 Int_t xbins = GetHistoXBins(); Float_t xmax = GetHistoXMax(); Float_t xmin = GetHistoXMin();
201 Int_t ybins = GetHistoYBins(); Float_t ymax = GetHistoYMax(); Float_t ymin = GetHistoYMin();
202 Int_t zbins = GetHistoZBins(); Float_t zmax = GetHistoZMax(); Float_t zmin = GetHistoZMin();
203 Int_t ssbins = GetHistoShowerShapeBins(); Float_t ssmax = GetHistoShowerShapeMax(); Float_t ssmin = GetHistoShowerShapeMin();
204 Int_t tdbins = GetHistoDiffTimeBins() ; Float_t tdmax = GetHistoDiffTimeMax(); Float_t tdmin = GetHistoDiffTimeMin();
206 Int_t nv0sbins = GetHistoV0SignalBins(); Int_t nv0smax = GetHistoV0SignalMax(); Int_t nv0smin = GetHistoV0SignalMin();
207 Int_t nv0mbins = GetHistoV0MultiplicityBins(); Int_t nv0mmax = GetHistoV0MultiplicityMax(); Int_t nv0mmin = GetHistoV0MultiplicityMin();
208 Int_t ntrmbins = GetHistoTrackMultiplicityBins(); Int_t ntrmmax = GetHistoTrackMultiplicityMax(); Int_t ntrmmin = GetHistoTrackMultiplicityMin();
216 if(fCalorimeter=="PHOS"){
223 fhE = new TH1F ("hE","E reconstructed clusters ", nptbins*5,ptmin,ptmax*5);
224 fhE->SetXTitle("E (GeV)");
225 outputContainer->Add(fhE);
228 fhPt = new TH1F ("hPt","p_{T} reconstructed clusters", nptbins,ptmin,ptmax);
229 fhPt->SetXTitle("p_{T} (GeV/c)");
230 outputContainer->Add(fhPt);
232 fhPhi = new TH1F ("hPhi","#phi reconstructed clusters ",nphibins,phimin,phimax);
233 fhPhi->SetXTitle("#phi (rad)");
234 outputContainer->Add(fhPhi);
236 fhEta = new TH1F ("hEta","#eta reconstructed clusters ",netabins,etamin,etamax);
237 fhEta->SetXTitle("#eta ");
238 outputContainer->Add(fhEta);
241 fhEtaPhiE = new TH3F ("hEtaPhiE","#eta vs #phi vs energy, reconstructed clusters",
242 netabins,etamin,etamax,nphibins,phimin,phimax,nptbins,ptmin,ptmax);
243 fhEtaPhiE->SetXTitle("#eta ");
244 fhEtaPhiE->SetYTitle("#phi (rad)");
245 fhEtaPhiE->SetZTitle("E (GeV) ");
246 outputContainer->Add(fhEtaPhiE);
248 fhClusterTimeEnergy = new TH2F ("hClusterTimeEnergy","energy vs TOF, reconstructed clusters",
249 nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
250 fhClusterTimeEnergy->SetXTitle("E (GeV) ");
251 fhClusterTimeEnergy->SetYTitle("TOF (ns)");
252 outputContainer->Add(fhClusterTimeEnergy);
254 fhClusterPairDiffTimeE = new TH2F("hClusterPairDiffTimeE","cluster pair time difference vs E, only good clusters",
255 nptbins,ptmin,ptmax, tdbins,tdmin,tdmax);
256 fhClusterPairDiffTimeE->SetXTitle("E_{cluster} (GeV)");
257 fhClusterPairDiffTimeE->SetYTitle("#Delta t (ns)");
258 outputContainer->Add(fhClusterPairDiffTimeE);
261 fhClusterMaxCellCloseCellRatio = new TH2F ("hClusterMaxCellCloseCellRatio","energy vs ratio of max cell / neighbour cell, reconstructed clusters",
262 nptbins,ptmin,ptmax, 100,0,1.);
263 fhClusterMaxCellCloseCellRatio->SetXTitle("E_{cluster} (GeV) ");
264 fhClusterMaxCellCloseCellRatio->SetYTitle("E_{cell i}/E_{cell max}");
265 outputContainer->Add(fhClusterMaxCellCloseCellRatio);
267 fhClusterMaxCellCloseCellDiff = new TH2F ("hClusterMaxCellCloseCellDiff","energy vs ratio of max cell / neighbour cell, reconstructed clusters",
268 nptbins,ptmin,ptmax, 200,0,10.);
269 fhClusterMaxCellCloseCellDiff->SetXTitle("E_{cluster} (GeV) ");
270 fhClusterMaxCellCloseCellDiff->SetYTitle("E_{cell max}-E_{cell i}");
271 outputContainer->Add(fhClusterMaxCellCloseCellDiff);
273 fhClusterMaxCellDiff = new TH2F ("hClusterMaxCellDiff","energy vs difference of cluster energy - max cell energy / cluster energy, good clusters",
274 nptbins,ptmin,ptmax, 500,0,1.);
275 fhClusterMaxCellDiff->SetXTitle("E_{cluster} (GeV) ");
276 fhClusterMaxCellDiff->SetYTitle("(E_{cluster} - E_{cell max})/ E_{cluster}");
277 outputContainer->Add(fhClusterMaxCellDiff);
279 fhClusterMaxCellDiffNoCut = new TH2F ("hClusterMaxCellDiffNoCut","energy vs difference of cluster energy - max cell energy / cluster energy",
280 nptbins,ptmin,ptmax, 500,0,1.);
281 fhClusterMaxCellDiffNoCut->SetXTitle("E_{cluster} (GeV) ");
282 fhClusterMaxCellDiffNoCut->SetYTitle("(E_{cluster} - E_{cell max})/ E_{cluster}");
283 outputContainer->Add(fhClusterMaxCellDiffNoCut);
285 // fhClusterMaxCellDiffDivLambda0 = new TH2F ("hClusterMaxCellDiffDivLambda0;","",
286 // nptbins,ptmin,ptmax, 500,0,5.);
287 // fhClusterMaxCellDiffDivLambda0->SetXTitle("E_{cluster} (GeV) ");
288 // fhClusterMaxCellDiffDivLambda0->SetYTitle("(E_{cluster} - E_{cell max})/ E_{cluster} / #lambda_{0}");
289 // outputContainer->Add(fhClusterMaxCellDiffDivLambda0);
291 fhLambda0vsClusterMaxCellDiffE0 = new TH2F ("hLambda0vsClusterMaxCellDiffE0","shower shape, #lambda^{2}_{0} vs fraction of energy carried by max cell, E < 2 GeV ",
292 ssbins,ssmin,ssmax,500,0,1.);
293 fhLambda0vsClusterMaxCellDiffE0->SetYTitle("(E_{cluster} - E_{cell max})/ E_{cluster}");
294 fhLambda0vsClusterMaxCellDiffE0->SetXTitle("#lambda^{2}_{0}");
295 outputContainer->Add(fhLambda0vsClusterMaxCellDiffE0);
297 fhLambda0vsClusterMaxCellDiffE2 = new TH2F ("hLambda0vsClusterMaxCellDiffE2","shower shape, #lambda^{2}_{0} vs fraction of energy carried by max cell, 2 < E < 6 GeV ",
298 ssbins,ssmin,ssmax,500,0,1.);
299 fhLambda0vsClusterMaxCellDiffE2->SetYTitle("(E_{cluster} - E_{cell max})/ E_{cluster}");
300 fhLambda0vsClusterMaxCellDiffE2->SetXTitle("#lambda^{2}_{0}");
301 outputContainer->Add(fhLambda0vsClusterMaxCellDiffE2);
303 fhLambda0vsClusterMaxCellDiffE6 = new TH2F ("hLambda0vsClusterMaxCellDiffE6","shower shape, #lambda^{2}_{0} vs fraction of energy carried by max cell, E > 6 ",
304 ssbins,ssmin,ssmax,500,0,1.);
305 fhLambda0vsClusterMaxCellDiffE6->SetYTitle("(E_{cluster} - E_{cell max})/ E_{cluster}");
306 fhLambda0vsClusterMaxCellDiffE6->SetXTitle("#lambda^{2}_{0}");
307 outputContainer->Add(fhLambda0vsClusterMaxCellDiffE6);
309 fhNCellsvsClusterMaxCellDiffE0 = new TH2F ("hNCellsvsClusterMaxCellDiffE0","N cells per cluster vs fraction of energy carried by max cell, E < 2 GeV ",
310 nbins/5,nmin,nmax/5,500,0,1.);
311 fhNCellsvsClusterMaxCellDiffE0->SetYTitle("(E_{cluster} - E_{cell max})/ E_{cluster}");
312 fhNCellsvsClusterMaxCellDiffE0->SetXTitle("N cells per cluster");
313 outputContainer->Add(fhNCellsvsClusterMaxCellDiffE0);
315 fhNCellsvsClusterMaxCellDiffE2 = new TH2F ("hNCellsvsClusterMaxCellDiffE2","N cells per cluster vs fraction of energy carried by max cell, 2 < E < 6 GeV ",
316 nbins/5,nmin,nmax/5,500,0,1.);
317 fhNCellsvsClusterMaxCellDiffE2->SetYTitle("(E_{cluster} - E_{cell max})/ E_{cluster}");
318 fhNCellsvsClusterMaxCellDiffE2->SetXTitle("N cells per cluster");
319 outputContainer->Add(fhNCellsvsClusterMaxCellDiffE2);
321 fhNCellsvsClusterMaxCellDiffE6 = new TH2F ("hNCellsvsClusterMaxCellDiffE6","N cells per cluster vs fraction of energy carried by max cell, E > 6 ",
322 nbins/5,nmin,nmax/5,500,0,1.);
323 fhNCellsvsClusterMaxCellDiffE6->SetYTitle("(E_{cluster} - E_{cell max})/ E_{cluster}");
324 fhNCellsvsClusterMaxCellDiffE6->SetXTitle("N cells per cluster");
325 outputContainer->Add(fhNCellsvsClusterMaxCellDiffE6);
328 if(fCalorimeter=="EMCAL" && !GetCaloUtils()->GetEMCALRecoUtils()->IsRejectExoticCluster()){
330 fhBadClusterEnergy = new TH1F ("hBadClusterEnergy","Bad cluster energy", nptbins,ptmin,ptmax);
331 fhBadClusterEnergy->SetXTitle("E_{cluster} (GeV) ");
332 outputContainer->Add(fhBadClusterEnergy);
334 fhBadClusterMaxCellCloseCellRatio = new TH2F ("hBadClusterMaxCellCloseCell","energy vs ratio of max cell / neighbour cell constributing cell, reconstructed bad clusters",
335 nptbins,ptmin,ptmax, 100,0,1.);
336 fhBadClusterMaxCellCloseCellRatio->SetXTitle("E_{cluster} (GeV) ");
337 fhBadClusterMaxCellCloseCellRatio->SetYTitle("ratio");
338 outputContainer->Add(fhBadClusterMaxCellCloseCellRatio);
340 fhBadClusterMaxCellDiff = new TH2F ("hBadClusterMaxCellDiff","energy vs difference of cluster energy - max cell energy / cluster energy for bad clusters",
341 nptbins,ptmin,ptmax, 500,0,1.);
342 fhBadClusterMaxCellDiff->SetXTitle("E_{cluster} (GeV) ");
343 fhBadClusterMaxCellDiff->SetYTitle("(E_{cluster} - E_{cell max}) / E_{cluster}");
344 outputContainer->Add(fhBadClusterMaxCellDiff);
346 fhBadClusterTimeEnergy = new TH2F ("hBadClusterTimeEnergy","energy vs TOF of reconstructed bad clusters",
347 nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
348 fhBadClusterTimeEnergy->SetXTitle("E_{cluster} (GeV) ");
349 fhBadClusterTimeEnergy->SetYTitle("TOF (ns)");
350 outputContainer->Add(fhBadClusterTimeEnergy);
352 fhBadClusterPairDiffTimeE = new TH2F("hBadClusterPairDiffTimeE","cluster pair time difference (bad - good) vs E from bad cluster",nptbins,ptmin,ptmax, tdbins,tdmin,tdmax);
353 fhBadClusterPairDiffTimeE->SetXTitle("E_{bad cluster} (GeV)");
354 fhBadClusterPairDiffTimeE->SetYTitle("#Delta t (ns)");
355 outputContainer->Add(fhBadClusterPairDiffTimeE);
362 fhECharged = new TH1F ("hECharged","E reconstructed clusters, matched with track", nptbins,ptmin,ptmax);
363 fhECharged->SetXTitle("E (GeV)");
364 outputContainer->Add(fhECharged);
366 fhPtCharged = new TH1F ("hPtCharged","p_{T} reconstructed clusters, matched with track", nptbins,ptmin,ptmax);
367 fhPtCharged->SetXTitle("p_{T} (GeV/c)");
368 outputContainer->Add(fhPtCharged);
370 fhPhiCharged = new TH1F ("hPhiCharged","#phi reconstructed clusters, matched with track",nphibins,phimin,phimax);
371 fhPhiCharged->SetXTitle("#phi (rad)");
372 outputContainer->Add(fhPhiCharged);
374 fhEtaCharged = new TH1F ("hEtaCharged","#eta reconstructed clusters, matched with track",netabins,etamin,etamax);
375 fhEtaCharged->SetXTitle("#eta ");
376 outputContainer->Add(fhEtaCharged);
379 fhEtaPhiECharged = new TH3F ("hEtaPhiECharged","#eta vs #phi, reconstructed clusters, matched with track",
380 netabins,etamin,etamax,nphibins,phimin,phimax,nptbins,ptmin,ptmax);
381 fhEtaPhiECharged->SetXTitle("#eta ");
382 fhEtaPhiECharged->SetYTitle("#phi ");
383 fhEtaPhiECharged->SetZTitle("E (GeV) ");
384 outputContainer->Add(fhEtaPhiECharged);
387 fh1pOverE = new TH2F("h1pOverE","TRACK matches p/E",nptbins,ptmin,ptmax, nPoverEbins,pOverEmin,pOverEmax);
388 fh1pOverE->SetYTitle("p/E");
389 fh1pOverE->SetXTitle("p_{T} (GeV/c)");
390 outputContainer->Add(fh1pOverE);
392 fh1dR = new TH1F("h1dR","TRACK matches dR",ndRbins,dRmin,dRmax);
393 fh1dR->SetXTitle("#Delta R (rad)");
394 outputContainer->Add(fh1dR) ;
396 fh2MatchdEdx = new TH2F("h2MatchdEdx","dE/dx vs. p for all matches",nptbins,ptmin,ptmax,ndedxbins,dedxmin,dedxmax);
397 fh2MatchdEdx->SetXTitle("p (GeV/c)");
398 fh2MatchdEdx->SetYTitle("<dE/dx>");
399 outputContainer->Add(fh2MatchdEdx);
401 fh2EledEdx = new TH2F("h2EledEdx","dE/dx vs. p for electrons",nptbins,ptmin,ptmax,ndedxbins,dedxmin,dedxmax);
402 fh2EledEdx->SetXTitle("p (GeV/c)");
403 fh2EledEdx->SetYTitle("<dE/dx>");
404 outputContainer->Add(fh2EledEdx) ;
406 fh1pOverER02 = new TH2F("h1pOverER02","TRACK matches p/E, all",nptbins,ptmin,ptmax, nPoverEbins,pOverEmin,pOverEmax);
407 fh1pOverER02->SetYTitle("p/E");
408 fh1pOverER02->SetXTitle("p_{T} (GeV/c)");
409 outputContainer->Add(fh1pOverER02);
412 if(fFillAllPi0Histo){
413 fhIM = new TH2F ("hIM","Cluster pairs Invariant mass vs reconstructed pair energy",nptbins,ptmin,ptmax,nmassbins,massmin,massmax);
414 fhIM->SetXTitle("p_{T, cluster pairs} (GeV) ");
415 fhIM->SetYTitle("M_{cluster pairs} (GeV/c^{2})");
416 outputContainer->Add(fhIM);
418 fhIMCellCut = new TH2F ("hIMCellCut","Cluster (n cell > 1) pairs Invariant mass vs reconstructed pair energy",nptbins,ptmin,ptmax,nmassbins,massmin,massmax);
419 fhIMCellCut->SetXTitle("p_{T, cluster pairs} (GeV) ");
420 fhIMCellCut->SetYTitle("M_{cluster pairs} (GeV/c^{2})");
421 outputContainer->Add(fhIMCellCut);
423 fhAsym = new TH2F ("hAssym","Cluster pairs Asymmetry vs reconstructed pair energy",nptbins,ptmin,ptmax,nasymbins,asymmin,asymmax);
424 fhAsym->SetXTitle("p_{T, cluster pairs} (GeV) ");
425 fhAsym->SetYTitle("Asymmetry");
426 outputContainer->Add(fhAsym);
430 fhNCellsPerClusterNoCut = new TH2F ("hNCellsPerClusterNoCut","# cells per cluster vs energy vs #eta, no bad clusters cut",nptbins,ptmin,ptmax, nbins/5,nmin,nmax/5);
431 fhNCellsPerClusterNoCut->SetXTitle("E (GeV)");
432 fhNCellsPerClusterNoCut->SetYTitle("n cells");
433 outputContainer->Add(fhNCellsPerClusterNoCut);
435 fhNCellsPerCluster = new TH2F ("hNCellsPerCluster","# cells per cluster vs energy vs #eta",nptbins,ptmin,ptmax, nbins/5,nmin,nmax/5);
436 fhNCellsPerCluster->SetXTitle("E (GeV)");
437 fhNCellsPerCluster->SetYTitle("n cells");
438 outputContainer->Add(fhNCellsPerCluster);
440 if((fCalorimeter=="EMCAL" && GetReader()->GetEMCALPtMin() < 0.3) ||
441 (fCalorimeter=="PHOS" && GetReader()->GetPHOSPtMin() < 0.3)) {
442 fhNCellsPerClusterMIP = new TH2F ("hNCellsPerClusterMIP","# cells per cluster vs energy vs #eta, smaller bin for MIP search",
444 fhNCellsPerClusterMIP->SetXTitle("E (GeV)");
445 fhNCellsPerClusterMIP->SetYTitle("n cells");
446 outputContainer->Add(fhNCellsPerClusterMIP);
450 fhNCellsPerClusterMIPCharged = new TH2F ("hNCellsPerClusterMIPCharged","# cells per track-matched cluster vs energy vs #eta, smaller bin for MIP search",
452 fhNCellsPerClusterMIPCharged->SetXTitle("E (GeV)");
453 fhNCellsPerClusterMIPCharged->SetYTitle("n cells");
454 outputContainer->Add(fhNCellsPerClusterMIPCharged);
458 fhNClusters = new TH1F ("hNClusters","# clusters", nbins,nmin,nmax);
459 fhNClusters->SetXTitle("number of clusters");
460 outputContainer->Add(fhNClusters);
462 if(fFillAllPosHisto2){
465 fhXYZ = new TH3F ("hXYZ","Cluster: x vs y vs z",xbins,xmin,xmax,ybins,ymin,ymax,zbins,zmin,zmax);
466 fhXYZ->SetXTitle("x (cm)");
467 fhXYZ->SetYTitle("y (cm)");
468 fhXYZ->SetZTitle("z (cm) ");
469 outputContainer->Add(fhXYZ);
472 fhXNCells = new TH2F ("hXNCells","Cluster X position vs N Clusters per Cell",xbins,xmin,xmax,nbins,nmin,nmax);
473 fhXNCells->SetXTitle("x (cm)");
474 fhXNCells->SetYTitle("N cells per cluster");
475 outputContainer->Add(fhXNCells);
477 fhZNCells = new TH2F ("hZNCells","Cluster Z position vs N Clusters per Cell",zbins,zmin,zmax,nbins,nmin,nmax);
478 fhZNCells->SetXTitle("z (cm)");
479 fhZNCells->SetYTitle("N cells per cluster");
480 outputContainer->Add(fhZNCells);
482 fhXE = new TH2F ("hXE","Cluster X position vs cluster energy",xbins,xmin,xmax,nptbins,ptmin,ptmax);
483 fhXE->SetXTitle("x (cm)");
484 fhXE->SetYTitle("E (GeV)");
485 outputContainer->Add(fhXE);
487 fhZE = new TH2F ("hZE","Cluster Z position vs cluster energy",zbins,zmin,zmax,nptbins,ptmin,ptmax);
488 fhZE->SetXTitle("z (cm)");
489 fhZE->SetYTitle("E (GeV)");
490 outputContainer->Add(fhZE);
493 fhRNCells = new TH2F ("hRNCells","Cluster R position vs N Clusters per Cell",rbins,rmin,rmax,nbins,nmin,nmax);
494 fhRNCells->SetXTitle("r = #sqrt{x^{2}+y^{2}} (cm)");
495 fhRNCells->SetYTitle("N cells per cluster");
496 outputContainer->Add(fhRNCells);
499 fhYNCells = new TH2F ("hYNCells","Cluster Y position vs N Clusters per Cell",ybins,ymin,ymax,nbins,nmin,nmax);
500 fhYNCells->SetXTitle("y (cm)");
501 fhYNCells->SetYTitle("N cells per cluster");
502 outputContainer->Add(fhYNCells);
504 fhRE = new TH2F ("hRE","Cluster R position vs cluster energy",rbins,rmin,rmax,nptbins,ptmin,ptmax);
505 fhRE->SetXTitle("r = #sqrt{x^{2}+y^{2}} (cm)");
506 fhRE->SetYTitle("E (GeV)");
507 outputContainer->Add(fhRE);
509 fhYE = new TH2F ("hYE","Cluster Y position vs cluster energy",ybins,ymin,ymax,nptbins,ptmin,ptmax);
510 fhYE->SetXTitle("y (cm)");
511 fhYE->SetYTitle("E (GeV)");
512 outputContainer->Add(fhYE);
514 if(fFillAllPosHisto){
516 fhRCellE = new TH2F ("hRCellE","Cell R position vs cell energy",rbins,rmin,rmax,nptbins,ptmin,ptmax);
517 fhRCellE->SetXTitle("r = #sqrt{x^{2}+y^{2}} (cm)");
518 fhRCellE->SetYTitle("E (GeV)");
519 outputContainer->Add(fhRCellE);
521 fhXCellE = new TH2F ("hXCellE","Cell X position vs cell energy",xbins,xmin,xmax,nptbins,ptmin,ptmax);
522 fhXCellE->SetXTitle("x (cm)");
523 fhXCellE->SetYTitle("E (GeV)");
524 outputContainer->Add(fhXCellE);
526 fhYCellE = new TH2F ("hYCellE","Cell Y position vs cell energy",ybins,ymin,ymax,nptbins,ptmin,ptmax);
527 fhYCellE->SetXTitle("y (cm)");
528 fhYCellE->SetYTitle("E (GeV)");
529 outputContainer->Add(fhYCellE);
531 fhZCellE = new TH2F ("hZCellE","Cell Z position vs cell energy",zbins,zmin,zmax,nptbins,ptmin,ptmax);
532 fhZCellE->SetXTitle("z (cm)");
533 fhZCellE->SetYTitle("E (GeV)");
534 outputContainer->Add(fhZCellE);
536 fhXYZCell = new TH3F ("hXYZCell","Cell : x vs y vs z",xbins,xmin,xmax,ybins,ymin,ymax,zbins,zmin,zmax);
537 fhXYZCell->SetXTitle("x (cm)");
538 fhXYZCell->SetYTitle("y (cm)");
539 fhXYZCell->SetZTitle("z (cm)");
540 outputContainer->Add(fhXYZCell);
543 Float_t dx = TMath::Abs(xmin)+TMath::Abs(xmax);
544 Float_t dy = TMath::Abs(ymin)+TMath::Abs(ymax);
545 Float_t dz = TMath::Abs(zmin)+TMath::Abs(zmax);
546 Float_t dr = TMath::Abs(rmin)+TMath::Abs(rmax);
548 fhDeltaCellClusterRNCells = new TH2F ("hDeltaCellClusterRNCells","Cluster-Cell R position vs N Clusters per Cell",rbins*2,-dr,dr,nbins,nmin,nmax);
549 fhDeltaCellClusterRNCells->SetXTitle("r = #sqrt{x^{2}+y^{2}} (cm)");
550 fhDeltaCellClusterRNCells->SetYTitle("N cells per cluster");
551 outputContainer->Add(fhDeltaCellClusterRNCells);
553 fhDeltaCellClusterXNCells = new TH2F ("hDeltaCellClusterXNCells","Cluster-Cell X position vs N Clusters per Cell",xbins*2,-dx,dx,nbins,nmin,nmax);
554 fhDeltaCellClusterXNCells->SetXTitle("x (cm)");
555 fhDeltaCellClusterXNCells->SetYTitle("N cells per cluster");
556 outputContainer->Add(fhDeltaCellClusterXNCells);
558 fhDeltaCellClusterYNCells = new TH2F ("hDeltaCellClusterYNCells","Cluster-Cell Y position vs N Clusters per Cell",ybins*2,-dy,dy,nbins,nmin,nmax);
559 fhDeltaCellClusterYNCells->SetXTitle("y (cm)");
560 fhDeltaCellClusterYNCells->SetYTitle("N cells per cluster");
561 outputContainer->Add(fhDeltaCellClusterYNCells);
563 fhDeltaCellClusterZNCells = new TH2F ("hDeltaCellClusterZNCells","Cluster-Cell Z position vs N Clusters per Cell",zbins*2,-dz,dz,nbins,nmin,nmax);
564 fhDeltaCellClusterZNCells->SetXTitle("z (cm)");
565 fhDeltaCellClusterZNCells->SetYTitle("N cells per cluster");
566 outputContainer->Add(fhDeltaCellClusterZNCells);
568 fhDeltaCellClusterRE = new TH2F ("hDeltaCellClusterRE","Cluster-Cell R position vs cluster energy",rbins*2,-dr,dr,nptbins,ptmin,ptmax);
569 fhDeltaCellClusterRE->SetXTitle("r = #sqrt{x^{2}+y^{2}} (cm)");
570 fhDeltaCellClusterRE->SetYTitle("E (GeV)");
571 outputContainer->Add(fhDeltaCellClusterRE);
573 fhDeltaCellClusterXE = new TH2F ("hDeltaCellClusterXE","Cluster-Cell X position vs cluster energy",xbins*2,-dx,dx,nptbins,ptmin,ptmax);
574 fhDeltaCellClusterXE->SetXTitle("x (cm)");
575 fhDeltaCellClusterXE->SetYTitle("E (GeV)");
576 outputContainer->Add(fhDeltaCellClusterXE);
578 fhDeltaCellClusterYE = new TH2F ("hDeltaCellClusterYE","Cluster-Cell Y position vs cluster energy",ybins*2,-dy,dy,nptbins,ptmin,ptmax);
579 fhDeltaCellClusterYE->SetXTitle("y (cm)");
580 fhDeltaCellClusterYE->SetYTitle("E (GeV)");
581 outputContainer->Add(fhDeltaCellClusterYE);
583 fhDeltaCellClusterZE = new TH2F ("hDeltaCellClusterZE","Cluster-Cell Z position vs cluster energy",zbins*2,-dz,dz,nptbins,ptmin,ptmax);
584 fhDeltaCellClusterZE->SetXTitle("z (cm)");
585 fhDeltaCellClusterZE->SetYTitle("E (GeV)");
586 outputContainer->Add(fhDeltaCellClusterZE);
588 fhEtaPhiAmp = new TH3F ("hEtaPhiAmp","Cell #eta vs cell #phi vs cell energy",netabins,etamin,etamax,nphibins,phimin,phimax,nptbins,ptmin,ptmax);
589 fhEtaPhiAmp->SetXTitle("#eta ");
590 fhEtaPhiAmp->SetYTitle("#phi (rad)");
591 fhEtaPhiAmp->SetZTitle("E (GeV) ");
592 outputContainer->Add(fhEtaPhiAmp);
597 fhNCells = new TH1F ("hNCells","# cells", colmax*rowmax*fNModules,0,colmax*rowmax*fNModules);
598 fhNCells->SetXTitle("n cells");
599 outputContainer->Add(fhNCells);
601 fhAmplitude = new TH1F ("hAmplitude","Cell Energy", nptbins*2,ptmin,ptmax);
602 fhAmplitude->SetXTitle("Cell Energy (GeV)");
603 outputContainer->Add(fhAmplitude);
605 fhAmpId = new TH2F ("hAmpId","Cell Energy", nfineptbins,ptfinemin,ptfinemax,rowmax*colmax*fNModules,0,rowmax*colmax*fNModules);
606 fhAmpId->SetXTitle("Cell Energy (GeV)");
607 outputContainer->Add(fhAmpId);
609 //Cell Time histograms, time only available in ESDs
610 if(GetReader()->GetDataType()==AliCaloTrackReader::kESD) {
612 fhCellTimeSpreadRespectToCellMax = new TH2F ("hCellTimeSpreadRespectToCellMax","t_{cell max}-t_{cell i} per cluster", nptbins,ptmin,ptmax, 250,-500,500);
613 fhCellTimeSpreadRespectToCellMax->SetXTitle("E (GeV)");
614 fhCellTimeSpreadRespectToCellMax->SetYTitle("#Delta t_{cell max-i} (ns)");
615 outputContainer->Add(fhCellTimeSpreadRespectToCellMax);
617 fhCellIdCellLargeTimeSpread= new TH1F ("hCellIdCellLargeTimeSpread","Cells with time 100 ns larger than cell max in cluster ",
618 colmax*rowmax*fNModules,0,colmax*rowmax*fNModules);
619 fhCellIdCellLargeTimeSpread->SetXTitle("Absolute Cell Id");
620 outputContainer->Add(fhCellIdCellLargeTimeSpread);
622 fhTime = new TH1F ("hTime","Cell Time",ntimebins,timemin,timemax);
623 fhTime->SetXTitle("Cell Time (ns)");
624 outputContainer->Add(fhTime);
626 fhTimeId = new TH2F ("hTimeId","Cell Time vs Absolute Id",ntimebins,timemin,timemax,rowmax*colmax*fNModules,0,rowmax*colmax*fNModules);
627 fhTimeId->SetXTitle("Cell Time (ns)");
628 fhTimeId->SetYTitle("Cell Absolute Id");
629 outputContainer->Add(fhTimeId);
631 fhTimeAmp = new TH2F ("hTimeAmp","Cell Time vs Cell Energy",nptbins*2,ptmin,ptmax,ntimebins,timemin,timemax);
632 fhTimeAmp->SetYTitle("Cell Time (ns)");
633 fhTimeAmp->SetXTitle("Cell Energy (GeV)");
634 outputContainer->Add(fhTimeAmp);
636 // fhT0Time = new TH1F ("hT0Time","Cell Time",ntimebins,timemin,timemax);
637 // fhT0Time->SetXTitle("T_{0} - T_{EMCal} (ns)");
638 // outputContainer->Add(fhT0Time);
640 // fhT0TimeId = new TH2F ("hT0TimeId","Cell Time vs Absolute Id",ntimebins,timemin,timemax,rowmax*colmax*fNModules,0,rowmax*colmax*fNModules);
641 // fhT0TimeId->SetXTitle("T_{0} - T_{EMCal} (ns)");
642 // fhT0TimeId->SetYTitle("Cell Absolute Id");
643 // outputContainer->Add(fhT0TimeId);
645 // fhT0TimeAmp = new TH2F ("hT0TimeAmp","Cell Time vs Cell Energy",nptbins*2,ptmin,ptmax,ntimebins,timemin,timemax);
646 // fhT0TimeAmp->SetYTitle("T_{0} - T_{EMCal} (ns)");
647 // fhT0TimeAmp->SetXTitle("Cell Energy (GeV)");
648 // outputContainer->Add(fhT0TimeAmp);
653 fhCaloCorrNClusters = new TH2F ("hCaloCorrNClusters","# clusters in EMCAL vs PHOS", nbins,nmin,nmax,nbins,nmin,nmax);
654 fhCaloCorrNClusters->SetXTitle("number of clusters in EMCAL");
655 fhCaloCorrNClusters->SetYTitle("number of clusters in PHOS");
656 outputContainer->Add(fhCaloCorrNClusters);
658 fhCaloCorrEClusters = new TH2F ("hCaloCorrEClusters","summed energy of clusters in EMCAL vs PHOS", nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
659 fhCaloCorrEClusters->SetXTitle("#Sigma E of clusters in EMCAL (GeV)");
660 fhCaloCorrEClusters->SetYTitle("#Sigma E of clusters in PHOS (GeV)");
661 outputContainer->Add(fhCaloCorrEClusters);
663 fhCaloCorrNCells = new TH2F ("hCaloCorrNCells","# Cells in EMCAL vs PHOS", nbins,nmin,nmax, nbins,nmin,nmax);
664 fhCaloCorrNCells->SetXTitle("number of Cells in EMCAL");
665 fhCaloCorrNCells->SetYTitle("number of Cells in PHOS");
666 outputContainer->Add(fhCaloCorrNCells);
668 fhCaloCorrECells = new TH2F ("hCaloCorrECells","summed energy of Cells in EMCAL vs PHOS", nptbins*2,ptmin,ptmax*2,nptbins*2,ptmin,ptmax*2);
669 fhCaloCorrECells->SetXTitle("#Sigma E of Cells in EMCAL (GeV)");
670 fhCaloCorrECells->SetYTitle("#Sigma E of Cells in PHOS (GeV)");
671 outputContainer->Add(fhCaloCorrECells);
673 //Calorimeter VS V0 signal
674 fhCaloV0SCorrNClusters = new TH2F ("hCaloV0SNClusters",Form("# clusters in %s vs V0 signal",fCalorimeter.Data()), nv0sbins,nv0smin,nv0smax,nbins,nmin,nmax);
675 fhCaloV0SCorrNClusters->SetXTitle("V0 signal");
676 fhCaloV0SCorrNClusters->SetYTitle(Form("number of clusters in %s",fCalorimeter.Data()));
677 outputContainer->Add(fhCaloV0SCorrNClusters);
679 fhCaloV0SCorrEClusters = new TH2F ("hCaloV0SEClusters",Form("summed energy of clusters in %s vs V0 signal",fCalorimeter.Data()), nv0sbins,nv0smin,nv0smax,nptbins,ptmin,ptmax);
680 fhCaloV0SCorrEClusters->SetXTitle("V0 signal");
681 fhCaloV0SCorrEClusters->SetYTitle(Form("#Sigma E of clusters in %s (GeV)",fCalorimeter.Data()));
682 outputContainer->Add(fhCaloV0SCorrEClusters);
684 fhCaloV0SCorrNCells = new TH2F ("hCaloV0SNCells",Form("# Cells in %s vs V0 signal",fCalorimeter.Data()), nv0sbins,nv0smin,nv0smax, nbins,nmin,nmax);
685 fhCaloV0SCorrNCells->SetXTitle("V0 signal");
686 fhCaloV0SCorrNCells->SetYTitle(Form("number of Cells in %s",fCalorimeter.Data()));
687 outputContainer->Add(fhCaloV0SCorrNCells);
689 fhCaloV0SCorrECells = new TH2F ("hCaloV0SECells",Form("summed energy of Cells in %s vs V0 signal",fCalorimeter.Data()), nv0sbins,nv0smin,nv0smax,nptbins,ptmin,ptmax);
690 fhCaloV0SCorrECells->SetXTitle("V0 signal");
691 fhCaloV0SCorrECells->SetYTitle(Form("#Sigma E of Cells in %s (GeV)",fCalorimeter.Data()));
692 outputContainer->Add(fhCaloV0SCorrECells);
694 //Calorimeter VS V0 multiplicity
695 fhCaloV0MCorrNClusters = new TH2F ("hCaloV0MNClusters",Form("# clusters in %s vs V0 signal",fCalorimeter.Data()), nv0mbins,nv0mmin,nv0mmax,nbins,nmin,nmax);
696 fhCaloV0MCorrNClusters->SetXTitle("V0 signal");
697 fhCaloV0MCorrNClusters->SetYTitle(Form("number of clusters in %s",fCalorimeter.Data()));
698 outputContainer->Add(fhCaloV0MCorrNClusters);
700 fhCaloV0MCorrEClusters = new TH2F ("hCaloV0MEClusters",Form("summed energy of clusters in %s vs V0 signal",fCalorimeter.Data()), nv0mbins,nv0mmin,nv0mmax,nptbins,ptmin,ptmax);
701 fhCaloV0MCorrEClusters->SetXTitle("V0 signal");
702 fhCaloV0MCorrEClusters->SetYTitle(Form("#Sigma E of clusters in %s (GeV)",fCalorimeter.Data()));
703 outputContainer->Add(fhCaloV0MCorrEClusters);
705 fhCaloV0MCorrNCells = new TH2F ("hCaloV0MNCells",Form("# Cells in %s vs V0 signal",fCalorimeter.Data()), nv0mbins,nv0mmin,nv0mmax, nbins,nmin,nmax);
706 fhCaloV0MCorrNCells->SetXTitle("V0 signal");
707 fhCaloV0MCorrNCells->SetYTitle(Form("number of Cells in %s",fCalorimeter.Data()));
708 outputContainer->Add(fhCaloV0MCorrNCells);
710 fhCaloV0MCorrECells = new TH2F ("hCaloV0MECells",Form("summed energy of Cells in %s vs V0 signal",fCalorimeter.Data()), nv0mbins,nv0mmin,nv0mmax,nptbins,ptmin,ptmax);
711 fhCaloV0MCorrECells->SetXTitle("V0 signal");
712 fhCaloV0MCorrECells->SetYTitle(Form("#Sigma E of Cells in %s (GeV)",fCalorimeter.Data()));
713 outputContainer->Add(fhCaloV0MCorrECells);
715 //Calorimeter VS Track multiplicity
716 fhCaloTrackMCorrNClusters = new TH2F ("hCaloTrackMNClusters",Form("# clusters in %s vs # tracks",fCalorimeter.Data()), ntrmbins,ntrmmin,ntrmmax,nbins,nmin,nmax);
717 fhCaloTrackMCorrNClusters->SetXTitle("# tracks");
718 fhCaloTrackMCorrNClusters->SetYTitle(Form("number of clusters in %s",fCalorimeter.Data()));
719 outputContainer->Add(fhCaloTrackMCorrNClusters);
721 fhCaloTrackMCorrEClusters = new TH2F ("hCaloTrackMEClusters",Form("summed energy of clusters in %s vs # tracks",fCalorimeter.Data()), ntrmbins,ntrmmin,ntrmmax,nptbins,ptmin,ptmax);
722 fhCaloTrackMCorrEClusters->SetXTitle("# tracks");
723 fhCaloTrackMCorrEClusters->SetYTitle(Form("#Sigma E of clusters in %s (GeV)",fCalorimeter.Data()));
724 outputContainer->Add(fhCaloTrackMCorrEClusters);
726 fhCaloTrackMCorrNCells = new TH2F ("hCaloTrackMNCells",Form("# Cells in %s vs # tracks",fCalorimeter.Data()), ntrmbins,ntrmmin,ntrmmax, nbins,nmin,nmax);
727 fhCaloTrackMCorrNCells->SetXTitle("# tracks");
728 fhCaloTrackMCorrNCells->SetYTitle(Form("number of Cells in %s",fCalorimeter.Data()));
729 outputContainer->Add(fhCaloTrackMCorrNCells);
731 fhCaloTrackMCorrECells = new TH2F ("hCaloTrackMECells",Form("summed energy of Cells in %s vs # tracks",fCalorimeter.Data()), ntrmbins,ntrmmin,ntrmmax,nptbins,ptmin,ptmax);
732 fhCaloTrackMCorrECells->SetXTitle("# tracks");
733 fhCaloTrackMCorrECells->SetYTitle(Form("#Sigma E of Cells in %s (GeV)",fCalorimeter.Data()));
734 outputContainer->Add(fhCaloTrackMCorrECells);
737 }//correlate calorimeters
740 fhEMod = new TH1F*[fNModules];
741 fhNClustersMod = new TH1F*[fNModules];
742 fhNCellsPerClusterMod = new TH2F*[fNModules];
743 fhNCellsPerClusterModNoCut = new TH2F*[fNModules];
744 fhNCellsMod = new TH1F*[fNModules];
745 fhGridCellsMod = new TH2F*[fNModules];
746 fhGridCellsEMod = new TH2F*[fNModules];
747 fhGridCellsTimeMod = new TH2F*[fNModules];
748 fhAmplitudeMod = new TH1F*[fNModules];
749 if(fCalorimeter=="EMCAL")
750 fhAmplitudeModFraction = new TH1F*[fNModules*3];
752 fhTimeAmpPerRCU = new TH2F*[fNModules*fNRCU];
753 //fhT0TimeAmpPerRCU = new TH2F*[fNModules*fNRCU];
754 //fhTimeCorrRCU = new TH2F*[fNModules*fNRCU*fNModules*fNRCU];
756 fhIMMod = new TH2F*[fNModules];
757 fhIMCellCutMod = new TH2F*[fNModules];
759 for(Int_t imod = 0; imod < fNModules; imod++){
761 fhEMod[imod] = new TH1F (Form("hE_Mod%d",imod),Form("Cluster reconstructed Energy in Module %d ",imod), nptbins,ptmin,ptmax);
762 fhEMod[imod]->SetXTitle("E (GeV)");
763 outputContainer->Add(fhEMod[imod]);
765 fhNClustersMod[imod] = new TH1F (Form("hNClusters_Mod%d",imod),Form("# clusters in Module %d",imod), nbins,nmin,nmax);
766 fhNClustersMod[imod]->SetXTitle("number of clusters");
767 outputContainer->Add(fhNClustersMod[imod]);
769 fhNCellsPerClusterMod[imod] = new TH2F (Form("hNCellsPerCluster_Mod%d",imod),
770 Form("# cells per cluster vs cluster energy in Module %d",imod),
771 nptbins,ptmin,ptmax, nbins,nmin,nmax);
772 fhNCellsPerClusterMod[imod]->SetXTitle("E (GeV)");
773 fhNCellsPerClusterMod[imod]->SetYTitle("n cells");
774 outputContainer->Add(fhNCellsPerClusterMod[imod]);
776 fhNCellsPerClusterModNoCut[imod] = new TH2F (Form("hNCellsPerClusterNoCut_Mod%d",imod),
777 Form("# cells per cluster vs cluster energy in Module %d, no cut",imod),
778 nptbins,ptmin,ptmax, nbins,nmin,nmax);
779 fhNCellsPerClusterModNoCut[imod]->SetXTitle("E (GeV)");
780 fhNCellsPerClusterModNoCut[imod]->SetYTitle("n cells");
781 outputContainer->Add(fhNCellsPerClusterModNoCut[imod]);
784 fhNCellsMod[imod] = new TH1F (Form("hNCells_Mod%d",imod),Form("# cells in Module %d",imod), colmax*rowmax,0,colmax*rowmax);
785 fhNCellsMod[imod]->SetXTitle("n cells");
786 outputContainer->Add(fhNCellsMod[imod]);
787 fhGridCellsMod[imod] = new TH2F (Form("hGridCells_Mod%d",imod),Form("Entries in grid of cells in Module %d",imod),
788 colmax+2,-1.5,colmax+0.5, rowmax+2,-1.5,rowmax+0.5);
789 fhGridCellsMod[imod]->SetYTitle("row (phi direction)");
790 fhGridCellsMod[imod]->SetXTitle("column (eta direction)");
791 outputContainer->Add(fhGridCellsMod[imod]);
793 fhGridCellsEMod[imod] = new TH2F (Form("hGridCellsE_Mod%d",imod),Form("Accumulated energy in grid of cells in Module %d",imod),
794 colmax+2,-1.5,colmax+0.5, rowmax+2,-1.5,rowmax+0.5);
795 fhGridCellsEMod[imod]->SetYTitle("row (phi direction)");
796 fhGridCellsEMod[imod]->SetXTitle("column (eta direction)");
797 outputContainer->Add(fhGridCellsEMod[imod]);
799 fhGridCellsTimeMod[imod] = new TH2F (Form("hGridCellsTime_Mod%d",imod),Form("Accumulated time in grid of cells in Module %d, with E > 0.5 GeV",imod),
800 colmax+2,-1.5,colmax+0.5, rowmax+2,-1.5,rowmax+0.5);
801 fhGridCellsTimeMod[imod]->SetYTitle("row (phi direction)");
802 fhGridCellsTimeMod[imod]->SetXTitle("column (eta direction)");
803 outputContainer->Add(fhGridCellsTimeMod[imod]);
805 fhAmplitudeMod[imod] = new TH1F (Form("hAmplitude_Mod%d",imod),Form("Cell Energy in Module %d",imod), nptbins*2,ptmin,ptmax);
806 fhAmplitudeMod[imod]->SetXTitle("Cell Energy (GeV)");
807 outputContainer->Add(fhAmplitudeMod[imod]);
809 if(fCalorimeter == "EMCAL"){
810 for(Int_t ifrac = 0; ifrac < 3; ifrac++){
811 fhAmplitudeModFraction[imod*3+ifrac] = new TH1F (Form("hAmplitude_Mod%d_Frac%d",imod,ifrac),Form("Cell reconstructed Energy in Module %d, Fraction %d ",imod,ifrac), nptbins,ptmin,ptmax);
812 fhAmplitudeModFraction[imod*3+ifrac]->SetXTitle("E (GeV)");
813 outputContainer->Add(fhAmplitudeModFraction[imod*3+ifrac]);
817 if(GetReader()->GetDataType()==AliCaloTrackReader::kESD) {
819 for(Int_t ircu = 0; ircu < fNRCU; ircu++){
820 fhTimeAmpPerRCU[imod*fNRCU+ircu] = new TH2F (Form("hTimeAmp_Mod%d_RCU%d",imod,ircu),
821 Form("Cell Energy vs Cell Time in Module %d, RCU %d ",imod,ircu),
822 nptbins,ptmin,ptmax,ntimebins,timemin,timemax);
823 fhTimeAmpPerRCU[imod*fNRCU+ircu]->SetXTitle("E (GeV)");
824 fhTimeAmpPerRCU[imod*fNRCU+ircu]->SetYTitle("time (ns)");
825 outputContainer->Add(fhTimeAmpPerRCU[imod*fNRCU+ircu]);
827 // fhT0TimeAmpPerRCU[imod*fNRCU+ircu] = new TH2F (Form("hT0TimeAmp_Mod%d_RCU%d",imod,ircu),
828 // Form("Cell Energy vs T0-Cell Time in Module %d, RCU %d ",imod,ircu),
829 // nptbins,ptmin,ptmax,ntimebins,timemin,timemax);
830 // fhT0TimeAmpPerRCU[imod*fNRCU+ircu]->SetXTitle("E (GeV)");
831 // fhT0TimeAmpPerRCU[imod*fNRCU+ircu]->SetYTitle("T_{0} - T_{EMCal} (ns)");
832 // outputContainer->Add(fhT0TimeAmpPerRCU[imod*fNRCU+ircu]);
835 // for(Int_t imod2 = 0; imod2 < fNModules; imod2++){
836 // for(Int_t ircu2 = 0; ircu2 < fNModules; ircu2++){
837 // Int_t index = (imod2*fNRCU+ircu2)+(fNModules*fNRCU)*(ircu+imod)+fNRCU*fNModules*imod;
838 // fhTimeCorrRCU[index] = new TH2F (Form("hTimeCorrRCU_Mod%d_RCU%d_CompareTo_Mod%d_RCU%d",imod, ircu,imod2, ircu2),
839 // Form("Cell Energy > 0.3, Correlate cell Time in Module %d, RCU %d to Module %d, RCU %d",imod,ircu,imod2, ircu2),
840 // ntimebins,timemin,timemax,ntimebins,timemin,timemax);
841 // fhTimeCorrRCU[index]->SetXTitle("Trigger Cell Time (ns)");
842 // fhTimeCorrRCU[index]->SetYTitle("Cell Time (ns)");
843 // outputContainer->Add(fhTimeCorrRCU[index]);
848 if(fFillAllPi0Histo){
849 fhIMMod[imod] = new TH2F (Form("hIM_Mod%d",imod),
850 Form("Cluster pairs Invariant mass vs reconstructed pair energy in Module %d",imod),
851 nptbins,ptmin,ptmax,nmassbins,massmin,massmax);
852 fhIMMod[imod]->SetXTitle("p_{T, cluster pairs} (GeV) ");
853 fhIMMod[imod]->SetYTitle("M_{cluster pairs} (GeV/c^{2})");
854 outputContainer->Add(fhIMMod[imod]);
856 fhIMCellCutMod[imod] = new TH2F (Form("hIMCellCut_Mod%d",imod),
857 Form("Cluster (n cells > 1) pairs Invariant mass vs reconstructed pair energy in Module %d",imod),
858 nptbins,ptmin,ptmax,nmassbins,massmin,massmax);
859 fhIMCellCutMod[imod]->SetXTitle("p_{T, cluster pairs} (GeV) ");
860 fhIMCellCutMod[imod]->SetYTitle("M_{cluster pairs} (GeV/c^{2})");
861 outputContainer->Add(fhIMCellCutMod[imod]);
866 //Monte Carlo Histograms
869 fhDeltaE = new TH1F ("hDeltaE","MC - Reco E ", nptbins*2,-ptmax,ptmax);
870 fhDeltaE->SetXTitle("#Delta E (GeV)");
871 outputContainer->Add(fhDeltaE);
873 fhDeltaPt = new TH1F ("hDeltaPt","MC - Reco p_{T} ", nptbins*2,-ptmax,ptmax);
874 fhDeltaPt->SetXTitle("#Delta p_{T} (GeV/c)");
875 outputContainer->Add(fhDeltaPt);
877 fhDeltaPhi = new TH1F ("hDeltaPhi","MC - Reco #phi ",nphibins*2,-phimax,phimax);
878 fhDeltaPhi->SetXTitle("#Delta #phi (rad)");
879 outputContainer->Add(fhDeltaPhi);
881 fhDeltaEta = new TH1F ("hDeltaEta","MC- Reco #eta",netabins*2,-etamax,etamax);
882 fhDeltaEta->SetXTitle("#Delta #eta ");
883 outputContainer->Add(fhDeltaEta);
885 fhRatioE = new TH1F ("hRatioE","Reco/MC E ", nratiobins,ratiomin,ratiomax);
886 fhRatioE->SetXTitle("E_{reco}/E_{gen}");
887 outputContainer->Add(fhRatioE);
889 fhRatioPt = new TH1F ("hRatioPt","Reco/MC p_{T} ", nratiobins,ratiomin,ratiomax);
890 fhRatioPt->SetXTitle("p_{T, reco}/p_{T, gen}");
891 outputContainer->Add(fhRatioPt);
893 fhRatioPhi = new TH1F ("hRatioPhi","Reco/MC #phi ",nratiobins,ratiomin,ratiomax);
894 fhRatioPhi->SetXTitle("#phi_{reco}/#phi_{gen}");
895 outputContainer->Add(fhRatioPhi);
897 fhRatioEta = new TH1F ("hRatioEta","Reco/MC #eta",nratiobins,ratiomin,ratiomax);
898 fhRatioEta->SetXTitle("#eta_{reco}/#eta_{gen} ");
899 outputContainer->Add(fhRatioEta);
901 fh2E = new TH2F ("h2E","E distribution, reconstructed vs generated", nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
902 fh2E->SetXTitle("E_{rec} (GeV)");
903 fh2E->SetYTitle("E_{gen} (GeV)");
904 outputContainer->Add(fh2E);
906 fh2Pt = new TH2F ("h2Pt","p_T distribution, reconstructed vs generated", nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
907 fh2Pt->SetXTitle("p_{T,rec} (GeV/c)");
908 fh2Pt->SetYTitle("p_{T,gen} (GeV/c)");
909 outputContainer->Add(fh2Pt);
911 fh2Phi = new TH2F ("h2Phi","#phi distribution, reconstructed vs generated", nphibins,phimin,phimax, nphibins,phimin,phimax);
912 fh2Phi->SetXTitle("#phi_{rec} (rad)");
913 fh2Phi->SetYTitle("#phi_{gen} (rad)");
914 outputContainer->Add(fh2Phi);
916 fh2Eta = new TH2F ("h2Eta","#eta distribution, reconstructed vs generated", netabins,etamin,etamax,netabins,etamin,etamax);
917 fh2Eta->SetXTitle("#eta_{rec} ");
918 fh2Eta->SetYTitle("#eta_{gen} ");
919 outputContainer->Add(fh2Eta);
921 //Fill histos depending on origin of cluster
922 fhGamE = new TH2F ("hGamE","E reconstructed vs E generated from #gamma", nptbins,ptmin,ptmax, nptbins,ptmin,ptmax);
923 fhGamE->SetXTitle("E_{rec} (GeV)");
924 fhGamE->SetXTitle("E_{gen} (GeV)");
925 outputContainer->Add(fhGamE);
927 fhGamPt = new TH2F ("hGamPt","p_{T} reconstructed vs E generated from #gamma", nptbins,ptmin,ptmax, nptbins,ptmin,ptmax);
928 fhGamPt->SetXTitle("p_{T rec} (GeV/c)");
929 fhGamPt->SetYTitle("p_{T gen} (GeV/c)");
930 outputContainer->Add(fhGamPt);
932 fhGamPhi = new TH2F ("hGamPhi","#phi reconstructed vs E generated from #gamma",nphibins,phimin,phimax,nphibins,phimin,phimax);
933 fhGamPhi->SetXTitle("#phi_{rec} (rad)");
934 fhGamPhi->SetYTitle("#phi_{gen} (rad)");
935 outputContainer->Add(fhGamPhi);
937 fhGamEta = new TH2F ("hGamEta","#eta reconstructed vs E generated from #gamma",netabins,etamin,etamax,netabins,etamin,etamax);
938 fhGamEta->SetXTitle("#eta_{rec} ");
939 fhGamEta->SetYTitle("#eta_{gen} ");
940 outputContainer->Add(fhGamEta);
942 fhGamDeltaE = new TH1F ("hGamDeltaE","#gamma MC - Reco E ", nptbins*2,-ptmax,ptmax);
943 fhGamDeltaE->SetXTitle("#Delta E (GeV)");
944 outputContainer->Add(fhGamDeltaE);
946 fhGamDeltaPt = new TH1F ("hGamDeltaPt","#gamma MC - Reco p_{T} ", nptbins*2,-ptmax,ptmax);
947 fhGamDeltaPt->SetXTitle("#Delta p_{T} (GeV/c)");
948 outputContainer->Add(fhGamDeltaPt);
950 fhGamDeltaPhi = new TH1F ("hGamDeltaPhi","#gamma MC - Reco #phi ",nphibins*2,-phimax,phimax);
951 fhGamDeltaPhi->SetXTitle("#Delta #phi (rad)");
952 outputContainer->Add(fhGamDeltaPhi);
954 fhGamDeltaEta = new TH1F ("hGamDeltaEta","#gamma MC- Reco #eta",netabins*2,-etamax,etamax);
955 fhGamDeltaEta->SetXTitle("#Delta #eta ");
956 outputContainer->Add(fhGamDeltaEta);
958 fhGamRatioE = new TH1F ("hGamRatioE","#gamma Reco/MC E ", nratiobins,ratiomin,ratiomax);
959 fhGamRatioE->SetXTitle("E_{reco}/E_{gen}");
960 outputContainer->Add(fhGamRatioE);
962 fhGamRatioPt = new TH1F ("hGamRatioPt","#gamma Reco/MC p_{T} ", nratiobins,ratiomin,ratiomax);
963 fhGamRatioPt->SetXTitle("p_{T, reco}/p_{T, gen}");
964 outputContainer->Add(fhGamRatioPt);
966 fhGamRatioPhi = new TH1F ("hGamRatioPhi","#gamma Reco/MC #phi ",nratiobins,ratiomin,ratiomax);
967 fhGamRatioPhi->SetXTitle("#phi_{reco}/#phi_{gen}");
968 outputContainer->Add(fhGamRatioPhi);
970 fhGamRatioEta = new TH1F ("hGamRatioEta","#gamma Reco/MC #eta",nratiobins,ratiomin,ratiomax);
971 fhGamRatioEta->SetXTitle("#eta_{reco}/#eta_{gen} ");
972 outputContainer->Add(fhGamRatioEta);
974 fhPi0E = new TH2F ("hPi0E","E reconstructed vs E generated from #pi^{0}", nptbins,ptmin,ptmax, nptbins,ptmin,ptmax);
975 fhPi0E->SetXTitle("E_{rec} (GeV)");
976 fhPi0E->SetYTitle("E_{gen} (GeV)");
977 outputContainer->Add(fhPi0E);
979 fhPi0Pt = new TH2F ("hPi0Pt","p_{T} reconstructed vs E generated from #pi^{0}", nptbins,ptmin,ptmax, nptbins,ptmin,ptmax);
980 fhPi0Pt->SetXTitle("p_{T rec} (GeV/c)");
981 fhPi0Pt->SetYTitle("p_{T gen} (GeV/c)");
982 outputContainer->Add(fhPi0Pt);
984 fhPi0Phi = new TH2F ("hPi0Phi","#phi reconstructed vs E generated from #pi^{0}",nphibins,phimin,phimax,nphibins,phimin,phimax);
985 fhPi0Phi->SetXTitle("#phi_{rec} (rad)");
986 fhPi0Phi->SetYTitle("#phi_{gen} (rad)");
987 outputContainer->Add(fhPi0Phi);
989 fhPi0Eta = new TH2F ("hPi0Eta","#eta reconstructed vs E generated from #pi^{0}",netabins,etamin,etamax,netabins,etamin,etamax);
990 fhPi0Eta->SetXTitle("#eta_{rec} ");
991 fhPi0Eta->SetYTitle("#eta_{gen} ");
992 outputContainer->Add(fhPi0Eta);
994 fhEleE = new TH2F ("hEleE","E reconstructed vs E generated from e^{#pm}", nptbins,ptmin,ptmax, nptbins,ptmin,ptmax);
995 fhEleE->SetXTitle("E_{rec} (GeV)");
996 fhEleE->SetXTitle("E_{gen} (GeV)");
997 outputContainer->Add(fhEleE);
999 fhElePt = new TH2F ("hElePt","p_{T} reconstructed vs E generated from e^{#pm}", nptbins,ptmin,ptmax, nptbins,ptmin,ptmax);
1000 fhElePt->SetXTitle("p_{T rec} (GeV/c)");
1001 fhElePt->SetYTitle("p_{T gen} (GeV/c)");
1002 outputContainer->Add(fhElePt);
1004 fhElePhi = new TH2F ("hElePhi","#phi reconstructed vs E generated from e^{#pm}",nphibins,phimin,phimax,nphibins,phimin,phimax);
1005 fhElePhi->SetXTitle("#phi_{rec} (rad)");
1006 fhElePhi->SetYTitle("#phi_{gen} (rad)");
1007 outputContainer->Add(fhElePhi);
1009 fhEleEta = new TH2F ("hEleEta","#eta reconstructed vs E generated from e^{#pm}",netabins,etamin,etamax,netabins,etamin,etamax);
1010 fhEleEta->SetXTitle("#eta_{rec} ");
1011 fhEleEta->SetYTitle("#eta_{gen} ");
1012 outputContainer->Add(fhEleEta);
1014 fhNeHadE = new TH2F ("hNeHadE","E reconstructed vs E generated from neutral hadron", nptbins,ptmin,ptmax, nptbins,ptmin,ptmax);
1015 fhNeHadE->SetXTitle("E_{rec} (GeV)");
1016 fhNeHadE->SetYTitle("E_{gen} (GeV)");
1017 outputContainer->Add(fhNeHadE);
1019 fhNeHadPt = new TH2F ("hNeHadPt","p_{T} reconstructed vs E generated from neutral hadron", nptbins,ptmin,ptmax, nptbins,ptmin,ptmax);
1020 fhNeHadPt->SetXTitle("p_{T rec} (GeV/c)");
1021 fhNeHadPt->SetYTitle("p_{T gen} (GeV/c)");
1022 outputContainer->Add(fhNeHadPt);
1024 fhNeHadPhi = new TH2F ("hNeHadPhi","#phi reconstructed vs E generated from neutral hadron",nphibins,phimin,phimax,nphibins,phimin,phimax);
1025 fhNeHadPhi->SetXTitle("#phi_{rec} (rad)");
1026 fhNeHadPhi->SetYTitle("#phi_{gen} (rad)");
1027 outputContainer->Add(fhNeHadPhi);
1029 fhNeHadEta = new TH2F ("hNeHadEta","#eta reconstructed vs E generated from neutral hadron",netabins,etamin,etamax,netabins,etamin,etamax);
1030 fhNeHadEta->SetXTitle("#eta_{rec} ");
1031 fhNeHadEta->SetYTitle("#eta_{gen} ");
1032 outputContainer->Add(fhNeHadEta);
1034 fhChHadE = new TH2F ("hChHadE","E reconstructed vs E generated from charged hadron", nptbins,ptmin,ptmax, nptbins,ptmin,ptmax);
1035 fhChHadE->SetXTitle("E_{rec} (GeV)");
1036 fhChHadE->SetYTitle("E_{gen} (GeV)");
1037 outputContainer->Add(fhChHadE);
1039 fhChHadPt = new TH2F ("hChHadPt","p_{T} reconstructed vs E generated from charged hadron", nptbins,ptmin,ptmax, nptbins,ptmin,ptmax);
1040 fhChHadPt->SetXTitle("p_{T rec} (GeV/c)");
1041 fhChHadPt->SetYTitle("p_{T gen} (GeV/c)");
1042 outputContainer->Add(fhChHadPt);
1044 fhChHadPhi = new TH2F ("hChHadPhi","#phi reconstructed vs E generated from charged hadron",nphibins,phimin,phimax,nphibins,phimin,phimax);
1045 fhChHadPhi->SetXTitle("#phi_{rec} (rad)");
1046 fhChHadPhi->SetYTitle("#phi_{gen} (rad)");
1047 outputContainer->Add(fhChHadPhi);
1049 fhChHadEta = new TH2F ("hChHadEta","#eta reconstructed vs E generated from charged hadron",netabins,etamin,etamax,netabins,etamin,etamax);
1050 fhChHadEta->SetXTitle("#eta_{rec} ");
1051 fhChHadEta->SetYTitle("#eta_{gen} ");
1052 outputContainer->Add(fhChHadEta);
1056 fhGamECharged = new TH2F ("hGamECharged","E reconstructed vs E generated from #gamma, track matched cluster", nptbins,ptmin,ptmax, nptbins,ptmin,ptmax);
1057 fhGamECharged->SetXTitle("E_{rec} (GeV)");
1058 fhGamECharged->SetXTitle("E_{gen} (GeV)");
1059 outputContainer->Add(fhGamECharged);
1061 fhGamPtCharged = new TH2F ("hGamPtCharged","p_{T} reconstructed vs E generated from #gamma, track matched cluster", nptbins,ptmin,ptmax, nptbins,ptmin,ptmax);
1062 fhGamPtCharged->SetXTitle("p_{T rec} (GeV/c)");
1063 fhGamPtCharged->SetYTitle("p_{T gen} (GeV/c)");
1064 outputContainer->Add(fhGamPtCharged);
1066 fhGamPhiCharged = new TH2F ("hGamPhiCharged","#phi reconstructed vs E generated from #gamma, track matched cluster",nphibins,phimin,phimax,nphibins,phimin,phimax);
1067 fhGamPhiCharged->SetXTitle("#phi_{rec} (rad)");
1068 fhGamPhiCharged->SetYTitle("#phi_{gen} (rad)");
1069 outputContainer->Add(fhGamPhiCharged);
1071 fhGamEtaCharged = new TH2F ("hGamEtaCharged","#eta reconstructed vs E generated from #gamma, track matched cluster",netabins,etamin,etamax,netabins,etamin,etamax);
1072 fhGamEtaCharged->SetXTitle("#eta_{rec} ");
1073 fhGamEtaCharged->SetYTitle("#eta_{gen} ");
1074 outputContainer->Add(fhGamEtaCharged);
1076 fhPi0ECharged = new TH2F ("hPi0ECharged","E reconstructed vs E generated from #pi^{0}, track matched cluster", nptbins,ptmin,ptmax, nptbins,ptmin,ptmax);
1077 fhPi0ECharged->SetXTitle("E_{rec} (GeV)");
1078 fhPi0ECharged->SetYTitle("E_{gen} (GeV)");
1079 outputContainer->Add(fhPi0ECharged);
1081 fhPi0PtCharged = new TH2F ("hPi0PtCharged","p_{T} reconstructed vs E generated from #pi^{0}, track matched cluster", nptbins,ptmin,ptmax, nptbins,ptmin,ptmax);
1082 fhPi0PtCharged->SetXTitle("p_{T rec} (GeV/c)");
1083 fhPi0PtCharged->SetYTitle("p_{T gen} (GeV/c)");
1084 outputContainer->Add(fhPi0PtCharged);
1086 fhPi0PhiCharged = new TH2F ("hPi0PhiCharged","#phi reconstructed vs E generated from #pi^{0}, track matched cluster",nphibins,phimin,phimax,nphibins,phimin,phimax);
1087 fhPi0PhiCharged->SetXTitle("#phi_{rec} (rad)");
1088 fhPi0PhiCharged->SetYTitle("#phi_{gen} (rad)");
1089 outputContainer->Add(fhPi0PhiCharged);
1091 fhPi0EtaCharged = new TH2F ("hPi0EtaCharged","#eta reconstructed vs E generated from #pi^{0}, track matched cluster",netabins,etamin,etamax,netabins,etamin,etamax);
1092 fhPi0EtaCharged->SetXTitle("#eta_{rec} ");
1093 fhPi0EtaCharged->SetYTitle("#eta_{gen} ");
1094 outputContainer->Add(fhPi0EtaCharged);
1096 fhEleECharged = new TH2F ("hEleECharged","E reconstructed vs E generated from e^{#pm}, track matched cluster", nptbins,ptmin,ptmax, nptbins,ptmin,ptmax);
1097 fhEleECharged->SetXTitle("E_{rec} (GeV)");
1098 fhEleECharged->SetXTitle("E_{gen} (GeV)");
1099 outputContainer->Add(fhEleECharged);
1101 fhElePtCharged = new TH2F ("hElePtCharged","p_{T} reconstructed vs E generated from e^{#pm}, track matched cluster", nptbins,ptmin,ptmax, nptbins,ptmin,ptmax);
1102 fhElePtCharged->SetXTitle("p_{T rec} (GeV/c)");
1103 fhElePtCharged->SetYTitle("p_{T gen} (GeV/c)");
1104 outputContainer->Add(fhElePtCharged);
1106 fhElePhiCharged = new TH2F ("hElePhiCharged","#phi reconstructed vs E generated from e^{#pm}, track matched cluster",nphibins,phimin,phimax,nphibins,phimin,phimax);
1107 fhElePhiCharged->SetXTitle("#phi_{rec} (rad)");
1108 fhElePhiCharged->SetYTitle("#phi_{gen} (rad)");
1109 outputContainer->Add(fhElePhiCharged);
1111 fhEleEtaCharged = new TH2F ("hEleEtaCharged","#eta reconstructed vs E generated from e^{#pm}, track matched cluster",netabins,etamin,etamax,netabins,etamin,etamax);
1112 fhEleEtaCharged->SetXTitle("#eta_{rec} ");
1113 fhEleEtaCharged->SetYTitle("#eta_{gen} ");
1114 outputContainer->Add(fhEleEtaCharged);
1116 fhNeHadECharged = new TH2F ("hNeHadECharged","E reconstructed vs E generated from neutral hadron, track matched cluster", nptbins,ptmin,ptmax, nptbins,ptmin,ptmax);
1117 fhNeHadECharged->SetXTitle("E_{rec} (GeV)");
1118 fhNeHadECharged->SetYTitle("E_{gen} (GeV)");
1119 outputContainer->Add(fhNeHadECharged);
1121 fhNeHadPtCharged = new TH2F ("hNeHadPtCharged","p_{T} reconstructed vs E generated from neutral hadron, track matched cluster", nptbins,ptmin,ptmax, nptbins,ptmin,ptmax);
1122 fhNeHadPtCharged->SetXTitle("p_{T rec} (GeV/c)");
1123 fhNeHadPtCharged->SetYTitle("p_{T gen} (GeV/c)");
1124 outputContainer->Add(fhNeHadPtCharged);
1126 fhNeHadPhiCharged = new TH2F ("hNeHadPhiCharged","#phi reconstructed vs E generated from neutral hadron, track matched cluster",nphibins,phimin,phimax,nphibins,phimin,phimax);
1127 fhNeHadPhiCharged->SetXTitle("#phi_{rec} (rad)");
1128 fhNeHadPhiCharged->SetYTitle("#phi_{gen} (rad)");
1129 outputContainer->Add(fhNeHadPhiCharged);
1131 fhNeHadEtaCharged = new TH2F ("hNeHadEtaCharged","#eta reconstructed vs E generated from neutral hadron, track matched cluster",netabins,etamin,etamax,netabins,etamin,etamax);
1132 fhNeHadEtaCharged->SetXTitle("#eta_{rec} ");
1133 fhNeHadEtaCharged->SetYTitle("#eta_{gen} ");
1134 outputContainer->Add(fhNeHadEtaCharged);
1136 fhChHadECharged = new TH2F ("hChHadECharged","E reconstructed vs E generated from charged hadron, track matched cluster", nptbins,ptmin,ptmax, nptbins,ptmin,ptmax);
1137 fhChHadECharged->SetXTitle("E_{rec} (GeV)");
1138 fhChHadECharged->SetYTitle("E_{gen} (GeV)");
1139 outputContainer->Add(fhChHadECharged);
1141 fhChHadPtCharged = new TH2F ("hChHadPtCharged","p_{T} reconstructed vs E generated from charged hadron, track matched cluster", nptbins,ptmin,ptmax, nptbins,ptmin,ptmax);
1142 fhChHadPtCharged->SetXTitle("p_{T rec} (GeV/c)");
1143 fhChHadPtCharged->SetYTitle("p_{T gen} (GeV/c)");
1144 outputContainer->Add(fhChHadPtCharged);
1146 fhChHadPhiCharged = new TH2F ("hChHadPhiCharged","#phi reconstructed vs E generated from charged hadron, track matched cluster",nphibins,phimin,phimax,nphibins,phimin,phimax);
1147 fhChHadPhiCharged->SetXTitle("#phi (rad)");
1148 fhChHadPhiCharged->SetXTitle("#phi_{rec} (rad)");
1149 fhChHadPhiCharged->SetYTitle("#phi_{gen} (rad)");
1150 outputContainer->Add(fhChHadPhiCharged);
1152 fhChHadEtaCharged = new TH2F ("hChHadEtaCharged","#eta reconstructed vs E generated from charged hadron, track matched cluster",netabins,etamin,etamax,netabins,etamin,etamax);
1153 fhChHadEtaCharged->SetXTitle("#eta_{rec} ");
1154 fhChHadEtaCharged->SetYTitle("#eta_{gen} ");
1155 outputContainer->Add(fhChHadEtaCharged);
1157 //Vertex of generated particles
1159 fhEMVxyz = new TH2F ("hEMVxyz","Production vertex of reconstructed ElectroMagnetic particles",nvdistbins,vdistmin,vdistmax,nvdistbins,vdistmin,vdistmax);//,100,0,500);
1160 fhEMVxyz->SetXTitle("v_{x}");
1161 fhEMVxyz->SetYTitle("v_{y}");
1162 //fhEMVxyz->SetZTitle("v_{z}");
1163 outputContainer->Add(fhEMVxyz);
1165 fhHaVxyz = new TH2F ("hHaVxyz","Production vertex of reconstructed hadrons",nvdistbins,vdistmin,vdistmax,nvdistbins,vdistmin,vdistmax);//,100,0,500);
1166 fhHaVxyz->SetXTitle("v_{x}");
1167 fhHaVxyz->SetYTitle("v_{y}");
1168 //fhHaVxyz->SetZTitle("v_{z}");
1169 outputContainer->Add(fhHaVxyz);
1171 fhEMR = new TH2F ("hEMR","Distance to production vertex of reconstructed ElectroMagnetic particles vs E rec",nptbins,ptmin,ptmax,nvdistbins,vdistmin,vdistmax);
1172 fhEMR->SetXTitle("E (GeV)");
1173 fhEMR->SetYTitle("TMath::Sqrt(v_{x}^{2}+v_{y}^{2})");
1174 outputContainer->Add(fhEMR);
1176 fhHaR = new TH2F ("hHaR","Distance to production vertex of reconstructed Hadrons vs E rec",nptbins,ptmin,ptmax,nvdistbins,vdistmin,vdistmax);
1177 fhHaR->SetXTitle("E (GeV)");
1178 fhHaR->SetYTitle("TMath::Sqrt(v_{x}^{2}+v_{y}^{2})");
1179 outputContainer->Add(fhHaR);
1184 fhGenGamPt = new TH1F("hGenGamPt" ,"p_{T} of generated #gamma",nptbins,ptmin,ptmax);
1185 fhGenGamEta = new TH1F("hGenGamEta","Y of generated #gamma",netabins,etamin,etamax);
1186 fhGenGamPhi = new TH1F("hGenGamPhi","#phi of generated #gamma",nphibins,phimin,phimax);
1188 fhGenPi0Pt = new TH1F("hGenPi0Pt" ,"p_{T} of generated #pi^{0}",nptbins,ptmin,ptmax);
1189 fhGenPi0Eta = new TH1F("hGenPi0Eta","Y of generated #pi^{0}",netabins,etamin,etamax);
1190 fhGenPi0Phi = new TH1F("hGenPi0Phi","#phi of generated #pi^{0}",nphibins,phimin,phimax);
1192 fhGenEtaPt = new TH1F("hGenEtaPt" ,"p_{T} of generated #eta",nptbins,ptmin,ptmax);
1193 fhGenEtaEta = new TH1F("hGenEtaEta","Y of generated #eta",netabins,etamin,etamax);
1194 fhGenEtaPhi = new TH1F("hGenEtaPhi","#phi of generated #eta",nphibins,phimin,phimax);
1196 fhGenOmegaPt = new TH1F("hGenOmegaPt" ,"p_{T} of generated #omega",nptbins,ptmin,ptmax);
1197 fhGenOmegaEta = new TH1F("hGenOmegaEta","Y of generated #omega",netabins,etamin,etamax);
1198 fhGenOmegaPhi = new TH1F("hGenOmegaPhi","#phi of generated #omega",nphibins,phimin,phimax);
1200 fhGenElePt = new TH1F("hGenElePt" ,"p_{T} of generated e^{#pm}",nptbins,ptmin,ptmax);
1201 fhGenEleEta = new TH1F("hGenEleEta","Y of generated e^{#pm}",netabins,etamin,etamax);
1202 fhGenElePhi = new TH1F("hGenElePhi","#phi of generated e^{#pm}",nphibins,phimin,phimax);
1204 fhGenGamPt->SetXTitle("p_{T} (GeV/c)");
1205 fhGenGamEta->SetXTitle("#eta");
1206 fhGenGamPhi->SetXTitle("#phi (rad)");
1207 outputContainer->Add(fhGenGamPt);
1208 outputContainer->Add(fhGenGamEta);
1209 outputContainer->Add(fhGenGamPhi);
1211 fhGenPi0Pt->SetXTitle("p_{T} (GeV/c)");
1212 fhGenPi0Eta->SetXTitle("#eta");
1213 fhGenPi0Phi->SetXTitle("#phi (rad)");
1214 outputContainer->Add(fhGenPi0Pt);
1215 outputContainer->Add(fhGenPi0Eta);
1216 outputContainer->Add(fhGenPi0Phi);
1218 fhGenEtaPt->SetXTitle("p_{T} (GeV/c)");
1219 fhGenEtaEta->SetXTitle("#eta");
1220 fhGenEtaPhi->SetXTitle("#phi (rad)");
1221 outputContainer->Add(fhGenEtaPt);
1222 outputContainer->Add(fhGenEtaEta);
1223 outputContainer->Add(fhGenEtaPhi);
1225 fhGenOmegaPt->SetXTitle("p_{T} (GeV/c)");
1226 fhGenOmegaEta->SetXTitle("#eta");
1227 fhGenOmegaPhi->SetXTitle("#phi (rad)");
1228 outputContainer->Add(fhGenOmegaPt);
1229 outputContainer->Add(fhGenOmegaEta);
1230 outputContainer->Add(fhGenOmegaPhi);
1232 fhGenElePt->SetXTitle("p_{T} (GeV/c)");
1233 fhGenEleEta->SetXTitle("#eta");
1234 fhGenElePhi->SetXTitle("#phi (rad)");
1235 outputContainer->Add(fhGenElePt);
1236 outputContainer->Add(fhGenEleEta);
1237 outputContainer->Add(fhGenElePhi);
1239 fhGenGamAccE = new TH1F("hGenGamAccE" ,"E of generated #gamma in calorimeter acceptance",nptbins,ptmin,ptmax);
1240 fhGenGamAccPt = new TH1F("hGenGamAccPt" ,"p_{T} of generated #gamma in calorimeter acceptance",nptbins,ptmin,ptmax);
1241 fhGenGamAccEta = new TH1F("hGenGamAccEta","Y of generated #gamma in calorimeter acceptance",netabins,etamin,etamax);
1242 fhGenGamAccPhi = new TH1F("hGenGamAccPhi","#phi of generated #gamma in calorimeter acceptance",nphibins,phimin,phimax);
1244 fhGenPi0AccE = new TH1F("hGenPi0AccE" ,"E of generated #pi^{0} in calorimeter acceptance",nptbins,ptmin,ptmax);
1245 fhGenPi0AccPt = new TH1F("hGenPi0AccPt" ,"p_{T} of generated #pi^{0} in calorimeter acceptance",nptbins,ptmin,ptmax);
1246 fhGenPi0AccEta = new TH1F("hGenPi0AccEta","Y of generated #pi^{0} in calorimeter acceptance",netabins,etamin,etamax);
1247 fhGenPi0AccPhi = new TH1F("hGenPi0AccPhi","#phi of generated #pi^{0} in calorimeter acceptance",nphibins,phimin,phimax);
1249 fhGenGamAccE ->SetXTitle("E (GeV)");
1250 fhGenGamAccPt ->SetXTitle("p_{T} (GeV/c)");
1251 fhGenGamAccEta->SetXTitle("#eta");
1252 fhGenGamAccPhi->SetXTitle("#phi (rad)");
1253 outputContainer->Add(fhGenGamAccE);
1254 outputContainer->Add(fhGenGamAccPt);
1255 outputContainer->Add(fhGenGamAccEta);
1256 outputContainer->Add(fhGenGamAccPhi);
1258 fhGenPi0AccE ->SetXTitle("E (GeV)");
1259 fhGenPi0AccPt ->SetXTitle("p_{T} (GeV/c)");
1260 fhGenPi0AccEta->SetXTitle("#eta");
1261 fhGenPi0AccPhi->SetXTitle("#phi (rad)");
1262 outputContainer->Add(fhGenPi0AccE);
1263 outputContainer->Add(fhGenPi0AccPt);
1264 outputContainer->Add(fhGenPi0AccEta);
1265 outputContainer->Add(fhGenPi0AccPhi);
1269 fhMCEle1pOverE = new TH2F("hMCEle1pOverE","TRACK matches p/E, MC electrons",nptbins,ptmin,ptmax, nPoverEbins,pOverEmin,pOverEmax);
1270 fhMCEle1pOverE->SetYTitle("p/E");
1271 fhMCEle1pOverE->SetXTitle("p_{T} (GeV/c)");
1272 outputContainer->Add(fhMCEle1pOverE);
1274 fhMCEle1dR = new TH1F("hMCEle1dR","TRACK matches dR, MC electrons",ndRbins,dRmin,dRmax);
1275 fhMCEle1dR->SetXTitle("#Delta R (rad)");
1276 outputContainer->Add(fhMCEle1dR) ;
1278 fhMCEle2MatchdEdx = new TH2F("hMCEle2MatchdEdx","dE/dx vs. p for all matches, MC electrons",nptbins,ptmin,ptmax,ndedxbins,dedxmin,dedxmax);
1279 fhMCEle2MatchdEdx->SetXTitle("p (GeV/c)");
1280 fhMCEle2MatchdEdx->SetYTitle("<dE/dx>");
1281 outputContainer->Add(fhMCEle2MatchdEdx);
1283 fhMCChHad1pOverE = new TH2F("hMCChHad1pOverE","TRACK matches p/E, MC charged hadrons",nptbins,ptmin,ptmax, nPoverEbins,pOverEmin,pOverEmax);
1284 fhMCChHad1pOverE->SetYTitle("p/E");
1285 fhMCChHad1pOverE->SetXTitle("p_{T} (GeV/c)");
1286 outputContainer->Add(fhMCChHad1pOverE);
1288 fhMCChHad1dR = new TH1F("hMCChHad1dR","TRACK matches dR, MC charged hadrons",ndRbins,dRmin,dRmax);
1289 fhMCChHad1dR->SetXTitle("#Delta R (rad)");
1290 outputContainer->Add(fhMCChHad1dR) ;
1292 fhMCChHad2MatchdEdx = new TH2F("hMCChHad2MatchdEdx","dE/dx vs. p for all matches, MC charged hadrons",nptbins,ptmin,ptmax,ndedxbins,dedxmin,dedxmax);
1293 fhMCChHad2MatchdEdx->SetXTitle("p (GeV/c)");
1294 fhMCChHad2MatchdEdx->SetYTitle("<dE/dx>");
1295 outputContainer->Add(fhMCChHad2MatchdEdx);
1297 fhMCNeutral1pOverE = new TH2F("hMCNeutral1pOverE","TRACK matches p/E, MC neutrals",nptbins,ptmin,ptmax, nPoverEbins,pOverEmin,pOverEmax);
1298 fhMCNeutral1pOverE->SetYTitle("p/E");
1299 fhMCNeutral1pOverE->SetXTitle("p_{T} (GeV/c)");
1300 outputContainer->Add(fhMCNeutral1pOverE);
1302 fhMCNeutral1dR = new TH1F("hMCNeutral1dR","TRACK matches dR, MC neutrals",ndRbins,dRmin,dRmax);
1303 fhMCNeutral1dR->SetXTitle("#Delta R (rad)");
1304 outputContainer->Add(fhMCNeutral1dR) ;
1306 fhMCNeutral2MatchdEdx = new TH2F("hMCNeutral2MatchdEdx","dE/dx vs. p for all matches, MC neutrals",nptbins,ptmin,ptmax,ndedxbins,dedxmin,dedxmax);
1307 fhMCNeutral2MatchdEdx->SetXTitle("p (GeV/c)");
1308 fhMCNeutral2MatchdEdx->SetYTitle("<dE/dx>");
1309 outputContainer->Add(fhMCNeutral2MatchdEdx);
1311 fhMCEle1pOverER02 = new TH2F("hMCEle1pOverER02","TRACK matches p/E, MC electrons",nptbins,ptmin,ptmax, nPoverEbins,pOverEmin,pOverEmax);
1312 fhMCEle1pOverER02->SetYTitle("p/E");
1313 fhMCEle1pOverER02->SetXTitle("p_{T} (GeV/c)");
1314 outputContainer->Add(fhMCEle1pOverER02);
1316 fhMCChHad1pOverER02 = new TH2F("hMCChHad1pOverER02","TRACK matches p/E, MC charged hadrons",nptbins,ptmin,ptmax, nPoverEbins,pOverEmin,pOverEmax);
1317 fhMCChHad1pOverER02->SetYTitle("p/E");
1318 fhMCChHad1pOverER02->SetXTitle("p_{T} (GeV/c)");
1319 outputContainer->Add(fhMCChHad1pOverER02);
1321 fhMCNeutral1pOverER02 = new TH2F("hMCNeutral1pOverER02","TRACK matches p/E, MC neutrals",nptbins,ptmin,ptmax, nPoverEbins,pOverEmin,pOverEmax);
1322 fhMCNeutral1pOverER02->SetYTitle("p/E");
1323 fhMCNeutral1pOverER02->SetXTitle("p_{T} (GeV/c)");
1324 outputContainer->Add(fhMCNeutral1pOverER02);
1327 // for(Int_t i = 0; i < outputContainer->GetEntries() ; i++)
1328 // printf("i=%d, name= %s\n",i,outputContainer->At(i)->GetName());
1330 return outputContainer;
1333 //_______________________________________________________________________________________________________________________________________
1334 Int_t AliAnaCalorimeterQA::GetNewRebinForRePlotting(TH1D* histo, const Float_t newXmin, const Float_t newXmax,const Int_t newXnbins) const
1336 //Calculate the rebinning for the new requested bin size, only used when replotting executing the Terminte
1338 Float_t oldbinsize = histo->GetBinWidth(0);
1339 Float_t newbinsize = TMath::Abs(newXmax-newXmin) / newXnbins;
1341 //printf("bin size, old %f, new %f\n",oldbinsize,newbinsize);
1342 if(newbinsize > oldbinsize) return (Int_t) (newbinsize/oldbinsize);
1347 //__________________________________________________
1348 void AliAnaCalorimeterQA::Init()
1350 //Check if the data or settings are ok
1352 if(fCalorimeter != "PHOS" && fCalorimeter !="EMCAL")
1353 AliFatal(Form("Wrong calorimeter name <%s>", fCalorimeter.Data()));
1355 if(GetReader()->GetDataType()== AliCaloTrackReader::kMC)
1356 AliFatal("Analysis of reconstructed data, MC reader not aplicable");
1361 //__________________________________________________
1362 void AliAnaCalorimeterQA::InitParameters()
1364 //Initialize the parameters of the analysis.
1365 AddToHistogramsName("AnaCaloQA_");
1367 fCalorimeter = "EMCAL"; //or PHOS
1369 fNModules = 12; // set maximum to maximum number of EMCAL modules
1370 fNRCU = 2; // set maximum number of RCU in EMCAL per SM
1372 fTimeCutMax = 9999999;
1373 fEMCALCellAmpMin = 0.0;
1374 fPHOSCellAmpMin = 0.0;
1378 //__________________________________________________________________
1379 void AliAnaCalorimeterQA::Print(const Option_t * opt) const
1381 //Print some relevant parameters set for the analysis
1385 printf("**** Print %s %s ****\n", GetName(), GetTitle() ) ;
1386 AliAnaPartCorrBaseClass::Print(" ");
1388 printf("Select Calorimeter %s \n",fCalorimeter.Data());
1389 printf("Plots style macro %s \n",fStyleMacro.Data());
1390 printf("Time Cut: %3.1f < TOF < %3.1f\n", fTimeCutMin, fTimeCutMax);
1391 printf("EMCAL Min Amplitude : %2.1f GeV/c\n", fEMCALCellAmpMin) ;
1392 printf("PHOS Min Amplitude : %2.1f GeV/c\n", fPHOSCellAmpMin) ;
1396 //__________________________________________________________________
1397 void AliAnaCalorimeterQA::MakeAnalysisFillHistograms()
1399 //Fill Calorimeter QA histograms
1400 TLorentzVector mom ;
1401 TLorentzVector mom2 ;
1402 TObjArray * caloClusters = NULL;
1405 Int_t nCaloClusters = 0;
1406 Int_t nCaloClustersAccepted = 0;
1407 Int_t nCaloCellsPerCluster = 0;
1408 Int_t nTracksMatched = 0;
1409 Int_t trackIndex = 0;
1412 //Get vertex for photon momentum calculation and event selection
1413 Double_t v[3] = {0,0,0}; //vertex ;
1414 GetReader()->GetVertex(v);
1415 if (TMath::Abs(v[2]) > GetZvertexCut()) return ;
1417 //Play with the MC stack if available
1418 //Get the MC arrays and do some checks
1420 if(GetReader()->ReadStack()){
1423 AliFatal("Stack not available, is the MC handler called?\n");
1425 //Fill some pure MC histograms, only primaries.
1426 for(Int_t i=0 ; i<GetMCStack()->GetNprimary(); i++){//Only primary particles, for all MC transport put GetNtrack()
1427 TParticle *primary = GetMCStack()->Particle(i) ;
1428 //printf("i %d, %s: status = %d, primary? %d\n",i, primary->GetName(), primary->GetStatusCode(), primary->IsPrimary());
1429 if (primary->GetStatusCode() > 11) continue; //Working for PYTHIA and simple generators, check for HERWIG
1430 primary->Momentum(mom);
1431 MCHistograms(mom,TMath::Abs(primary->GetPdgCode()));
1434 else if(GetReader()->ReadAODMCParticles()){
1436 if(!GetReader()->GetAODMCParticles(0))
1437 AliFatal("AODMCParticles not available!");
1439 //Fill some pure MC histograms, only primaries.
1440 for(Int_t i=0 ; i < (GetReader()->GetAODMCParticles(0))->GetEntriesFast(); i++){
1441 AliAODMCParticle *aodprimary = (AliAODMCParticle*) (GetReader()->GetAODMCParticles(0))->At(i) ;
1442 //printf("i %d, %s: primary? %d physical primary? %d, flag %d\n",
1443 // i,(TDatabasePDG::Instance()->GetParticle(aodprimary->GetPdgCode()))->GetName(),
1444 // aodprimary->IsPrimary(), aodprimary->IsPhysicalPrimary(), aodprimary->GetFlag());
1445 if (!aodprimary->IsPrimary()) continue; //accept all which is not MC transport generated. Don't know how to avoid partons
1446 //aodprimary->Momentum(mom);
1447 mom.SetPxPyPzE(aodprimary->Px(),aodprimary->Py(),aodprimary->Pz(),aodprimary->E());
1448 MCHistograms(mom,TMath::Abs(aodprimary->GetPdgCode()));
1455 //Get List with CaloClusters
1456 if (fCalorimeter == "PHOS") caloClusters = GetPHOSClusters();
1457 else if (fCalorimeter == "EMCAL") caloClusters = GetEMCALClusters();
1459 AliFatal(Form("AliAnaCalorimeterQA::MakeAnalysisFillHistograms() - Wrong calorimeter name <%s>, END\n", fCalorimeter.Data()));
1461 // if (fCalorimeter == "EMCAL") GetReader()->GetInputEvent()->GetEMCALClusters(caloClusters);//GetEMCALClusters();
1462 // else if(fCalorimeter == "PHOS") GetReader()->GetInputEvent()->GetPHOSClusters (caloClusters);//GetPHOSClusters();
1464 // AliFatal(Form("AliAnaCalorimeterQA::MakeAnalysisFillHistograms() - Wrong calorimeter name <%s>, END\n", fCalorimeter.Data()));
1467 AliFatal(Form("AliAnaCalorimeterQA::MakeAnalysisFillHistograms() - No CaloClusters available\n"));
1470 //----------------------------------------------------------
1471 //Correlate Calorimeters and V0 and track Multiplicity
1472 //----------------------------------------------------------
1473 if(fCorrelate) Correlate();
1475 //----------------------------------------------------------
1477 //----------------------------------------------------------
1479 nCaloClusters = caloClusters->GetEntriesFast() ;
1480 Int_t *nClustersInModule = new Int_t[fNModules];
1481 for(Int_t imod = 0; imod < fNModules; imod++ ) nClustersInModule[imod] = 0;
1484 printf("AliAnaCalorimeterQA::MakeAnalysisFillHistograms() - In %s there are %d clusters \n", fCalorimeter.Data(), nCaloClusters);
1486 AliVTrack * track = 0x0;
1489 //Loop over CaloClusters
1490 //if(nCaloClusters > 0)printf("QA : Vertex Cut passed %f, cut %f, entries %d, %s\n",v[2], 40., nCaloClusters, fCalorimeter.Data());
1491 for(Int_t iclus = 0; iclus < nCaloClusters; iclus++){
1493 if(GetDebug() > 0) printf("AliAnaCalorimeterQA::MakeAnalysisFillHistograms() - cluster: %d/%d, data %d \n",
1494 iclus+1,nCaloClusters,GetReader()->GetDataType());
1496 AliVCluster* clus = (AliVCluster*)caloClusters->At(iclus);
1497 AliVCaloCells * cell = 0x0;
1498 if(fCalorimeter == "PHOS") cell = GetPHOSCells();
1499 else cell = GetEMCALCells();
1501 //Get cluster kinematics
1502 clus->GetPosition(pos);
1503 clus->GetMomentum(mom,v);
1504 tof = clus->GetTOF()*1e9;
1505 if(tof < fTimeCutMin || tof > fTimeCutMax) continue;
1507 //Check only certain regions
1509 if(IsFiducialCutOn()) in = GetFiducialCut()->IsInFiducialCut(mom,fCalorimeter) ;
1513 nLabel = clus->GetNLabels();
1514 labels = clus->GetLabels();
1517 nCaloCellsPerCluster = clus->GetNCells();
1518 //if(mom.E() > 10 && nCaloCellsPerCluster == 1 ) printf("%s:************** E = %f ********** ncells = %d\n",fCalorimeter.Data(), mom.E(),nCaloCellsPerCluster);
1520 //matched cluster with tracks
1521 nTracksMatched = clus->GetNTracksMatched();
1522 if(GetReader()->GetDataType() == AliCaloTrackReader::kESD){
1523 trackIndex = clus->GetTrackMatchedIndex();
1524 if(trackIndex >= 0){
1525 track = (AliVTrack*)GetReader()->GetInputEvent()->GetTrack(trackIndex);
1528 if(nTracksMatched == 1) nTracksMatched = 0;
1533 if(nTracksMatched > 0) track = (AliVTrack*)clus->GetTrackMatched(0);
1536 //======================
1538 //======================
1540 //Get list of contributors
1541 UShort_t * indexList = clus->GetCellsAbsId() ;
1542 // check time of cells respect to max energy cell
1543 //Get maximum energy cell
1545 //printf("nCaloCellsPerCluster %d\n",nCaloCellsPerCluster);
1546 if(fFillAllPosHisto){
1547 //Loop on cluster cells
1548 for (Int_t ipos = 0; ipos < nCaloCellsPerCluster; ipos++) {
1549 // printf("Index %d\n",ipos);
1550 absId = indexList[ipos];
1552 //Get position of cell compare to cluster
1554 if(fCalorimeter=="EMCAL" && GetCaloUtils()->IsEMCALGeoMatrixSet()){
1556 Double_t cellpos[] = {0, 0, 0};
1557 GetEMCALGeometry()->GetGlobal(absId, cellpos);
1559 fhDeltaCellClusterXNCells->Fill(pos[0]-cellpos[0],nCaloCellsPerCluster) ;
1560 fhDeltaCellClusterYNCells->Fill(pos[1]-cellpos[1],nCaloCellsPerCluster) ;
1561 fhDeltaCellClusterZNCells->Fill(pos[2]-cellpos[2],nCaloCellsPerCluster) ;
1563 fhDeltaCellClusterXE->Fill(pos[0]-cellpos[0],mom.E()) ;
1564 fhDeltaCellClusterYE->Fill(pos[1]-cellpos[1],mom.E()) ;
1565 fhDeltaCellClusterZE->Fill(pos[2]-cellpos[2],mom.E()) ;
1567 Float_t r = TMath::Sqrt(pos[0]*pos[0] +pos[1]*pos[1]);// +pos[2]*pos[2]);
1568 Float_t rcell = TMath::Sqrt(cellpos[0]*cellpos[0]+cellpos[1]*cellpos[1]);//+cellpos[2]*cellpos[2]);
1569 fhDeltaCellClusterRNCells->Fill(r-rcell, nCaloCellsPerCluster) ;
1570 fhDeltaCellClusterRE ->Fill(r-rcell, mom.E()) ;
1572 // Float_t celleta = 0, cellphi = 0;
1573 // GetEMCALGeometry()->EtaPhiFromIndex(absId, celleta, cellphi);
1574 // Int_t imod = -1, iTower = -1, iIphi = -1, iIeta = -1, iphi = -1, ieta = -1;
1575 // GetEMCALGeometry()->GetCellIndex(absId,imod,iTower,iIphi,iIeta);
1576 // GetEMCALGeometry()->GetCellPhiEtaIndexInSModule(imod,iTower,
1577 // iIphi, iIeta,iphi,ieta);
1578 // printf("AbsId %d, SM %d, Index eta %d, phi %d\n", absId, imod, ieta, iphi);
1579 // printf("Cluster E %f, eta %f, phi %f; Cell: Amp %f, eta %f, phi%f\n", mom.E(),mom.Eta(), mom.Phi()*TMath::RadToDeg(), cell->GetCellAmplitude(absId),celleta, cellphi*TMath::RadToDeg());
1580 // printf("x cluster %f, x cell %f, cluster-cell %f\n",pos[0], cellpos[0],pos[0]-cellpos[0]);
1581 // printf("y cluster %f, y cell %f, cluster-cell %f\n",pos[1], cellpos[1],pos[1]-cellpos[1]);
1582 // printf("z cluster %f, z cell %f, cluster-cell %f\n",pos[2], cellpos[2],pos[2]-cellpos[2]);
1583 // printf("r cluster %f, r cell %f, cluster-cell %f\n",r, rcell, r-rcell);
1586 }//EMCAL and its matrices are available
1587 else if(fCalorimeter=="PHOS" && GetCaloUtils()->IsPHOSGeoMatrixSet()){
1589 Int_t relId[4], module;
1590 Float_t xCell, zCell;
1592 GetPHOSGeometry()->AbsToRelNumbering(absId,relId);
1594 GetPHOSGeometry()->RelPosInModule(relId,xCell,zCell);
1595 GetPHOSGeometry()->Local2Global(module,xCell,zCell,xyz);
1597 fhDeltaCellClusterXNCells->Fill(pos[0]-xyz.X(),nCaloCellsPerCluster) ;
1598 fhDeltaCellClusterYNCells->Fill(pos[1]-xyz.Y(),nCaloCellsPerCluster) ;
1599 fhDeltaCellClusterZNCells->Fill(pos[2]-xyz.Z(),nCaloCellsPerCluster) ;
1601 fhDeltaCellClusterXE->Fill(pos[0]-xyz.X(),mom.E()) ;
1602 fhDeltaCellClusterYE->Fill(pos[1]-xyz.Y(),mom.E()) ;
1603 fhDeltaCellClusterZE->Fill(pos[2]-xyz.Z(),mom.E()) ;
1605 Float_t r = TMath::Sqrt(pos[0]*pos[0] +pos[1]*pos[1]);// +pos[2]*pos[2]);
1606 Float_t rcell = TMath::Sqrt(xyz.X()*xyz.X()+xyz.Y()*xyz.Y());//+xyz.Z()*xyz.Z());
1607 fhDeltaCellClusterRNCells->Fill(r-rcell, nCaloCellsPerCluster) ;
1608 fhDeltaCellClusterRE ->Fill(r-rcell, mom.E()) ;
1610 // printf("x cluster %f, x cell %f, cluster-cell %f\n",pos[0], cellpos[0],pos[0]-cellpos[0]);
1611 // printf("y cluster %f, y cell %f, cluster-cell %f\n",pos[1], cellpos[1],pos[1]-cellpos[1]);
1612 // printf("z cluster %f, z cell %f, cluster-cell %f\n",pos[2], cellpos[2],pos[2]-cellpos[2]);
1613 // printf("r cluster %f, r cell %f, cluster-cell %f\n",r, rcell, r-rcell);
1614 }//PHOS and its matrices are available
1617 }// cluster cell loop
1618 }//Fill all position histograms
1620 // Get the fraction of the cluster energy that carries the cell with highest energy
1621 Float_t maxCellFraction = 0.;
1622 Int_t absIdMax = GetCaloUtils()->GetMaxEnergyCell(cell, clus,maxCellFraction);
1623 Double_t tmax = cell->GetCellTime(absIdMax)*1e9;
1625 if (clus->E() < 2.){
1626 fhLambda0vsClusterMaxCellDiffE0->Fill(clus->GetM02(), maxCellFraction);
1627 fhNCellsvsClusterMaxCellDiffE0 ->Fill(nCaloCellsPerCluster,maxCellFraction);
1629 else if(clus->E() < 6.){
1630 fhLambda0vsClusterMaxCellDiffE2->Fill(clus->GetM02(), maxCellFraction);
1631 fhNCellsvsClusterMaxCellDiffE2 ->Fill(nCaloCellsPerCluster,maxCellFraction);
1634 fhLambda0vsClusterMaxCellDiffE6->Fill(clus->GetM02(), maxCellFraction);
1635 fhNCellsvsClusterMaxCellDiffE6 ->Fill(nCaloCellsPerCluster,maxCellFraction);
1638 fhNCellsPerClusterNoCut ->Fill(clus->E(), nCaloCellsPerCluster);
1639 nModule = GetModuleNumber(clus);
1640 if(nModule >=0 && nModule < fNModules) fhNCellsPerClusterModNoCut[nModule]->Fill(clus->E(), nCaloCellsPerCluster);
1642 fhClusterMaxCellDiffNoCut->Fill(clus->E(),maxCellFraction);
1643 //fhClusterMaxCellDiffDivLambda0->Fill(clus->E(),maxCellFraction / clus->GetNCells());
1645 //Check bad clusters if rejection was not on
1646 Bool_t badCluster = kFALSE;
1647 if(fCalorimeter=="EMCAL" && !GetCaloUtils()->GetEMCALRecoUtils()->IsRejectExoticCluster()){
1648 //Bad clusters histograms
1649 //Float_t minNCells = TMath::Max(1,TMath::Nint(1 + TMath::Log(clus->E() - 5 )*1.5 ));
1650 //if(nCaloCellsPerCluster <= minNCells) {
1651 if(clus->GetM02() < 0.05) {
1653 //if(clus->GetM02() > 0 || TMath::Abs(clus->GetM20()) > 0 || clus->GetDispersion() > 0)
1655 Int_t sm =0; Int_t ietaa=-1; Int_t iphii = 0; Int_t rcu = 0;
1656 sm = GetModuleNumberCellIndexes(absId,fCalorimeter, ietaa, iphii, rcu);
1657 // if(clus->GetNCells() > 3){
1658 // printf("Bad : E %f, ncells %d, nclusters %d, dist to bad %f, l0 %f, l1 %f, d %f, cell max t %f, cluster TOF %f, sm %d, icol %d, irow %d, rcu %d\n",
1659 // clus->E(), clus->GetNCells(),nCaloClusters, clus->GetDistanceToBadChannel(),
1660 // clus->GetM02(), clus->GetM20(), clus->GetDispersion(),tmax, tof,sm,ietaa,iphii,rcu);
1665 fhBadClusterEnergy ->Fill(clus->E());
1666 fhBadClusterMaxCellDiff->Fill(clus->E(),maxCellFraction);
1667 fhBadClusterTimeEnergy ->Fill(clus->E(),tof);
1668 //printf("bad tof : %2.3f\n",tof);
1669 //if(clus->E() - emax < 0)printf("What?\n");
1671 //Clusters in event time difference
1673 for(Int_t iclus2 = 0; iclus2 < nCaloClusters; iclus2++ ){
1675 AliVCluster* clus2 = (AliVCluster*)caloClusters->At(iclus2);
1677 if(clus->GetID()==clus2->GetID()) continue;
1679 if(clus->GetM02() > 0.01) {
1680 fhBadClusterPairDiffTimeE ->Fill(clus->E(), tof-clus2->GetTOF()*1.e9);
1681 // if(clus->GetNCells()>3) printf("\t i %d, E %f, nCells %d, dist to bad %f, good tof %f, bad tof %f, diff %f \n",
1682 // iclus2, clus2->E(), clus2->GetNCells(), clus->GetDistanceToBadChannel(), tof,clus2->GetTOF()*1.e9,tof-clus2->GetTOF()*1.e9);
1685 // Int_t absId2 = clus2->GetCellsAbsId()[0];
1686 // Int_t sm2 =0; Int_t ietaa2=-1; Int_t iphii2 = 0; Int_t rcu2 = 0;
1687 // sm2 = GetModuleNumberCellIndexes(absId2,fCalorimeter, ietaa2, iphii2, rcu2);
1688 // if(clus->GetNCells()>3) printf("i %d, E %f, nCells %d, bad tof %f, bad tof %f, diff %f, sm %d, icol %d, irow %d, rcu %d\n",
1689 // iclus2, clus2->E(), clus2->GetNCells(), tof,clus2->GetTOF()*1.e9,tof-clus2->GetTOF()*1.e9,sm2,ietaa2,iphii2,rcu2);
1696 for (Int_t ipos = 0; ipos < nCaloCellsPerCluster; ipos++) {
1697 // printf("Index %d\n",ipos);
1698 absId = indexList[ipos];
1699 if(absId!=absIdMax){
1700 Float_t frac = cell->GetCellAmplitude(absId)/cell->GetCellAmplitude(absIdMax);
1701 //printf("bad frac : %2.3f, e %2.2f, ncells %d, min %2.1f\n",frac,mom.E(),nCaloCellsPerCluster,minNCells);
1702 fhBadClusterMaxCellCloseCellRatio->Fill(mom.E(),frac);
1710 // if(TMath::Abs(clus->GetM20()) < 0.0001 && clus->GetNCells() > 3){
1711 // Int_t sm =0; Int_t ietaa=-1; Int_t iphii = 0; Int_t rcu = 0;
1712 // sm = GetModuleNumberCellIndexes(absId,fCalorimeter, ietaa, iphii, rcu);
1713 // printf("Good : E %f, mcells %d, l0 %f, l1 %f, d %f, cell max t %f, cluster TOF %f, sm %d, icol %d, irow %d \n",
1714 // clus->E(), clus->GetNCells(),clus->GetM02(), clus->GetM20(), clus->GetDispersion(),tmax, tof,sm,ietaa,iphii);
1718 fhClusterMaxCellDiff->Fill(clus->E(),maxCellFraction);
1719 fhClusterTimeEnergy ->Fill(mom.E(),tof);
1721 //Clusters in event time difference
1722 for(Int_t iclus2 = 0; iclus2 < nCaloClusters; iclus2++ ){
1724 AliVCluster* clus2 = (AliVCluster*) caloClusters->At(iclus2);
1726 if(clus->GetID()==clus2->GetID()) continue;
1728 if(clus->GetM02() > 0.01) {
1729 fhClusterPairDiffTimeE ->Fill(clus->E(), tof-clus2->GetTOF()*1.e9);
1733 if(nCaloCellsPerCluster > 1){
1735 // check time of cells respect to max energy cell
1737 for (Int_t ipos = 0; ipos < nCaloCellsPerCluster; ipos++) {
1739 absId = indexList[ipos];
1740 if(absId == absIdMax) continue;
1742 Float_t frac = cell->GetCellAmplitude(absId)/cell->GetCellAmplitude(absIdMax);
1743 fhClusterMaxCellCloseCellRatio->Fill(mom.E(),frac);
1744 fhClusterMaxCellCloseCellDiff ->Fill(mom.E(),cell->GetCellAmplitude(absIdMax)-cell->GetCellAmplitude(absId));
1746 if(GetReader()->GetDataType()==AliCaloTrackReader::kESD) {
1747 Float_t diff = (tmax-cell->GetCellTime(absId)*1e9);
1748 fhCellTimeSpreadRespectToCellMax->Fill(mom.E(), diff);
1749 if(TMath::Abs(TMath::Abs(diff) > 100) && mom.E() > 1 ) fhCellIdCellLargeTimeSpread->Fill(absId);
1752 }// fill cell-cluster histogram loop
1753 }//check time and energy of cells respect to max energy cell if cluster of more than 1 cell
1756 //Get module of cluster
1757 nCaloClustersAccepted++;
1758 if(nModule >=0 && nModule < fNModules) nClustersInModule[nModule]++;
1760 //-----------------------------------------------------------
1761 //Fill histograms related to single cluster or track matching
1762 //-----------------------------------------------------------
1763 ClusterHistograms(mom, pos, nCaloCellsPerCluster, nModule, nTracksMatched, track, labels, nLabel);
1766 //-----------------------------------------------------------
1768 //-----------------------------------------------------------
1769 if(fFillAllPi0Histo){
1770 if(GetDebug()>1) printf("Invariant mass \n");
1772 //do not do for bad vertex
1773 // Float_t fZvtxCut = 40. ;
1774 if(v[2]<-GetZvertexCut() || v[2]> GetZvertexCut()) continue ; //Event can not be used (vertex, centrality,... cuts not fulfilled)
1776 Int_t nModule2 = -1;
1777 Int_t nCaloCellsPerCluster2=0;
1778 if (nCaloClusters > 1 ) {
1779 for(Int_t jclus = iclus + 1 ; jclus < nCaloClusters ; jclus++) {
1780 AliVCluster* clus2 = (AliVCluster*)caloClusters->At(jclus);
1782 //Get cluster kinematics
1783 clus2->GetMomentum(mom2,v);
1784 //Check only certain regions
1786 if(IsFiducialCutOn()) in2 = GetFiducialCut()->IsInFiducialCut(mom2,fCalorimeter) ;
1788 //Get module of cluster
1789 nModule2 = GetModuleNumber(clus2);
1791 nCaloCellsPerCluster2 = clus2->GetNCells();
1793 //Fill invariant mass histograms
1796 //printf("QA : Fill inv mass histo: pt1 %f, pt2 %f, pt12 %f, mass %f, calo %s \n",mom.Pt(),mom2.Pt(),(mom+mom2).Pt(),(mom+mom2).M(), fCalorimeter.Data());
1797 fhIM ->Fill((mom+mom2).Pt(),(mom+mom2).M());
1799 if(nModule == nModule2 && nModule >=0 && nModule < fNModules)
1800 fhIMMod[nModule]->Fill((mom+mom2).Pt(),(mom+mom2).M());
1802 //Select only clusters with at least 2 cells
1803 if(nCaloCellsPerCluster > 1 && nCaloCellsPerCluster2 > 1) {
1805 fhIMCellCut ->Fill((mom+mom2).Pt(),(mom+mom2).M());
1807 if(nModule == nModule2 && nModule >=0 && nModule < fNModules)
1808 fhIMCellCutMod[nModule]->Fill((mom+mom2).Pt(),(mom+mom2).M());
1811 //Asymetry histograms
1812 fhAsym->Fill((mom+mom2).Pt(),TMath::Abs((mom.E()-mom2.E())/(mom.E()+mom2.E())));
1814 }// 2nd cluster loop
1821 //Number of clusters histograms
1822 if(nCaloClustersAccepted > 0) fhNClusters->Fill(nCaloClustersAccepted);
1823 // Number of clusters per module
1824 for(Int_t imod = 0; imod < fNModules; imod++ ){
1826 printf("AliAnaCalorimeterQA::MakeAnalysisFillHistograms() - module %d calo %s clusters %d\n", imod, fCalorimeter.Data(), nClustersInModule[imod]);
1827 fhNClustersMod[imod]->Fill(nClustersInModule[imod]);
1829 delete [] nClustersInModule;
1830 //delete caloClusters;
1831 }// calo clusters array exists
1833 //----------------------------------------------------------
1835 //----------------------------------------------------------
1837 AliVCaloCells * cell = 0x0;
1839 if(fCalorimeter == "PHOS")
1840 cell = GetPHOSCells();
1842 cell = GetEMCALCells();
1845 AliFatal(Form("No %s CELLS available for analysis",fCalorimeter.Data()));
1846 return; // just to trick coverity
1850 printf("AliAnaCalorimeterQA::MakeAnalysisFillHistograms() - %s cell entries %d\n", fCalorimeter.Data(), cell->GetNumberOfCells());
1852 //Init arrays and used variables
1853 Int_t *nCellsInModule = new Int_t[fNModules];
1854 for(Int_t imod = 0; imod < fNModules; imod++ ) nCellsInModule[imod] = 0;
1861 Float_t recalF = 1.;
1863 for (Int_t iCell = 0; iCell < cell->GetNumberOfCells(); iCell++) {
1865 printf("AliAnaCalorimeterQA::MakeAnalysisFillHistograms() - Cell : amp %f, absId %d \n", cell->GetAmplitude(iCell), cell->GetCellNumber(iCell));
1866 nModule = GetModuleNumberCellIndexes(cell->GetCellNumber(iCell),fCalorimeter, icol, irow, iRCU);
1868 printf("\t module %d, column %d, row %d \n", nModule,icol,irow);
1870 if(nModule < fNModules) {
1872 //Check if the cell is a bad channel
1873 if(GetCaloUtils()->IsBadChannelsRemovalSwitchedOn()){
1874 if(fCalorimeter=="EMCAL"){
1875 if(GetCaloUtils()->GetEMCALChannelStatus(nModule,icol,irow)) continue;
1878 if(GetCaloUtils()->GetPHOSChannelStatus(nModule,icol,irow)) {
1879 printf("PHOS bad channel\n");
1883 } // use bad channel map
1885 //Get Recalibration factor if set
1886 if (GetCaloUtils()->IsRecalibrationOn()) {
1887 if(fCalorimeter == "PHOS") recalF = GetCaloUtils()->GetPHOSChannelRecalibrationFactor(nModule,icol,irow);
1888 else recalF = GetCaloUtils()->GetEMCALChannelRecalibrationFactor(nModule,icol,irow);
1889 //if(fCalorimeter == "PHOS")printf("Recalibration factor (sm,row,col)=(%d,%d,%d) - %f\n",nModule,icol,irow,recalF);
1892 amp = cell->GetAmplitude(iCell)*recalF;
1893 time = cell->GetTime(iCell)*1e9;//transform time to ns
1895 //Remove noisy channels, only possible in ESDs
1896 if(GetReader()->GetDataType() == AliCaloTrackReader::kESD){
1897 if(time < fTimeCutMin || time > fTimeCutMax) continue;
1899 //if(amp > 3 && fCalorimeter=="EMCAL") printf("Amp = %f, time = %f, (mod, col, row)= (%d,%d,%d)\n",
1900 // amp,time,nModule,icol,irow);
1902 id = cell->GetCellNumber(iCell);
1903 fhAmplitude->Fill(amp);
1904 fhAmpId ->Fill(amp,id);
1906 fhAmplitudeMod[nModule]->Fill(amp);
1907 if(fCalorimeter=="EMCAL"){
1909 if(icol > 15 && icol < 32) ifrac = 1;
1910 else if(icol > 31) ifrac = 2;
1911 fhAmplitudeModFraction[nModule*3+ifrac]->Fill(amp);
1914 nCellsInModule[nModule]++;
1915 fhGridCellsMod[nModule] ->Fill(icol,irow);
1916 fhGridCellsEMod[nModule] ->Fill(icol,irow,amp);
1918 if(GetReader()->GetDataType() == AliCaloTrackReader::kESD){
1919 //printf("%s: time %g\n",fCalorimeter.Data(), time);
1920 fhTime ->Fill(time);
1921 fhTimeId ->Fill(time,id);
1922 fhTimeAmp ->Fill(amp,time);
1924 //Double_t t0 = GetReader()->GetInputEvent()->GetT0();
1925 //printf("---->>> Time EMCal %e, T0 %e, T0 vertex %e, T0 clock %e, T0 trig %d \n",time,t0,
1926 // GetReader()->GetInputEvent()->GetT0zVertex(),
1927 // GetReader()->GetInputEvent()->GetT0clock(),
1928 // GetReader()->GetInputEvent()->GetT0Trig());
1929 //fhT0Time ->Fill(time-t0);
1930 //fhT0TimeId ->Fill(time-t0,id);
1931 //fhT0TimeAmp ->Fill(amp,time-t0);
1933 //printf("id %d, nModule %d, iRCU %d: Histo Name %s\n",id, nModule,iRCU, fhTimeAmpPerRCU[nModule*fNRCU+iRCU]->GetName());
1934 //fhT0TimeAmpPerRCU[nModule*fNRCU+iRCU]->Fill(amp, time-t0);
1936 fhTimeAmpPerRCU [nModule*fNRCU+iRCU]->Fill(amp, time);
1939 fhGridCellsTimeMod[nModule]->Fill(icol,irow,time);
1941 // AliESDCaloCells * cell2 = 0x0;
1942 // if(fCalorimeter == "PHOS") cell2 = GetReader()->GetInputEvent()->GetPHOSCells();
1943 // else cell2 = GetReader()->GetInputEvent()->GetEMCALCells();
1944 // Int_t icol2 = -1;
1945 // Int_t irow2 = -1;
1946 // Int_t iRCU2 = -1;
1947 // Float_t amp2 = 0.;
1948 // Float_t time2 = 0.;
1950 // Int_t nModule2 = -1;
1951 // for (Int_t iCell2 = 0; iCell2 < ncells; iCell2++) {
1952 // amp2 = cell2->GetAmplitude(iCell2);
1953 // if(amp2 < 0.3) continue;
1954 // if(iCell2 == iCell) continue;
1955 // time2 = cell2->GetTime(iCell2)*1e9;//transform time to ns
1956 // //printf("%s: time %g\n",fCalorimeter.Data(), time);
1957 // id2 = cell2->GetCellNumber(iCell2);
1958 // nModule2 = GetModuleNumberCellIndexes(cell2->GetCellNumber(iCell2), fCalorimeter, icol2, irow2, iRCU2);
1959 // Int_t index = (nModule2*fNRCU+iRCU2)+(fNModules*fNRCU)*(iRCU+fNRCU*nModule);
1960 // //printf("id %d, nModule %d, iRCU %d, id2 %d, nModule2 %d, iRCU2 %d, index %d: Histo Name %s\n",id, nModule,iRCU,cell2->GetCellNumber(iCell2),nModule2,iRCU2,index, fhTimeCorrRCU[index]->GetName());
1961 // fhTimeCorrRCU[index]->Fill(time,time2);
1963 // }// second cell loop
1969 //Get Eta-Phi position of Cell
1970 if(fFillAllPosHisto)
1972 if(fCalorimeter=="EMCAL" && GetCaloUtils()->IsEMCALGeoMatrixSet()){
1973 Float_t celleta = 0.;
1974 Float_t cellphi = 0.;
1975 GetEMCALGeometry()->EtaPhiFromIndex(id, celleta, cellphi);
1977 fhEtaPhiAmp->Fill(celleta,cellphi,amp);
1978 Double_t cellpos[] = {0, 0, 0};
1979 GetEMCALGeometry()->GetGlobal(id, cellpos);
1980 fhXCellE->Fill(cellpos[0],amp) ;
1981 fhYCellE->Fill(cellpos[1],amp) ;
1982 fhZCellE->Fill(cellpos[2],amp) ;
1983 Float_t rcell = TMath::Sqrt(cellpos[0]*cellpos[0]+cellpos[1]*cellpos[1]);//+cellpos[2]*cellpos[2]);
1984 fhRCellE->Fill(rcell,amp) ;
1985 fhXYZCell->Fill(cellpos[0],cellpos[1],cellpos[2]) ;
1987 else if(fCalorimeter=="PHOS" && GetCaloUtils()->IsPHOSGeoMatrixSet()){
1989 Int_t relId[4], module;
1990 Float_t xCell, zCell;
1992 GetPHOSGeometry()->AbsToRelNumbering(id,relId);
1994 GetPHOSGeometry()->RelPosInModule(relId,xCell,zCell);
1995 GetPHOSGeometry()->Local2Global(module,xCell,zCell,xyz);
1996 Float_t rcell = TMath::Sqrt(xyz.X()*xyz.X()+xyz.Y()*xyz.Y());
1997 fhXCellE ->Fill(xyz.X(),amp) ;
1998 fhYCellE ->Fill(xyz.Y(),amp) ;
1999 fhZCellE ->Fill(xyz.Z(),amp) ;
2000 fhRCellE ->Fill(rcell ,amp) ;
2001 fhXYZCell->Fill(xyz.X(),xyz.Y(),xyz.Z()) ;
2003 }//fill cell position histograms
2005 if (fCalorimeter=="EMCAL" && amp > fEMCALCellAmpMin) ncells ++ ;
2006 else if(fCalorimeter=="PHOS" && amp > fPHOSCellAmpMin) ncells ++ ;
2008 // printf("AliAnaCalorimeterQA::MakeAnalysisFillHistograms() - no %s CELLS passed the analysis cut\n",fCalorimeter.Data());
2011 if(ncells > 0 )fhNCells->Fill(ncells) ; //fill the cells after the cut
2013 //Number of cells per module
2014 for(Int_t imod = 0; imod < fNModules; imod++ ) {
2016 printf("AliAnaCalorimeterQA::MakeAnalysisFillHistograms() - module %d calo %s cells %d\n", imod, fCalorimeter.Data(), nCellsInModule[imod]);
2017 fhNCellsMod[imod]->Fill(nCellsInModule[imod]) ;
2019 delete [] nCellsInModule;
2022 printf("AliAnaCalorimeterQA::MakeAnalysisFillHistograms() - End \n");
2026 //_____________________________________________________________________________________________
2027 void AliAnaCalorimeterQA::ClusterHistograms(const TLorentzVector mom,
2028 Float_t *pos, const Int_t nCaloCellsPerCluster,const Int_t nModule,
2029 const Int_t nTracksMatched, const AliVTrack * track,
2030 const Int_t * labels, const Int_t nLabels){
2031 //Fill CaloCluster related histograms
2033 AliAODMCParticle * aodprimary = 0x0;
2034 TParticle * primary = 0x0;
2037 Float_t e = mom.E();
2038 Float_t pt = mom.Pt();
2039 Float_t eta = mom.Eta();
2040 Float_t phi = mom.Phi();
2041 if(phi < 0) phi +=TMath::TwoPi();
2042 if(GetDebug() > 0) {
2043 printf("AliAnaCalorimeterQA::ClusterHistograms() - cluster: E %2.3f, pT %2.3f, eta %2.3f, phi %2.3f \n",e,pt,eta,phi*TMath::RadToDeg());
2045 //printf("\t Primaries: nlabels %d, labels pointer %p\n",nLabels,labels);
2046 printf("\t Primaries: nlabels %d\n",nLabels);
2047 if(!nLabels || !labels) printf("\t Strange, no labels!!!\n");
2052 if(nModule >=0 && nModule < fNModules) fhEMod[nModule]->Fill(e);
2059 fhEtaPhiE->Fill(eta,phi,e);
2062 fhNCellsPerCluster ->Fill(e, nCaloCellsPerCluster);
2063 if((fCalorimeter=="EMCAL" && GetReader()->GetEMCALPtMin() < 0.3) ||
2064 (fCalorimeter=="PHOS" && GetReader()->GetPHOSPtMin() < 0.3)) fhNCellsPerClusterMIP->Fill(e, nCaloCellsPerCluster);
2067 if(fFillAllPosHisto2){
2068 fhXE ->Fill(pos[0],e);
2069 fhYE ->Fill(pos[1],e);
2070 fhZE ->Fill(pos[2],e);
2072 fhXYZ ->Fill(pos[0], pos[1],pos[2]);
2074 fhXNCells->Fill(pos[0],nCaloCellsPerCluster);
2075 fhYNCells->Fill(pos[1],nCaloCellsPerCluster);
2076 fhZNCells->Fill(pos[2],nCaloCellsPerCluster);
2077 Float_t rxyz = TMath::Sqrt(pos[0]*pos[0]+pos[1]*pos[1]);//+pos[2]*pos[2]);
2078 fhRE ->Fill(rxyz,e);
2079 fhRNCells->Fill(rxyz ,nCaloCellsPerCluster);
2082 if(nModule >=0 && nModule < fNModules) fhNCellsPerClusterMod[nModule]->Fill(e, nCaloCellsPerCluster);
2084 //Fill histograms only possible when simulation
2085 if(IsDataMC() && nLabels > 0 && labels){
2087 //Play with the MC stack if available
2088 Int_t label = labels[0];
2091 if(GetDebug() >= 0) printf("AliAnaCalorimeterQA::ClusterHistograms() *** bad label ***: label %d \n", label);
2095 Int_t pdg =-1; Int_t pdg0 =-1;Int_t status = -1; Int_t iMother = -1; Int_t iParent = -1;
2096 Float_t vxMC= 0; Float_t vyMC = 0;
2097 Float_t eMC = 0; Float_t ptMC= 0; Float_t phiMC =0; Float_t etaMC = 0;
2101 tag = GetMCAnalysisUtils()->CheckOrigin(labels,nLabels, GetReader(),0);
2103 if(GetReader()->ReadStack() && !GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCUnknown)){ //it MC stack and known tag
2105 if( label >= GetMCStack()->GetNtrack()) {
2106 if(GetDebug() >= 0) printf("AliAnaCalorimeterQA::ClusterHistograms() *** large label ***: label %d, n tracks %d \n", label, GetMCStack()->GetNtrack());
2110 primary = GetMCStack()->Particle(label);
2112 pdg0 = TMath::Abs(primary->GetPdgCode());
2114 status = primary->GetStatusCode();
2115 vxMC = primary->Vx();
2116 vyMC = primary->Vy();
2117 iParent = primary->GetFirstMother();
2119 if(GetDebug() > 1 ) {
2120 printf("AliAnaCalorimeterQA::ClusterHistograms() - Cluster most contributing mother: \n");
2121 printf("\t Mother label %d, pdg %d, %s, status %d, parent %d \n",iMother, pdg0, primary->GetName(),status, iParent);
2124 //Get final particle, no conversion products
2125 if(GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCConversion)){
2127 primary = GetMCStack()->Particle(iParent);
2128 pdg = TMath::Abs(primary->GetPdgCode());
2129 if(GetDebug() > 1 ) printf("AliAnaCalorimeterQA::ClusterHistograms() - Converted cluster!. Find before conversion: \n");
2130 while((pdg == 22 || pdg == 11) && status != 1){
2132 primary = GetMCStack()->Particle(iMother);
2133 status = primary->GetStatusCode();
2134 iParent = primary->GetFirstMother();
2135 pdg = TMath::Abs(primary->GetPdgCode());
2136 if(GetDebug() > 1 )printf("\t pdg %d, index %d, %s, status %d \n",pdg, iMother, primary->GetName(),status);
2139 if(GetDebug() > 1 ) {
2140 printf("AliAnaCalorimeterQA::ClusterHistograms() - Converted Cluster mother before conversion: \n");
2141 printf("\t Mother label %d, pdg %d, %s, status %d, parent %d \n",iMother, pdg, primary->GetName(), status, iParent);
2146 //Overlapped pi0 (or eta, there will be very few), get the meson
2147 if(GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCPi0) ||
2148 GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCEta)){
2149 if(GetDebug() > 1 ) printf("AliAnaCalorimeterQA::ClusterHistograms() - Overlapped Meson decay!, Find it: \n");
2150 while(pdg != 111 && pdg != 221){
2152 primary = GetMCStack()->Particle(iMother);
2153 status = primary->GetStatusCode();
2154 iParent = primary->GetFirstMother();
2155 pdg = TMath::Abs(primary->GetPdgCode());
2156 if(GetDebug() > 1 ) printf("\t pdg %d, %s, index %d\n",pdg, primary->GetName(),iMother);
2158 printf("AliAnaCalorimeterQA::ClusterHistograms() - Tagged as Overlapped photon but meson not found, why?\n");
2163 if(GetDebug() > 2 ) printf("AliAnaCalorimeterQA::ClusterHistograms() - Overlapped %s decay, label %d \n",
2164 primary->GetName(),iMother);
2167 eMC = primary->Energy();
2168 ptMC = primary->Pt();
2169 phiMC = primary->Phi();
2170 etaMC = primary->Eta();
2171 pdg = TMath::Abs(primary->GetPdgCode());
2172 charge = (Int_t) TDatabasePDG::Instance()->GetParticle(pdg)->Charge();
2175 else if(GetReader()->ReadAODMCParticles() && !GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCUnknown)){//it MC AOD and known tag
2176 //Get the list of MC particles
2177 if(!GetReader()->GetAODMCParticles(0))
2178 AliFatal("MCParticles not available!");
2180 aodprimary = (AliAODMCParticle*) (GetReader()->GetAODMCParticles(0))->At(label);
2182 pdg0 = TMath::Abs(aodprimary->GetPdgCode());
2184 status = aodprimary->IsPrimary();
2185 vxMC = aodprimary->Xv();
2186 vyMC = aodprimary->Yv();
2187 iParent = aodprimary->GetMother();
2189 if(GetDebug() > 1 ) {
2190 printf("AliAnaCalorimeterQA::ClusterHistograms() - Cluster most contributing mother: \n");
2191 printf("\t Mother label %d, pdg %d, Primary? %d, Physical Primary? %d, parent %d \n",
2192 iMother, pdg0, aodprimary->IsPrimary(), aodprimary->IsPhysicalPrimary(), iParent);
2195 //Get final particle, no conversion products
2196 if(GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCConversion)){
2198 printf("AliAnaCalorimeterQA::ClusterHistograms() - Converted cluster!. Find before conversion: \n");
2200 aodprimary = (AliAODMCParticle*)(GetReader()->GetAODMCParticles(0))->At(iParent);
2201 pdg = TMath::Abs(aodprimary->GetPdgCode());
2202 while ((pdg == 22 || pdg == 11) && !aodprimary->IsPhysicalPrimary()) {
2204 aodprimary = (AliAODMCParticle*)(GetReader()->GetAODMCParticles(0))->At(iMother);
2205 status = aodprimary->IsPrimary();
2206 iParent = aodprimary->GetMother();
2207 pdg = TMath::Abs(aodprimary->GetPdgCode());
2209 printf("\t pdg %d, index %d, Primary? %d, Physical Primary? %d \n",
2210 pdg, iMother, aodprimary->IsPrimary(), aodprimary->IsPhysicalPrimary());
2213 if(GetDebug() > 1 ) {
2214 printf("AliAnaCalorimeterQA::ClusterHistograms() - Converted Cluster mother before conversion: \n");
2215 printf("\t Mother label %d, pdg %d, parent %d, Primary? %d, Physical Primary? %d \n",
2216 iMother, pdg, iParent, aodprimary->IsPrimary(), aodprimary->IsPhysicalPrimary());
2221 //Overlapped pi0 (or eta, there will be very few), get the meson
2222 if(GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCPi0) ||
2223 GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCEta)){
2224 if(GetDebug() > 1 ) printf("AliAnaCalorimeterQA::ClusterHistograms() - Overlapped Meson decay!, Find it: PDG %d, mom %d \n",pdg, iMother);
2225 while(pdg != 111 && pdg != 221){
2228 aodprimary = (AliAODMCParticle*)(GetReader()->GetAODMCParticles(0))->At(iMother);
2229 status = aodprimary->IsPrimary();
2230 iParent = aodprimary->GetMother();
2231 pdg = TMath::Abs(aodprimary->GetPdgCode());
2233 if(GetDebug() > 1 ) printf("\t pdg %d, index %d\n",pdg, iMother);
2236 printf("AliAnaCalorimeterQA::ClusterHistograms() - Tagged as Overlapped photon but meson not found, why?\n");
2241 if(GetDebug() > 2 ) printf("AliAnaCalorimeterQA::ClusterHistograms() - Overlapped %s decay, label %d \n",
2242 aodprimary->GetName(),iMother);
2245 status = aodprimary->IsPrimary();
2246 eMC = aodprimary->E();
2247 ptMC = aodprimary->Pt();
2248 phiMC = aodprimary->Phi();
2249 etaMC = aodprimary->Eta();
2250 pdg = TMath::Abs(aodprimary->GetPdgCode());
2251 charge = aodprimary->Charge();
2255 //Float_t vz = primary->Vz();
2256 Float_t rVMC = TMath::Sqrt(vxMC*vxMC + vyMC*vyMC);
2257 if((pdg == 22 || TMath::Abs(pdg)==11) && status!=1) {
2258 fhEMVxyz ->Fill(vxMC,vyMC);//,vz);
2259 fhEMR ->Fill(e,rVMC);
2262 //printf("reco e %f, pt %f, phi %f, eta %f \n", e, pt, phi, eta);
2263 //printf("prim e %f, pt %f, phi %f, eta %f \n", eMC,ptMC,phiMC ,etaMC );
2264 //printf("vertex: vx %f, vy %f, vz %f, r %f \n", vxMC, vyMC, vz, r);
2267 fh2E ->Fill(e, eMC);
2268 fh2Pt ->Fill(pt, ptMC);
2269 fh2Phi ->Fill(phi, phiMC);
2270 fh2Eta ->Fill(eta, etaMC);
2271 fhDeltaE ->Fill(eMC-e);
2272 fhDeltaPt ->Fill(ptMC-pt);
2273 fhDeltaPhi->Fill(phiMC-phi);
2274 fhDeltaEta->Fill(etaMC-eta);
2275 if(eMC > 0) fhRatioE ->Fill(e/eMC);
2276 if(ptMC > 0) fhRatioPt ->Fill(pt/ptMC);
2277 if(phiMC > 0) fhRatioPhi->Fill(phi/phiMC);
2278 if(etaMC > 0) fhRatioEta->Fill(eta/etaMC);
2281 //Overlapped pi0 (or eta, there will be very few)
2282 if(GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCPi0) ||
2283 GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCEta)){
2284 fhPi0E ->Fill(e,eMC);
2285 fhPi0Pt ->Fill(pt,ptMC);
2286 fhPi0Eta ->Fill(eta,etaMC);
2287 fhPi0Phi ->Fill(phi,phiMC);
2288 if( nTracksMatched > 0){
2289 fhPi0ECharged ->Fill(e,eMC);
2290 fhPi0PtCharged ->Fill(pt,ptMC);
2291 fhPi0PhiCharged ->Fill(phi,phiMC);
2292 fhPi0EtaCharged ->Fill(eta,etaMC);
2294 }//Overlapped pizero decay
2295 else if(GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCPhoton)){
2296 fhGamE ->Fill(e,eMC);
2297 fhGamPt ->Fill(pt,ptMC);
2298 fhGamEta ->Fill(eta,etaMC);
2299 fhGamPhi ->Fill(phi,phiMC);
2300 fhGamDeltaE ->Fill(eMC-e);
2301 fhGamDeltaPt ->Fill(ptMC-pt);
2302 fhGamDeltaPhi->Fill(phiMC-phi);
2303 fhGamDeltaEta->Fill(etaMC-eta);
2304 if(eMC > 0) fhGamRatioE ->Fill(e/eMC);
2305 if(ptMC > 0) fhGamRatioPt ->Fill(pt/ptMC);
2306 if(phiMC > 0) fhGamRatioPhi->Fill(phi/phiMC);
2307 if(etaMC > 0) fhGamRatioEta->Fill(eta/etaMC);
2308 if( nTracksMatched > 0){
2309 fhGamECharged ->Fill(e,eMC);
2310 fhGamPtCharged ->Fill(pt,ptMC);
2311 fhGamPhiCharged ->Fill(phi,phiMC);
2312 fhGamEtaCharged ->Fill(eta,etaMC);
2315 else if(GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCElectron)){
2316 fhEleE ->Fill(e,eMC);
2317 fhElePt ->Fill(pt,ptMC);
2318 fhEleEta ->Fill(eta,etaMC);
2319 fhElePhi ->Fill(phi,phiMC);
2320 fhEMVxyz ->Fill(vxMC,vyMC);//,vz);
2321 fhEMR ->Fill(e,rVMC);
2322 if( nTracksMatched > 0){
2323 fhEleECharged ->Fill(e,eMC);
2324 fhElePtCharged ->Fill(pt,ptMC);
2325 fhElePhiCharged ->Fill(phi,phiMC);
2326 fhEleEtaCharged ->Fill(eta,etaMC);
2329 else if(charge == 0){
2330 fhNeHadE ->Fill(e,eMC);
2331 fhNeHadPt ->Fill(pt,ptMC);
2332 fhNeHadEta ->Fill(eta,etaMC);
2333 fhNeHadPhi ->Fill(phi,phiMC);
2334 fhHaVxyz ->Fill(vxMC,vyMC);//,vz);
2335 fhHaR ->Fill(e,rVMC);
2336 if( nTracksMatched > 0){
2337 fhNeHadECharged ->Fill(e,eMC);
2338 fhNeHadPtCharged ->Fill(pt,ptMC);
2339 fhNeHadPhiCharged ->Fill(phi,phiMC);
2340 fhNeHadEtaCharged ->Fill(eta,etaMC);
2344 fhChHadE ->Fill(e,eMC);
2345 fhChHadPt ->Fill(pt,ptMC);
2346 fhChHadEta ->Fill(eta,etaMC);
2347 fhChHadPhi ->Fill(phi,phiMC);
2348 fhHaVxyz ->Fill(vxMC,vyMC);//,vz);
2349 fhHaR ->Fill(e,rVMC);
2350 if( nTracksMatched > 0){
2351 fhChHadECharged ->Fill(e,eMC);
2352 fhChHadPtCharged ->Fill(pt,ptMC);
2353 fhChHadPhiCharged ->Fill(phi,phiMC);
2354 fhChHadEtaCharged ->Fill(eta,etaMC);
2360 //Match tracks and clusters
2361 //To be Modified in case of AODs
2363 if( nTracksMatched > 0 && fFillAllTMHisto){
2364 if(fFillAllTH12 && fFillAllTMHisto){
2365 fhECharged ->Fill(e);
2366 fhPtCharged ->Fill(pt);
2367 fhPhiCharged ->Fill(phi);
2368 fhEtaCharged ->Fill(eta);
2371 if(fFillAllTMHisto){
2372 if(fFillAllTH3)fhEtaPhiECharged->Fill(eta,phi,e);
2373 if((fCalorimeter=="EMCAL" && GetReader()->GetEMCALPtMin() < 0.3) ||
2374 (fCalorimeter=="PHOS" && GetReader()->GetPHOSPtMin() < 0.3)) fhNCellsPerClusterMIPCharged->Fill(e, nCaloCellsPerCluster);
2376 //printf("track index %d ntracks %d\n", esd->GetNumberOfTracks());
2377 //Study the track and matched cluster if track exists.
2379 Double_t emcpos[3] = {0.,0.,0.};
2380 Double_t emcmom[3] = {0.,0.,0.};
2381 Double_t radius = 441.0; //[cm] EMCAL radius +13cm
2382 Double_t bfield = 0.;
2388 Double_t tpcSignal = 0;
2389 Bool_t okpos = kFALSE;
2390 Bool_t okmom = kFALSE;
2391 Bool_t okout = kFALSE;
2395 //In case of ESDs get the parameters in this way
2396 if(GetReader()->GetDataType()==AliCaloTrackReader::kESD) {
2397 if (track->GetOuterParam() ) {
2400 bfield = GetReader()->GetInputEvent()->GetMagneticField();
2401 okpos = track->GetOuterParam()->GetXYZAt(radius,bfield,emcpos);
2402 okmom = track->GetOuterParam()->GetPxPyPzAt(radius,bfield,emcmom);
2403 if(!(okpos && okmom)) return;
2405 TVector3 position(emcpos[0],emcpos[1],emcpos[2]);
2406 TVector3 momentum(emcmom[0],emcmom[1],emcmom[2]);
2407 tphi = position.Phi();
2408 teta = position.Eta();
2409 tmom = momentum.Mag();
2411 //Double_t tphi = track->GetOuterParam()->Phi();
2412 //Double_t teta = track->GetOuterParam()->Eta();
2413 //Double_t tmom = track->GetOuterParam()->P();
2416 tpcSignal = track->GetTPCsignal();
2418 nITS = track->GetNcls(0);
2419 nTPC = track->GetNcls(1);
2420 }//Outer param available
2422 else if(GetReader()->GetDataType()==AliCaloTrackReader::kAOD) {
2423 AliAODPid* pid = (AliAODPid*) ((AliAODTrack *) track)->GetDetPid();
2426 pid->GetEMCALPosition(emcpos);
2427 pid->GetEMCALMomentum(emcmom);
2429 TVector3 position(emcpos[0],emcpos[1],emcpos[2]);
2430 TVector3 momentum(emcmom[0],emcmom[1],emcmom[2]);
2431 tphi = position.Phi();
2432 teta = position.Eta();
2433 tmom = momentum.Mag();
2437 tpcSignal = pid->GetTPCsignal();
2439 //nITS = ((AliAODTrack*)track)->GetNcls(0);
2440 //nTPC = ((AliAODTrack*)track)->GetNcls(1);
2445 Double_t deta = teta - eta;
2446 Double_t dphi = tphi - phi;
2447 if(dphi > TMath::Pi()) dphi -= 2*TMath::Pi();
2448 if(dphi < -TMath::Pi()) dphi += 2*TMath::Pi();
2449 Double_t dR = sqrt(dphi*dphi + deta*deta);
2451 Double_t pOverE = tmom/e;
2453 fh1pOverE->Fill(tpt, pOverE);
2454 if(dR < 0.02) fh1pOverER02->Fill(tpt,pOverE);
2457 fh2MatchdEdx->Fill(tmom2,tpcSignal);
2459 if(IsDataMC() && primary){
2460 Int_t pdg = primary->GetPdgCode();
2461 Double_t charge = TDatabasePDG::Instance()->GetParticle(pdg)->Charge();
2463 if(TMath::Abs(pdg) == 11){
2464 fhMCEle1pOverE->Fill(tpt,pOverE);
2465 fhMCEle1dR->Fill(dR);
2466 fhMCEle2MatchdEdx->Fill(tmom2,tpcSignal);
2467 if(dR < 0.02) fhMCEle1pOverER02->Fill(tpt,pOverE);
2470 fhMCChHad1pOverE->Fill(tpt,pOverE);
2471 fhMCChHad1dR->Fill(dR);
2472 fhMCChHad2MatchdEdx->Fill(tmom2,tpcSignal);
2473 if(dR < 0.02) fhMCChHad1pOverER02->Fill(tpt,pOverE);
2475 else if(charge == 0){
2476 fhMCNeutral1pOverE->Fill(tpt,pOverE);
2477 fhMCNeutral1dR->Fill(dR);
2478 fhMCNeutral2MatchdEdx->Fill(tmom2,tpcSignal);
2479 if(dR < 0.02) fhMCNeutral1pOverER02->Fill(tpt,pOverE);
2483 if(dR < 0.02 && pOverE > 0.5 && pOverE < 1.5
2484 && nCaloCellsPerCluster > 1 && nITS > 3 && nTPC > 20) {
2485 fh2EledEdx->Fill(tmom2,tpcSignal);
2488 else{//no ESD external param or AODPid
2489 // ULong_t status=AliESDtrack::kTPCrefit;
2490 // status|=AliESDtrack::kITSrefit;
2491 //printf("track status %d\n", track->GetStatus() );
2492 // fhEChargedNoOut ->Fill(e);
2493 // fhPtChargedNoOut ->Fill(pt);
2494 // fhPhiChargedNoOut ->Fill(phi);
2495 // fhEtaChargedNoOut ->Fill(eta);
2496 // fhEtaPhiChargedNoOut ->Fill(eta,phi);
2497 // if(GetDebug() >= 0 && ((track->GetStatus() & status) == status)) printf("ITS+TPC\n");
2498 if(GetDebug() >= 0) printf("No ESD external param or AliAODPid \n");
2501 }//matched clusters with tracks
2506 //__________________________________
2507 void AliAnaCalorimeterQA::Correlate(){
2508 // Correlate information from PHOS and EMCAL and with V0 and track multiplicity
2511 TObjArray * caloClustersEMCAL = GetEMCALClusters();
2512 TObjArray * caloClustersPHOS = GetPHOSClusters();
2514 Int_t nclEMCAL = caloClustersEMCAL->GetEntriesFast();
2515 Int_t nclPHOS = caloClustersPHOS ->GetEntriesFast();
2517 Float_t sumClusterEnergyEMCAL = 0;
2518 Float_t sumClusterEnergyPHOS = 0;
2520 for(iclus = 0 ; iclus < caloClustersEMCAL->GetEntriesFast() ; iclus++)
2521 sumClusterEnergyEMCAL += ((AliVCluster*)caloClustersEMCAL->At(iclus))->E();
2522 for(iclus = 0 ; iclus < caloClustersPHOS->GetEntriesFast(); iclus++)
2523 sumClusterEnergyPHOS += ((AliVCluster*)caloClustersPHOS->At(iclus))->E();
2528 AliVCaloCells * cellsEMCAL = GetEMCALCells();
2529 AliVCaloCells * cellsPHOS = GetPHOSCells();
2531 Int_t ncellsEMCAL = cellsEMCAL->GetNumberOfCells();
2532 Int_t ncellsPHOS = cellsPHOS ->GetNumberOfCells();
2534 Float_t sumCellEnergyEMCAL = 0;
2535 Float_t sumCellEnergyPHOS = 0;
2537 for(icell = 0 ; icell < cellsEMCAL->GetNumberOfCells() ; icell++)
2538 sumCellEnergyEMCAL += cellsEMCAL->GetAmplitude(icell);
2539 for(icell = 0 ; icell < cellsPHOS->GetNumberOfCells(); icell++)
2540 sumCellEnergyPHOS += cellsPHOS->GetAmplitude(icell);
2544 fhCaloCorrNClusters->Fill(nclEMCAL,nclPHOS);
2545 fhCaloCorrEClusters->Fill(sumClusterEnergyEMCAL,sumClusterEnergyPHOS);
2546 fhCaloCorrNCells ->Fill(ncellsEMCAL,ncellsPHOS);
2547 fhCaloCorrECells ->Fill(sumCellEnergyEMCAL,sumCellEnergyPHOS);
2549 Int_t v0S = GetV0Signal(0)+GetV0Signal(1);
2550 Int_t v0M = GetV0Multiplicity(0)+GetV0Multiplicity(1);
2551 Int_t trM = GetTrackMultiplicity();
2552 if(fCalorimeter=="PHOS"){
2553 fhCaloV0MCorrNClusters ->Fill(v0M,nclPHOS);
2554 fhCaloV0MCorrEClusters ->Fill(v0M,sumClusterEnergyPHOS);
2555 fhCaloV0MCorrNCells ->Fill(v0M,ncellsPHOS);
2556 fhCaloV0MCorrECells ->Fill(v0M,sumCellEnergyPHOS);
2558 fhCaloV0SCorrNClusters ->Fill(v0S,nclPHOS);
2559 fhCaloV0SCorrEClusters ->Fill(v0S,sumClusterEnergyPHOS);
2560 fhCaloV0SCorrNCells ->Fill(v0S,ncellsPHOS);
2561 fhCaloV0SCorrECells ->Fill(v0S,sumCellEnergyPHOS);
2563 fhCaloTrackMCorrNClusters->Fill(trM,nclPHOS);
2564 fhCaloTrackMCorrEClusters->Fill(trM,sumClusterEnergyPHOS);
2565 fhCaloTrackMCorrNCells ->Fill(trM,ncellsPHOS);
2566 fhCaloTrackMCorrECells ->Fill(trM,sumCellEnergyPHOS);
2569 fhCaloV0MCorrNClusters ->Fill(v0M,nclEMCAL);
2570 fhCaloV0MCorrEClusters ->Fill(v0M,sumClusterEnergyEMCAL);
2571 fhCaloV0MCorrNCells ->Fill(v0M,ncellsEMCAL);
2572 fhCaloV0MCorrECells ->Fill(v0M,sumCellEnergyEMCAL);
2574 fhCaloV0SCorrNClusters ->Fill(v0S,nclEMCAL);
2575 fhCaloV0SCorrEClusters ->Fill(v0S,sumClusterEnergyEMCAL);
2576 fhCaloV0SCorrNCells ->Fill(v0S,ncellsEMCAL);
2577 fhCaloV0SCorrECells ->Fill(v0S,sumCellEnergyEMCAL);
2579 fhCaloTrackMCorrNClusters->Fill(trM,nclEMCAL);
2580 fhCaloTrackMCorrEClusters->Fill(trM,sumClusterEnergyEMCAL);
2581 fhCaloTrackMCorrNCells ->Fill(trM,ncellsEMCAL);
2582 fhCaloTrackMCorrECells ->Fill(trM,sumCellEnergyEMCAL);
2587 printf("AliAnaCalorimeterQA::Correlate(): \n");
2588 printf("\t EMCAL: N cells %d, N clusters %d, summed E cells %f, summed E clusters %f \n",
2589 ncellsEMCAL,nclEMCAL, sumCellEnergyEMCAL,sumClusterEnergyEMCAL);
2590 printf("\t PHOS : N cells %d, N clusters %d, summed E cells %f, summed E clusters %f \n",
2591 ncellsPHOS,nclPHOS,sumCellEnergyPHOS,sumClusterEnergyPHOS);
2592 printf("\t V0 : Signal %d, Multiplicity %d, Track Multiplicity %d \n", v0S,v0M,trM);
2597 //______________________________________________________________________________
2598 void AliAnaCalorimeterQA::MCHistograms(const TLorentzVector mom, const Int_t pdg){
2599 //Fill pure monte carlo related histograms
2601 Float_t eMC = mom.E();
2602 Float_t ptMC = mom.Pt();
2603 Float_t phiMC = mom.Phi();
2605 phiMC += TMath::TwoPi();
2606 Float_t etaMC = mom.Eta();
2608 if (TMath::Abs(etaMC) > 1) return;
2611 if(IsFiducialCutOn()) in = GetFiducialCut()->IsInFiducialCut(mom,fCalorimeter) ;
2614 fhGenGamPt ->Fill(ptMC);
2615 fhGenGamEta->Fill(etaMC);
2616 fhGenGamPhi->Fill(phiMC);
2618 fhGenGamAccE ->Fill(eMC);
2619 fhGenGamAccPt ->Fill(ptMC);
2620 fhGenGamAccEta->Fill(etaMC);
2621 fhGenGamAccPhi->Fill(phiMC);
2624 else if (pdg==111) {
2625 fhGenPi0Pt ->Fill(ptMC);
2626 fhGenPi0Eta->Fill(etaMC);
2627 fhGenPi0Phi->Fill(phiMC);
2629 fhGenPi0AccE ->Fill(eMC);
2630 fhGenPi0AccPt ->Fill(ptMC);
2631 fhGenPi0AccEta->Fill(etaMC);
2632 fhGenPi0AccPhi->Fill(phiMC);
2635 else if (pdg==221) {
2636 fhGenEtaPt ->Fill(ptMC);
2637 fhGenEtaEta->Fill(etaMC);
2638 fhGenEtaPhi->Fill(phiMC);
2640 else if (pdg==223) {
2641 fhGenOmegaPt ->Fill(ptMC);
2642 fhGenOmegaEta->Fill(etaMC);
2643 fhGenOmegaPhi->Fill(phiMC);
2645 else if (TMath::Abs(pdg)==11) {
2646 fhGenElePt ->Fill(ptMC);
2647 fhGenEleEta->Fill(etaMC);
2648 fhGenElePhi->Fill(phiMC);
git://git.uio.no / u / mrichter / AliRoot.git / blob