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"
36 #include <TObjString.h>
38 //---- AliRoot system ----
39 #include "AliAnaCalorimeterQA.h"
40 #include "AliCaloTrackReader.h"
42 #include "AliVCaloCells.h"
43 #include "AliFiducialCut.h"
44 #include "AliAODTrack.h"
45 #include "AliVCluster.h"
46 #include "AliVEvent.h"
47 #include "AliVEventHandler.h"
48 #include "AliAnalysisManager.h"
49 #include "AliAODMCParticle.h"
50 #include "AliMCAnalysisUtils.h"
51 #include "AliAODPid.h"
52 #include "AliExternalTrackParam.h"
54 ClassImp(AliAnaCalorimeterQA)
56 //____________________________________________________________________________
57 AliAnaCalorimeterQA::AliAnaCalorimeterQA() :
58 AliAnaPartCorrBaseClass(),
60 fFillAllPosHisto(kFALSE), fFillAllPosHisto2(kTRUE),
61 fFillAllTH12(kFALSE), fFillAllTH3(kTRUE),
62 fFillAllTMHisto(kTRUE), fFillAllPi0Histo(kTRUE),
64 fNModules(12), fNRCU(2),
65 fNMaxCols(48), fNMaxRows(24),
66 fTimeCutMin(-1), fTimeCutMax(9999999),
67 fEMCALCellAmpMin(0), fPHOSCellAmpMin(0),
69 fhPhi(0), fhEta(0), fhEtaPhiE(0),
70 fhECharged(0), fhPtCharged(0),
71 fhPhiCharged(0), fhEtaCharged(0), fhEtaPhiECharged(0),
75 fhNCellsPerCluster(0), fhNCellsPerClusterNoCut(0),
76 fhNCellsPerClusterMIP(0), fhNCellsPerClusterMIPCharged(0),
77 fhNCellsvsClusterMaxCellDiffE0(0), fhNCellsvsClusterMaxCellDiffE2(0),
78 fhNCellsvsClusterMaxCellDiffE6(0), fhNClusters(0),
81 fhClusterTimeEnergy(0), fhCellTimeSpreadRespectToCellMax(0),
82 fhClusterMaxCellDiffAverageTime(0), fhClusterMaxCellDiffWeightTime(0),
83 fhClusterDiffWeightAverTime(0),
84 fhClusterMaxCellDiffAverageNoMaxTime(0), fhClusterMaxCellDiffWeightNoMaxTime(0),
85 fhClusterNoMaxCellWeight(0),
86 fhCellIdCellLargeTimeSpread(0), fhClusterPairDiffTimeE(0),
88 fhClusterMaxCellCloseCellRatio(0), fhClusterMaxCellCloseCellDiff(0),
89 fhClusterMaxCellDiff(0), fhClusterMaxCellDiffNoCut(0),
90 fhLambda0vsClusterMaxCellDiffE0(0), fhLambda0vsClusterMaxCellDiffE2(0), fhLambda0vsClusterMaxCellDiffE6(0),
93 fhBadClusterEnergy(0), fhBadClusterTimeEnergy(0), fhBadClusterPairDiffTimeE(0),
94 fhBadClusterMaxCellCloseCellRatio(0), fhBadClusterMaxCellCloseCellDiff(0), fhBadClusterMaxCellDiff(0),
95 fhBadClusterMaxCellDiffAverageTime(0), fhBadClusterMaxCellDiffWeightTime(0),
96 fhBadClusterDiffWeightAverTime(0),
97 fhBadClusterMaxCellDiffAverageNoMaxTime(0), fhBadClusterMaxCellDiffWeightNoMaxTime(0),
98 fhBadClusterNoMaxCellWeight(0), fhBadCellTimeSpreadRespectToCellMax(0),
99 fhBadClusterL0(0), fhBadClusterL1(0), fhBadClusterD(0),
102 fhDeltaIEtaDeltaIPhiE0(), fhDeltaIEtaDeltaIPhiE2(), fhDeltaIEtaDeltaIPhiE6(),
103 fhDeltaIA(), fhDeltaIAL0(),
104 fhDeltaIAL1(), fhDeltaIANCells(),
107 fhRNCells(0), fhXNCells(0),
108 fhYNCells(0), fhZNCells(0),
112 fhRCellE(0), fhXCellE(0),
113 fhYCellE(0), fhZCellE(0),
115 fhDeltaCellClusterRNCells(0), fhDeltaCellClusterXNCells(0),
116 fhDeltaCellClusterYNCells(0), fhDeltaCellClusterZNCells(0),
117 fhDeltaCellClusterRE(0), fhDeltaCellClusterXE(0),
118 fhDeltaCellClusterYE(0), fhDeltaCellClusterZE(0),
120 fhNCells(0), fhAmplitude(0),
121 fhAmpId(0), fhEtaPhiAmp(0),
122 fhTime(0), fhTimeId(0), fhTimeAmp(0),
123 fhCaloCorrNClusters(0), fhCaloCorrEClusters(0),
124 fhCaloCorrNCells(0), fhCaloCorrECells(0),
125 fhCaloV0SCorrNClusters(0), fhCaloV0SCorrEClusters(0),
126 fhCaloV0SCorrNCells(0), fhCaloV0SCorrECells(0),
127 fhCaloV0MCorrNClusters(0), fhCaloV0MCorrEClusters(0),
128 fhCaloV0MCorrNCells(0), fhCaloV0MCorrECells(0),
129 fhCaloTrackMCorrNClusters(0), fhCaloTrackMCorrEClusters(0),
130 fhCaloTrackMCorrNCells(0), fhCaloTrackMCorrECells(0),
131 //Super-Module dependent histgrams
132 fhEMod(0), fhNClustersMod(0),
133 fhNCellsPerClusterMod(0), fhNCellsPerClusterModNoCut(0),
135 fhGridCellsMod(0), fhGridCellsEMod(0), fhGridCellsTimeMod(0),
136 fhTimeAmpPerRCU(0), fhIMMod(0),
139 fhRecoMCE(), fhRecoMCPhi(), fhRecoMCEta(),
140 fhRecoMCDeltaE(), fhRecoMCDeltaPhi(), fhRecoMCDeltaEta(),
143 fhGenMCE(), fhGenMCEtaPhi(),
144 fhGenMCAccE(), fhGenMCAccEtaPhi(),
147 fhEMVxyz(0), fhEMR(0), fhHaVxyz(0), fhHaR(0),
148 fh1pOverE(0), fh1dR(0), fh2EledEdx(0), fh2MatchdEdx(0),
149 fhMCEle1pOverE(0), fhMCEle1dR(0), fhMCEle2MatchdEdx(0),
150 fhMCChHad1pOverE(0), fhMCChHad1dR(0), fhMCChHad2MatchdEdx(0),
151 fhMCNeutral1pOverE(0), fhMCNeutral1dR(0), fhMCNeutral2MatchdEdx(0),fh1pOverER02(0),
152 fhMCEle1pOverER02(0), fhMCChHad1pOverER02(0), fhMCNeutral1pOverER02(0)
156 //Initialize parameters
160 //________________________________________________________________________
161 TObjString * AliAnaCalorimeterQA::GetAnalysisCuts()
163 //Save parameters used for analysis
164 TString parList ; //this will be list of parameters used for this analysis.
165 const Int_t buffersize = 255;
166 char onePar[buffersize] ;
168 snprintf(onePar,buffersize,"--- AliAnaCalorimeterQA ---\n") ;
170 snprintf(onePar,buffersize,"Calorimeter: %s\n",fCalorimeter.Data()) ;
172 snprintf(onePar,buffersize,"Time Cut : %2.2f < T < %2.2f ns \n",fTimeCutMin, fTimeCutMax) ;
174 snprintf(onePar,buffersize,"PHOS Cell Amplitude > %2.2f GeV, EMCAL Cell Amplitude > %2.2f GeV \n",fPHOSCellAmpMin, fEMCALCellAmpMin) ;
176 //Get parameters set in base class.
177 //parList += GetBaseParametersList() ;
179 //Get parameters set in FiducialCut class (not available yet)
180 //parlist += GetFidCut()->GetFidCutParametersList()
182 return new TObjString(parList) ;
186 //________________________________________________________________________
187 TList * AliAnaCalorimeterQA::GetCreateOutputObjects()
189 // Create histograms to be saved in output file and
190 // store them in outputContainer
192 TList * outputContainer = new TList() ;
193 outputContainer->SetName("QAHistos") ;
196 Int_t nptbins = GetHistoPtBins(); Float_t ptmax = GetHistoPtMax(); Float_t ptmin = GetHistoPtMin();
197 Int_t nfineptbins = GetHistoFinePtBins(); Float_t ptfinemax = GetHistoFinePtMax(); Float_t ptfinemin = GetHistoFinePtMin();
198 Int_t nphibins = GetHistoPhiBins(); Float_t phimax = GetHistoPhiMax(); Float_t phimin = GetHistoPhiMin();
199 Int_t netabins = GetHistoEtaBins(); Float_t etamax = GetHistoEtaMax(); Float_t etamin = GetHistoEtaMin();
200 Int_t nmassbins = GetHistoMassBins(); Float_t massmax = GetHistoMassMax(); Float_t massmin = GetHistoMassMin();
201 Int_t nasymbins = GetHistoAsymmetryBins(); Float_t asymmax = GetHistoAsymmetryMax(); Float_t asymmin = GetHistoAsymmetryMin();
202 Int_t nPoverEbins = GetHistoPOverEBins(); Float_t pOverEmax = GetHistoPOverEMax(); Float_t pOverEmin = GetHistoPOverEMin();
203 Int_t ndedxbins = GetHistodEdxBins(); Float_t dedxmax = GetHistodEdxMax(); Float_t dedxmin = GetHistodEdxMin();
204 Int_t ndRbins = GetHistodRBins(); Float_t dRmax = GetHistodRMax(); Float_t dRmin = GetHistodRMin();
205 Int_t ntimebins = GetHistoTimeBins(); Float_t timemax = GetHistoTimeMax(); Float_t timemin = GetHistoTimeMin();
206 Int_t nclbins = GetHistoNClustersBins(); Int_t nclmax = GetHistoNClustersMax(); Int_t nclmin = GetHistoNClustersMin();
207 Int_t ncebins = GetHistoNCellsBins(); Int_t ncemax = GetHistoNCellsMax(); Int_t ncemin = GetHistoNCellsMin();
208 Int_t nceclbins = GetHistoNClusterCellBins(); Int_t nceclmax = GetHistoNClusterCellMax(); Int_t nceclmin = GetHistoNClusterCellMin();
209 Int_t nvdistbins = GetHistoVertexDistBins(); Float_t vdistmax = GetHistoVertexDistMax(); Float_t vdistmin = GetHistoVertexDistMin();
210 Int_t rbins = GetHistoRBins(); Float_t rmax = GetHistoRMax(); Float_t rmin = GetHistoRMin();
211 Int_t xbins = GetHistoXBins(); Float_t xmax = GetHistoXMax(); Float_t xmin = GetHistoXMin();
212 Int_t ybins = GetHistoYBins(); Float_t ymax = GetHistoYMax(); Float_t ymin = GetHistoYMin();
213 Int_t zbins = GetHistoZBins(); Float_t zmax = GetHistoZMax(); Float_t zmin = GetHistoZMin();
214 Int_t ssbins = GetHistoShowerShapeBins(); Float_t ssmax = GetHistoShowerShapeMax(); Float_t ssmin = GetHistoShowerShapeMin();
215 Int_t tdbins = GetHistoDiffTimeBins() ; Float_t tdmax = GetHistoDiffTimeMax(); Float_t tdmin = GetHistoDiffTimeMin();
217 Int_t nv0sbins = GetHistoV0SignalBins(); Int_t nv0smax = GetHistoV0SignalMax(); Int_t nv0smin = GetHistoV0SignalMin();
218 Int_t nv0mbins = GetHistoV0MultiplicityBins(); Int_t nv0mmax = GetHistoV0MultiplicityMax(); Int_t nv0mmin = GetHistoV0MultiplicityMin();
219 Int_t ntrmbins = GetHistoTrackMultiplicityBins(); Int_t ntrmmax = GetHistoTrackMultiplicityMax(); Int_t ntrmmin = GetHistoTrackMultiplicityMin();
226 if(fCalorimeter=="PHOS"){
232 fhE = new TH1F ("hE","E reconstructed clusters ", nptbins*5,ptmin,ptmax*5);
233 fhE->SetXTitle("E (GeV)");
234 outputContainer->Add(fhE);
237 fhPt = new TH1F ("hPt","p_{T} reconstructed clusters", nptbins,ptmin,ptmax);
238 fhPt->SetXTitle("p_{T} (GeV/c)");
239 outputContainer->Add(fhPt);
241 fhPhi = new TH1F ("hPhi","#phi reconstructed clusters ",nphibins,phimin,phimax);
242 fhPhi->SetXTitle("#phi (rad)");
243 outputContainer->Add(fhPhi);
245 fhEta = new TH1F ("hEta","#eta reconstructed clusters ",netabins,etamin,etamax);
246 fhEta->SetXTitle("#eta ");
247 outputContainer->Add(fhEta);
250 fhEtaPhiE = new TH3F ("hEtaPhiE","#eta vs #phi vs energy, reconstructed clusters",
251 netabins,etamin,etamax,nphibins,phimin,phimax,nptbins,ptmin,ptmax);
252 fhEtaPhiE->SetXTitle("#eta ");
253 fhEtaPhiE->SetYTitle("#phi (rad)");
254 fhEtaPhiE->SetZTitle("E (GeV) ");
255 outputContainer->Add(fhEtaPhiE);
257 fhClusterTimeEnergy = new TH2F ("hClusterTimeEnergy","energy vs TOF, reconstructed clusters",
258 nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
259 fhClusterTimeEnergy->SetXTitle("E (GeV) ");
260 fhClusterTimeEnergy->SetYTitle("TOF (ns)");
261 outputContainer->Add(fhClusterTimeEnergy);
263 fhClusterPairDiffTimeE = new TH2F("hClusterPairDiffTimeE","cluster pair time difference vs E, only good clusters",
264 nptbins,ptmin,ptmax, tdbins,tdmin,tdmax);
265 fhClusterPairDiffTimeE->SetXTitle("E_{cluster} (GeV)");
266 fhClusterPairDiffTimeE->SetYTitle("#Delta t (ns)");
267 outputContainer->Add(fhClusterPairDiffTimeE);
270 fhClusterMaxCellCloseCellRatio = new TH2F ("hClusterMaxCellCloseCellRatio","energy vs ratio of max cell / neighbour cell, reconstructed clusters",
271 nptbins,ptmin,ptmax, 100,0,1.);
272 fhClusterMaxCellCloseCellRatio->SetXTitle("E_{cluster} (GeV) ");
273 fhClusterMaxCellCloseCellRatio->SetYTitle("E_{cell i}/E_{cell max}");
274 outputContainer->Add(fhClusterMaxCellCloseCellRatio);
276 fhClusterMaxCellCloseCellDiff = new TH2F ("hClusterMaxCellCloseCellDiff","energy vs ratio of max cell / neighbour cell, reconstructed clusters",
277 nptbins,ptmin,ptmax, 500,0,100.);
278 fhClusterMaxCellCloseCellDiff->SetXTitle("E_{cluster} (GeV) ");
279 fhClusterMaxCellCloseCellDiff->SetYTitle("E_{cell max}-E_{cell i} (GeV)");
280 outputContainer->Add(fhClusterMaxCellCloseCellDiff);
282 fhClusterMaxCellDiff = new TH2F ("hClusterMaxCellDiff","energy vs difference of cluster energy - max cell energy / cluster energy, good clusters",
283 nptbins,ptmin,ptmax, 500,0,1.);
284 fhClusterMaxCellDiff->SetXTitle("E_{cluster} (GeV) ");
285 fhClusterMaxCellDiff->SetYTitle("(E_{cluster} - E_{cell max})/ E_{cluster}");
286 outputContainer->Add(fhClusterMaxCellDiff);
288 fhClusterMaxCellDiffNoCut = new TH2F ("hClusterMaxCellDiffNoCut","energy vs difference of cluster energy - max cell energy / cluster energy",
289 nptbins,ptmin,ptmax, 500,0,1.);
290 fhClusterMaxCellDiffNoCut->SetXTitle("E_{cluster} (GeV) ");
291 fhClusterMaxCellDiffNoCut->SetYTitle("(E_{cluster} - E_{cell max})/ E_{cluster}");
292 outputContainer->Add(fhClusterMaxCellDiffNoCut);
294 fhLambda0vsClusterMaxCellDiffE0 = new TH2F ("hLambda0vsClusterMaxCellDiffE0","shower shape, #lambda^{2}_{0} vs fraction of energy carried by max cell, E < 2 GeV ",
295 ssbins,ssmin,ssmax,500,0,1.);
296 fhLambda0vsClusterMaxCellDiffE0->SetYTitle("(E_{cluster} - E_{cell max})/ E_{cluster}");
297 fhLambda0vsClusterMaxCellDiffE0->SetXTitle("#lambda^{2}_{0}");
298 outputContainer->Add(fhLambda0vsClusterMaxCellDiffE0);
300 fhLambda0vsClusterMaxCellDiffE2 = new TH2F ("hLambda0vsClusterMaxCellDiffE2","shower shape, #lambda^{2}_{0} vs fraction of energy carried by max cell, 2 < E < 6 GeV ",
301 ssbins,ssmin,ssmax,500,0,1.);
302 fhLambda0vsClusterMaxCellDiffE2->SetYTitle("(E_{cluster} - E_{cell max})/ E_{cluster}");
303 fhLambda0vsClusterMaxCellDiffE2->SetXTitle("#lambda^{2}_{0}");
304 outputContainer->Add(fhLambda0vsClusterMaxCellDiffE2);
306 fhLambda0vsClusterMaxCellDiffE6 = new TH2F ("hLambda0vsClusterMaxCellDiffE6","shower shape, #lambda^{2}_{0} vs fraction of energy carried by max cell, E > 6 ",
307 ssbins,ssmin,ssmax,500,0,1.);
308 fhLambda0vsClusterMaxCellDiffE6->SetYTitle("(E_{cluster} - E_{cell max})/ E_{cluster}");
309 fhLambda0vsClusterMaxCellDiffE6->SetXTitle("#lambda^{2}_{0}");
310 outputContainer->Add(fhLambda0vsClusterMaxCellDiffE6);
312 fhNCellsvsClusterMaxCellDiffE0 = new TH2F ("hNCellsvsClusterMaxCellDiffE0","N cells per cluster vs fraction of energy carried by max cell, E < 2 GeV ",
313 nceclbins,nceclmin,nceclmax,500,0,1.);
314 fhNCellsvsClusterMaxCellDiffE0->SetYTitle("(E_{cluster} - E_{cell max})/ E_{cluster}");
315 fhNCellsvsClusterMaxCellDiffE0->SetXTitle("N cells per cluster");
316 outputContainer->Add(fhNCellsvsClusterMaxCellDiffE0);
318 fhNCellsvsClusterMaxCellDiffE2 = new TH2F ("hNCellsvsClusterMaxCellDiffE2","N cells per cluster vs fraction of energy carried by max cell, 2 < E < 6 GeV ",
319 nceclbins,nceclmin,nceclmax,500,0,1.);
320 fhNCellsvsClusterMaxCellDiffE2->SetYTitle("(E_{cluster} - E_{cell max})/ E_{cluster}");
321 fhNCellsvsClusterMaxCellDiffE2->SetXTitle("N cells per cluster");
322 outputContainer->Add(fhNCellsvsClusterMaxCellDiffE2);
324 fhNCellsvsClusterMaxCellDiffE6 = new TH2F ("hNCellsvsClusterMaxCellDiffE6","N cells per cluster vs fraction of energy carried by max cell, E > 6 ",
325 nceclbins,nceclmin,nceclmax,500,0,1.);
326 fhNCellsvsClusterMaxCellDiffE6->SetYTitle("(E_{cluster} - E_{cell max})/ E_{cluster}");
327 fhNCellsvsClusterMaxCellDiffE6->SetXTitle("N cells per cluster");
328 outputContainer->Add(fhNCellsvsClusterMaxCellDiffE6);
331 if(fCalorimeter=="EMCAL" && !GetCaloUtils()->GetEMCALRecoUtils()->IsRejectExoticCluster()){
333 fhBadClusterEnergy = new TH1F ("hBadClusterEnergy","Bad cluster energy", nptbins,ptmin,ptmax);
334 fhBadClusterEnergy->SetXTitle("E_{cluster} (GeV) ");
335 outputContainer->Add(fhBadClusterEnergy);
337 fhBadClusterMaxCellCloseCellRatio = new TH2F ("hBadClusterMaxCellCloseCellRatio","energy vs ratio of max cell / neighbour cell constributing cell, reconstructed bad clusters",
338 nptbins,ptmin,ptmax, 100,0,1.);
339 fhBadClusterMaxCellCloseCellRatio->SetXTitle("E_{cluster} (GeV) ");
340 fhBadClusterMaxCellCloseCellRatio->SetYTitle("ratio");
341 outputContainer->Add(fhBadClusterMaxCellCloseCellRatio);
343 fhBadClusterMaxCellCloseCellDiff = new TH2F ("hBadClusterMaxCellCloseCellDiff","energy vs ratio of max cell - neighbour cell constributing cell, reconstructed bad clusters",
344 nptbins,ptmin,ptmax, 500,0,100);
345 fhBadClusterMaxCellCloseCellDiff->SetXTitle("E_{cluster} (GeV) ");
346 fhBadClusterMaxCellCloseCellDiff->SetYTitle("E_{cell max} - E_{cell i} (GeV)");
347 outputContainer->Add(fhBadClusterMaxCellCloseCellDiff);
349 fhBadClusterMaxCellDiff = new TH2F ("hBadClusterMaxCellDiff","energy vs difference of cluster energy - max cell energy / cluster energy for bad clusters",
350 nptbins,ptmin,ptmax, 500,0,1.);
351 fhBadClusterMaxCellDiff->SetXTitle("E_{cluster} (GeV) ");
352 fhBadClusterMaxCellDiff->SetYTitle("(E_{cluster} - E_{cell max}) / E_{cluster}");
353 outputContainer->Add(fhBadClusterMaxCellDiff);
355 fhBadClusterTimeEnergy = new TH2F ("hBadClusterTimeEnergy","energy vs TOF of reconstructed bad clusters",
356 nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
357 fhBadClusterTimeEnergy->SetXTitle("E_{cluster} (GeV) ");
358 fhBadClusterTimeEnergy->SetYTitle("TOF (ns)");
359 outputContainer->Add(fhBadClusterTimeEnergy);
361 fhBadClusterPairDiffTimeE = new TH2F("hBadClusterPairDiffTimeE","cluster pair time difference (bad - good) vs E from bad cluster",nptbins,ptmin,ptmax, tdbins,tdmin,tdmax);
362 fhBadClusterPairDiffTimeE->SetXTitle("E_{bad cluster} (GeV)");
363 fhBadClusterPairDiffTimeE->SetYTitle("#Delta t (ns)");
364 outputContainer->Add(fhBadClusterPairDiffTimeE);
366 fhBadClusterL0 = new TH2F ("hBadClusterL0","shower shape, #lambda^{2}_{0} vs E for bad cluster ",
367 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
368 fhBadClusterL0->SetXTitle("E_{cluster}");
369 fhBadClusterL0->SetYTitle("#lambda^{2}_{0}");
370 outputContainer->Add(fhBadClusterL0);
372 fhBadClusterL1 = new TH2F ("hBadClusterL1","shower shape, #lambda^{2}_{1} vs E for bad cluster ",
373 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
374 fhBadClusterL1->SetXTitle("E_{cluster}");
375 fhBadClusterL1->SetYTitle("#lambda^{2}_{1}");
376 outputContainer->Add(fhBadClusterL1);
378 fhBadClusterD = new TH2F ("hBadClusterD","shower shape, Dispersion^{2} vs E for bad cluster ",
379 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
380 fhBadClusterD->SetXTitle("E_{cluster}");
381 fhBadClusterD->SetYTitle("Dispersion");
382 outputContainer->Add(fhBadClusterD);
384 if(GetReader()->GetDataType()==AliCaloTrackReader::kESD) {
385 fhBadCellTimeSpreadRespectToCellMax = new TH2F ("hBadCellTimeSpreadRespectToCellMax","t_{cell max}-t_{cell i} per cluster", nptbins,ptmin,ptmax, 250,-500,500);
386 fhBadCellTimeSpreadRespectToCellMax->SetXTitle("E (GeV)");
387 fhBadCellTimeSpreadRespectToCellMax->SetYTitle("#Delta t_{cell max - i} (ns)");
388 outputContainer->Add(fhBadCellTimeSpreadRespectToCellMax);
390 fhBadClusterMaxCellDiffAverageTime = new TH2F ("hBadClusterMaxCellDiffAverageTime","t_{cell max}-t_{average} per cluster", nptbins,ptmin,ptmax, 250,-500,500);
391 fhBadClusterMaxCellDiffAverageTime->SetXTitle("E (GeV)");
392 fhBadClusterMaxCellDiffAverageTime->SetYTitle("#Delta t_{cell max - average} (ns)");
393 outputContainer->Add(fhBadClusterMaxCellDiffAverageTime);
395 fhBadClusterMaxCellDiffWeightTime = new TH2F ("hBadClusterMaxCellDiffWeightTime","t_{cell max weighted}-t_{average weighted} per cluster", nptbins,ptmin,ptmax, 250,-500,500);
396 fhBadClusterMaxCellDiffWeightTime->SetXTitle("E (GeV)");
397 fhBadClusterMaxCellDiffWeightTime->SetYTitle("#Delta t_{cell max - average weighted} (ns)");
398 outputContainer->Add(fhBadClusterMaxCellDiffWeightTime);
400 fhBadClusterDiffWeightAverTime = new TH2F ("hBadClusterDiffWeightAverTime","Average Time in cluster - Weighted time in cluster", nptbins,ptmin,ptmax, 250,-500,500);
401 fhBadClusterDiffWeightAverTime->SetXTitle("E (GeV)");
402 fhBadClusterDiffWeightAverTime->SetYTitle("#bar{t}-wt");
403 outputContainer->Add(fhBadClusterDiffWeightAverTime);
405 fhBadClusterMaxCellDiffAverageNoMaxTime = new TH2F ("hBadClusterMaxCellDiffAverageNoMaxTime","t_{cell max}-t_{average} per cluster", nptbins,ptmin,ptmax, 250,-500,500);
406 fhBadClusterMaxCellDiffAverageNoMaxTime->SetXTitle("E (GeV)");
407 fhBadClusterMaxCellDiffAverageNoMaxTime->SetYTitle("#Delta t_{cell max - average} (ns)");
408 outputContainer->Add(fhBadClusterMaxCellDiffAverageNoMaxTime);
410 fhBadClusterMaxCellDiffWeightNoMaxTime = new TH2F ("hBadClusterMaxCellDiffWeightNoMaxTime","t_{cell max weighted}-t_{average weighted} per cluster", nptbins,ptmin,ptmax, 250,-500,500);
411 fhBadClusterMaxCellDiffWeightNoMaxTime->SetXTitle("E (GeV)");
412 fhBadClusterMaxCellDiffWeightNoMaxTime->SetYTitle("#Delta t_{cell max - average weighted} (ns)");
413 outputContainer->Add(fhBadClusterMaxCellDiffWeightNoMaxTime);
417 fhBadClusterNoMaxCellWeight = new TH2F ("hBadClusterNoMaxCellWeight"," #S weight_{no max} / #S weight_{tot} per cluster", nptbins,ptmin,ptmax, 100,0,1);
418 fhBadClusterNoMaxCellWeight->SetXTitle("E (GeV)");
419 fhBadClusterNoMaxCellWeight->SetYTitle("#S weight_{no max} / #S weight_{tot}");
420 outputContainer->Add(fhBadClusterNoMaxCellWeight);
424 // Cluster size in terms of cells
426 fhDeltaIEtaDeltaIPhiE0[0] = new TH2F ("hDeltaIEtaDeltaIPhiE0"," Cluster size in columns vs rows for E < 2 GeV, n cells > 3",
428 fhDeltaIEtaDeltaIPhiE0[0]->SetXTitle("#Delta Column");
429 fhDeltaIEtaDeltaIPhiE0[0]->SetYTitle("#Delta Row");
430 outputContainer->Add(fhDeltaIEtaDeltaIPhiE0[0]);
432 fhDeltaIEtaDeltaIPhiE2[0] = new TH2F ("hDeltaIEtaDeltaIPhiE2"," Cluster size in columns vs rows for 2 <E < 6 GeV, n cells > 3",
434 fhDeltaIEtaDeltaIPhiE2[0]->SetXTitle("#Delta Column");
435 fhDeltaIEtaDeltaIPhiE2[0]->SetYTitle("#Delta Row");
436 outputContainer->Add(fhDeltaIEtaDeltaIPhiE2[0]);
438 fhDeltaIEtaDeltaIPhiE6[0] = new TH2F ("hDeltaIEtaDeltaIPhiE6"," Cluster size in columns vs rows for E > 6 GeV, n cells > 3",
440 fhDeltaIEtaDeltaIPhiE6[0]->SetXTitle("#Delta Column");
441 fhDeltaIEtaDeltaIPhiE6[0]->SetYTitle("#Delta Row");
442 outputContainer->Add(fhDeltaIEtaDeltaIPhiE6[0]);
444 fhDeltaIA[0] = new TH2F ("hDeltaIA"," Cluster *asymmetry* in cell units vs E",
445 nptbins,ptmin,ptmax,21,-1.05,1.05);
446 fhDeltaIA[0]->SetXTitle("E_{cluster}");
447 fhDeltaIA[0]->SetYTitle("A_{cell in cluster}");
448 outputContainer->Add(fhDeltaIA[0]);
450 fhDeltaIAL0[0] = new TH2F ("hDeltaIAL0"," Cluster *asymmetry* in cell units vs #lambda^{2}_{0}",
451 ssbins,ssmin,ssmax,21,-1.05,1.05);
452 fhDeltaIAL0[0]->SetXTitle("#lambda^{2}_{0}");
453 fhDeltaIAL0[0]->SetYTitle("A_{cell in cluster}");
454 outputContainer->Add(fhDeltaIAL0[0]);
456 fhDeltaIAL1[0] = new TH2F ("hDeltaIAL1"," Cluster *asymmetry* in cell units vs #lambda^{2}_{1}",
457 ssbins,ssmin,ssmax,21,-1.05,1.05);
458 fhDeltaIAL1[0]->SetXTitle("#lambda^{2}_{1}");
459 fhDeltaIAL1[0]->SetYTitle("A_{cell in cluster}");
460 outputContainer->Add(fhDeltaIAL1[0]);
462 fhDeltaIANCells[0] = new TH2F ("hDeltaIANCells"," Cluster *asymmetry* in cell units vs N cells in cluster",
463 nceclbins,nceclmin,nceclmax,21,-1.05,1.05);
464 fhDeltaIANCells[0]->SetXTitle("N_{cell in cluster}");
465 fhDeltaIANCells[0]->SetYTitle("A_{cell in cluster}");
466 outputContainer->Add(fhDeltaIANCells[0]);
469 fhDeltaIEtaDeltaIPhiE0[1] = new TH2F ("hDeltaIEtaDeltaIPhiE0Charged"," Cluster size in columns vs rows for E < 2 GeV, n cells > 3, matched with track",
471 fhDeltaIEtaDeltaIPhiE0[1]->SetXTitle("#Delta Column");
472 fhDeltaIEtaDeltaIPhiE0[1]->SetYTitle("#Delta Row");
473 outputContainer->Add(fhDeltaIEtaDeltaIPhiE0[1]);
475 fhDeltaIEtaDeltaIPhiE2[1] = new TH2F ("hDeltaIEtaDeltaIPhiE2Charged"," Cluster size in columns vs rows for 2 <E < 6 GeV, n cells > 3, matched with track",
477 fhDeltaIEtaDeltaIPhiE2[1]->SetXTitle("#Delta Column");
478 fhDeltaIEtaDeltaIPhiE2[1]->SetYTitle("#Delta Row");
479 outputContainer->Add(fhDeltaIEtaDeltaIPhiE2[1]);
481 fhDeltaIEtaDeltaIPhiE6[1] = new TH2F ("hDeltaIEtaDeltaIPhiE6Charged"," Cluster size in columns vs rows for E > 6 GeV, n cells > 3, matched with track",
483 fhDeltaIEtaDeltaIPhiE6[1]->SetXTitle("#Delta Column");
484 fhDeltaIEtaDeltaIPhiE6[1]->SetYTitle("#Delta Row");
485 outputContainer->Add(fhDeltaIEtaDeltaIPhiE6[1]);
487 fhDeltaIA[1] = new TH2F ("hDeltaIACharged"," Cluster *asymmetry* in cell units vs E, matched with track",
488 nptbins,ptmin,ptmax,21,-1.05,1.05);
489 fhDeltaIA[1]->SetXTitle("E_{cluster}");
490 fhDeltaIA[1]->SetYTitle("A_{cell in cluster}");
491 outputContainer->Add(fhDeltaIA[1]);
493 fhDeltaIAL0[1] = new TH2F ("hDeltaIAL0Charged"," Cluster *asymmetry* in cell units vs #lambda^{2}_{0}, matched with track",
494 ssbins,ssmin,ssmax,21,-1.05,1.05);
495 fhDeltaIAL0[1]->SetXTitle("#lambda^{2}_{0}");
496 fhDeltaIAL0[1]->SetYTitle("A_{cell in cluster}");
497 outputContainer->Add(fhDeltaIAL0[1]);
499 fhDeltaIAL1[1] = new TH2F ("hDeltaIAL1Charged"," Cluster *asymmetry* in cell units vs #lambda^{2}_{1}, matched with track",
500 ssbins,ssmin,ssmax,21,-1.05,1.05);
501 fhDeltaIAL1[1]->SetXTitle("#lambda^{2}_{1}");
502 fhDeltaIAL1[1]->SetYTitle("A_{cell in cluster}");
503 outputContainer->Add(fhDeltaIAL1[1]);
505 fhDeltaIANCells[1] = new TH2F ("hDeltaIANCellsCharged"," Cluster *asymmetry* in cell units vs N cells in cluster, matched with track",
506 nceclbins,nceclmin,nceclmax,21,-1.05,1.05);
507 fhDeltaIANCells[1]->SetXTitle("N_{cell in cluster}");
508 fhDeltaIANCells[1]->SetYTitle("A_{cell in cluster}");
509 outputContainer->Add(fhDeltaIANCells[1]);
512 TString particle[]={"Photon","Electron","Conversion","Hadron"};
513 for (Int_t iPart = 0; iPart < 4; iPart++) {
515 fhDeltaIAMC[iPart] = new TH2F (Form("hDeltaIA_MC%s",particle[iPart].Data()),Form(" Cluster *asymmetry* in cell units vs E, from %s",particle[iPart].Data()),
516 nptbins,ptmin,ptmax,21,-1.05,1.05);
517 fhDeltaIAMC[iPart]->SetXTitle("E_{cluster}");
518 fhDeltaIAMC[iPart]->SetYTitle("A_{cell in cluster}");
519 outputContainer->Add(fhDeltaIAMC[iPart]);
526 fhECharged = new TH1F ("hECharged","E reconstructed clusters, matched with track", nptbins,ptmin,ptmax);
527 fhECharged->SetXTitle("E (GeV)");
528 outputContainer->Add(fhECharged);
530 fhPtCharged = new TH1F ("hPtCharged","p_{T} reconstructed clusters, matched with track", nptbins,ptmin,ptmax);
531 fhPtCharged->SetXTitle("p_{T} (GeV/c)");
532 outputContainer->Add(fhPtCharged);
534 fhPhiCharged = new TH1F ("hPhiCharged","#phi reconstructed clusters, matched with track",nphibins,phimin,phimax);
535 fhPhiCharged->SetXTitle("#phi (rad)");
536 outputContainer->Add(fhPhiCharged);
538 fhEtaCharged = new TH1F ("hEtaCharged","#eta reconstructed clusters, matched with track",netabins,etamin,etamax);
539 fhEtaCharged->SetXTitle("#eta ");
540 outputContainer->Add(fhEtaCharged);
543 fhEtaPhiECharged = new TH3F ("hEtaPhiECharged","#eta vs #phi, reconstructed clusters, matched with track",
544 netabins,etamin,etamax,nphibins,phimin,phimax,nptbins,ptmin,ptmax);
545 fhEtaPhiECharged->SetXTitle("#eta ");
546 fhEtaPhiECharged->SetYTitle("#phi ");
547 fhEtaPhiECharged->SetZTitle("E (GeV) ");
548 outputContainer->Add(fhEtaPhiECharged);
551 fh1pOverE = new TH2F("h1pOverE","TRACK matches p/E",nptbins,ptmin,ptmax, nPoverEbins,pOverEmin,pOverEmax);
552 fh1pOverE->SetYTitle("p/E");
553 fh1pOverE->SetXTitle("p_{T} (GeV/c)");
554 outputContainer->Add(fh1pOverE);
556 fh1dR = new TH1F("h1dR","TRACK matches dR",ndRbins,dRmin,dRmax);
557 fh1dR->SetXTitle("#Delta R (rad)");
558 outputContainer->Add(fh1dR) ;
560 fh2MatchdEdx = new TH2F("h2MatchdEdx","dE/dx vs. p for all matches",nptbins,ptmin,ptmax,ndedxbins,dedxmin,dedxmax);
561 fh2MatchdEdx->SetXTitle("p (GeV/c)");
562 fh2MatchdEdx->SetYTitle("<dE/dx>");
563 outputContainer->Add(fh2MatchdEdx);
565 fh2EledEdx = new TH2F("h2EledEdx","dE/dx vs. p for electrons",nptbins,ptmin,ptmax,ndedxbins,dedxmin,dedxmax);
566 fh2EledEdx->SetXTitle("p (GeV/c)");
567 fh2EledEdx->SetYTitle("<dE/dx>");
568 outputContainer->Add(fh2EledEdx) ;
570 fh1pOverER02 = new TH2F("h1pOverER02","TRACK matches p/E, all",nptbins,ptmin,ptmax, nPoverEbins,pOverEmin,pOverEmax);
571 fh1pOverER02->SetYTitle("p/E");
572 fh1pOverER02->SetXTitle("p_{T} (GeV/c)");
573 outputContainer->Add(fh1pOverER02);
576 if(fFillAllPi0Histo){
577 fhIM = new TH2F ("hIM","Cluster pairs Invariant mass vs reconstructed pair energy, ncell > 1",nptbins,ptmin,ptmax,nmassbins,massmin,massmax);
578 fhIM->SetXTitle("p_{T, cluster pairs} (GeV) ");
579 fhIM->SetYTitle("M_{cluster pairs} (GeV/c^{2})");
580 outputContainer->Add(fhIM);
582 fhAsym = new TH2F ("hAssym","Cluster pairs Asymmetry vs reconstructed pair energy",nptbins,ptmin,ptmax,nasymbins,asymmin,asymmax);
583 fhAsym->SetXTitle("p_{T, cluster pairs} (GeV) ");
584 fhAsym->SetYTitle("Asymmetry");
585 outputContainer->Add(fhAsym);
589 fhNCellsPerClusterNoCut = new TH2F ("hNCellsPerClusterNoCut","# cells per cluster vs energy vs #eta, no bad clusters cut",
590 nptbins,ptmin,ptmax, nceclbins,nceclmin,nceclmax);
591 fhNCellsPerClusterNoCut->SetXTitle("E (GeV)");
592 fhNCellsPerClusterNoCut->SetYTitle("n cells");
593 outputContainer->Add(fhNCellsPerClusterNoCut);
595 fhNCellsPerCluster = new TH2F ("hNCellsPerCluster","# cells per cluster vs energy vs #eta",nptbins,ptmin,ptmax, nceclbins,nceclmin,nceclmax);
596 fhNCellsPerCluster->SetXTitle("E (GeV)");
597 fhNCellsPerCluster->SetYTitle("n cells");
598 outputContainer->Add(fhNCellsPerCluster);
600 if((fCalorimeter=="EMCAL" && GetReader()->GetEMCALPtMin() < 0.3) ||
601 (fCalorimeter=="PHOS" && GetReader()->GetPHOSPtMin() < 0.3)) {
602 fhNCellsPerClusterMIP = new TH2F ("hNCellsPerClusterMIP","# cells per cluster vs energy vs #eta, smaller bin for MIP search",
604 fhNCellsPerClusterMIP->SetXTitle("E (GeV)");
605 fhNCellsPerClusterMIP->SetYTitle("n cells");
606 outputContainer->Add(fhNCellsPerClusterMIP);
610 fhNCellsPerClusterMIPCharged = new TH2F ("hNCellsPerClusterMIPCharged","# cells per track-matched cluster vs energy vs #eta, smaller bin for MIP search",
612 fhNCellsPerClusterMIPCharged->SetXTitle("E (GeV)");
613 fhNCellsPerClusterMIPCharged->SetYTitle("n cells");
614 outputContainer->Add(fhNCellsPerClusterMIPCharged);
618 fhNClusters = new TH1F ("hNClusters","# clusters", nclbins,nclmin,nclmax);
619 fhNClusters->SetXTitle("number of clusters");
620 outputContainer->Add(fhNClusters);
622 if(fFillAllPosHisto2){
625 fhXYZ = new TH3F ("hXYZ","Cluster: x vs y vs z",xbins,xmin,xmax,ybins,ymin,ymax,zbins,zmin,zmax);
626 fhXYZ->SetXTitle("x (cm)");
627 fhXYZ->SetYTitle("y (cm)");
628 fhXYZ->SetZTitle("z (cm) ");
629 outputContainer->Add(fhXYZ);
632 fhXNCells = new TH2F ("hXNCells","Cluster X position vs N Cells per Cluster",xbins,xmin,xmax,nceclbins,nceclmin,nceclmax);
633 fhXNCells->SetXTitle("x (cm)");
634 fhXNCells->SetYTitle("N cells per cluster");
635 outputContainer->Add(fhXNCells);
637 fhZNCells = new TH2F ("hZNCells","Cluster Z position vs N Cells per Cluster",zbins,zmin,zmax,nceclbins,nceclmin,nceclmax);
638 fhZNCells->SetXTitle("z (cm)");
639 fhZNCells->SetYTitle("N cells per cluster");
640 outputContainer->Add(fhZNCells);
642 fhXE = new TH2F ("hXE","Cluster X position vs cluster energy",xbins,xmin,xmax,nptbins,ptmin,ptmax);
643 fhXE->SetXTitle("x (cm)");
644 fhXE->SetYTitle("E (GeV)");
645 outputContainer->Add(fhXE);
647 fhZE = new TH2F ("hZE","Cluster Z position vs cluster energy",zbins,zmin,zmax,nptbins,ptmin,ptmax);
648 fhZE->SetXTitle("z (cm)");
649 fhZE->SetYTitle("E (GeV)");
650 outputContainer->Add(fhZE);
652 fhRNCells = new TH2F ("hRNCells","Cluster R position vs N Cells per Cluster",rbins,rmin,rmax,nceclbins,nceclmin,nceclmax);
653 fhRNCells->SetXTitle("r = #sqrt{x^{2}+y^{2}} (cm)");
654 fhRNCells->SetYTitle("N cells per cluster");
655 outputContainer->Add(fhRNCells);
658 fhYNCells = new TH2F ("hYNCells","Cluster Y position vs N Cells per Cluster",ybins,ymin,ymax,nceclbins,nceclmin,nceclmax);
659 fhYNCells->SetXTitle("y (cm)");
660 fhYNCells->SetYTitle("N cells per cluster");
661 outputContainer->Add(fhYNCells);
663 fhRE = new TH2F ("hRE","Cluster R position vs cluster energy",rbins,rmin,rmax,nptbins,ptmin,ptmax);
664 fhRE->SetXTitle("r = #sqrt{x^{2}+y^{2}} (cm)");
665 fhRE->SetYTitle("E (GeV)");
666 outputContainer->Add(fhRE);
668 fhYE = new TH2F ("hYE","Cluster Y position vs cluster energy",ybins,ymin,ymax,nptbins,ptmin,ptmax);
669 fhYE->SetXTitle("y (cm)");
670 fhYE->SetYTitle("E (GeV)");
671 outputContainer->Add(fhYE);
673 if(fFillAllPosHisto){
675 fhRCellE = new TH2F ("hRCellE","Cell R position vs cell energy",rbins,rmin,rmax,nptbins,ptmin,ptmax);
676 fhRCellE->SetXTitle("r = #sqrt{x^{2}+y^{2}} (cm)");
677 fhRCellE->SetYTitle("E (GeV)");
678 outputContainer->Add(fhRCellE);
680 fhXCellE = new TH2F ("hXCellE","Cell X position vs cell energy",xbins,xmin,xmax,nptbins,ptmin,ptmax);
681 fhXCellE->SetXTitle("x (cm)");
682 fhXCellE->SetYTitle("E (GeV)");
683 outputContainer->Add(fhXCellE);
685 fhYCellE = new TH2F ("hYCellE","Cell Y position vs cell energy",ybins,ymin,ymax,nptbins,ptmin,ptmax);
686 fhYCellE->SetXTitle("y (cm)");
687 fhYCellE->SetYTitle("E (GeV)");
688 outputContainer->Add(fhYCellE);
690 fhZCellE = new TH2F ("hZCellE","Cell Z position vs cell energy",zbins,zmin,zmax,nptbins,ptmin,ptmax);
691 fhZCellE->SetXTitle("z (cm)");
692 fhZCellE->SetYTitle("E (GeV)");
693 outputContainer->Add(fhZCellE);
695 fhXYZCell = new TH3F ("hXYZCell","Cell : x vs y vs z",xbins,xmin,xmax,ybins,ymin,ymax,zbins,zmin,zmax);
696 fhXYZCell->SetXTitle("x (cm)");
697 fhXYZCell->SetYTitle("y (cm)");
698 fhXYZCell->SetZTitle("z (cm)");
699 outputContainer->Add(fhXYZCell);
702 Float_t dx = TMath::Abs(xmin)+TMath::Abs(xmax);
703 Float_t dy = TMath::Abs(ymin)+TMath::Abs(ymax);
704 Float_t dz = TMath::Abs(zmin)+TMath::Abs(zmax);
705 Float_t dr = TMath::Abs(rmin)+TMath::Abs(rmax);
707 fhDeltaCellClusterRNCells = new TH2F ("hDeltaCellClusterRNCells","Cluster-Cell R position vs N Cells per Cluster",rbins*2,-dr,dr,nceclbins,nceclmin,nceclmax);
708 fhDeltaCellClusterRNCells->SetXTitle("r = #sqrt{x^{2}+y^{2}} (cm)");
709 fhDeltaCellClusterRNCells->SetYTitle("N cells per cluster");
710 outputContainer->Add(fhDeltaCellClusterRNCells);
712 fhDeltaCellClusterXNCells = new TH2F ("hDeltaCellClusterXNCells","Cluster-Cell X position vs N Cells per Cluster",xbins*2,-dx,dx,nceclbins,nceclmin,nceclmax);
713 fhDeltaCellClusterXNCells->SetXTitle("x (cm)");
714 fhDeltaCellClusterXNCells->SetYTitle("N cells per cluster");
715 outputContainer->Add(fhDeltaCellClusterXNCells);
717 fhDeltaCellClusterYNCells = new TH2F ("hDeltaCellClusterYNCells","Cluster-Cell Y position vs N Cells per Cluster",ybins*2,-dy,dy,nceclbins,nceclmin,nceclmax);
718 fhDeltaCellClusterYNCells->SetXTitle("y (cm)");
719 fhDeltaCellClusterYNCells->SetYTitle("N cells per cluster");
720 outputContainer->Add(fhDeltaCellClusterYNCells);
722 fhDeltaCellClusterZNCells = new TH2F ("hDeltaCellClusterZNCells","Cluster-Cell Z position vs N Cells per Cluster",zbins*2,-dz,dz,nceclbins,nceclmin,nceclmax);
723 fhDeltaCellClusterZNCells->SetXTitle("z (cm)");
724 fhDeltaCellClusterZNCells->SetYTitle("N cells per cluster");
725 outputContainer->Add(fhDeltaCellClusterZNCells);
727 fhDeltaCellClusterRE = new TH2F ("hDeltaCellClusterRE","Cluster-Cell R position vs cluster energy",rbins*2,-dr,dr,nptbins,ptmin,ptmax);
728 fhDeltaCellClusterRE->SetXTitle("r = #sqrt{x^{2}+y^{2}} (cm)");
729 fhDeltaCellClusterRE->SetYTitle("E (GeV)");
730 outputContainer->Add(fhDeltaCellClusterRE);
732 fhDeltaCellClusterXE = new TH2F ("hDeltaCellClusterXE","Cluster-Cell X position vs cluster energy",xbins*2,-dx,dx,nptbins,ptmin,ptmax);
733 fhDeltaCellClusterXE->SetXTitle("x (cm)");
734 fhDeltaCellClusterXE->SetYTitle("E (GeV)");
735 outputContainer->Add(fhDeltaCellClusterXE);
737 fhDeltaCellClusterYE = new TH2F ("hDeltaCellClusterYE","Cluster-Cell Y position vs cluster energy",ybins*2,-dy,dy,nptbins,ptmin,ptmax);
738 fhDeltaCellClusterYE->SetXTitle("y (cm)");
739 fhDeltaCellClusterYE->SetYTitle("E (GeV)");
740 outputContainer->Add(fhDeltaCellClusterYE);
742 fhDeltaCellClusterZE = new TH2F ("hDeltaCellClusterZE","Cluster-Cell Z position vs cluster energy",zbins*2,-dz,dz,nptbins,ptmin,ptmax);
743 fhDeltaCellClusterZE->SetXTitle("z (cm)");
744 fhDeltaCellClusterZE->SetYTitle("E (GeV)");
745 outputContainer->Add(fhDeltaCellClusterZE);
747 fhEtaPhiAmp = new TH3F ("hEtaPhiAmp","Cell #eta vs cell #phi vs cell energy",netabins,etamin,etamax,nphibins,phimin,phimax,nptbins,ptmin,ptmax);
748 fhEtaPhiAmp->SetXTitle("#eta ");
749 fhEtaPhiAmp->SetYTitle("#phi (rad)");
750 fhEtaPhiAmp->SetZTitle("E (GeV) ");
751 outputContainer->Add(fhEtaPhiAmp);
756 fhNCells = new TH1F ("hNCells","# cells", ncebins,ncemin,ncemax);
757 fhNCells->SetXTitle("n cells");
758 outputContainer->Add(fhNCells);
760 fhAmplitude = new TH1F ("hAmplitude","Cell Energy", nptbins*2,ptmin,ptmax);
761 fhAmplitude->SetXTitle("Cell Energy (GeV)");
762 outputContainer->Add(fhAmplitude);
764 fhAmpId = new TH2F ("hAmpId","Cell Energy", nfineptbins,ptfinemin,ptfinemax,fNMaxRows*fNMaxCols*fNModules,0,fNMaxRows*fNMaxCols*fNModules);
765 fhAmpId->SetXTitle("Cell Energy (GeV)");
766 outputContainer->Add(fhAmpId);
768 //Cell Time histograms, time only available in ESDs
769 if(GetReader()->GetDataType()==AliCaloTrackReader::kESD) {
771 fhCellTimeSpreadRespectToCellMax = new TH2F ("hCellTimeSpreadRespectToCellMax","t_{cell max}-t_{cell i} per cluster", nptbins,ptmin,ptmax, 250,-500,500);
772 fhCellTimeSpreadRespectToCellMax->SetXTitle("E (GeV)");
773 fhCellTimeSpreadRespectToCellMax->SetYTitle("#Delta t_{cell max-i} (ns)");
774 outputContainer->Add(fhCellTimeSpreadRespectToCellMax);
776 fhClusterMaxCellDiffAverageTime = new TH2F ("hClusterMaxCellDiffAverageTime","t_{cell max}-t_{average} per cluster", nptbins,ptmin,ptmax, 250,-500,500);
777 fhClusterMaxCellDiffAverageTime->SetXTitle("E (GeV)");
778 fhClusterMaxCellDiffAverageTime->SetYTitle("#Delta t_{cell max - average} (ns)");
779 outputContainer->Add(fhClusterMaxCellDiffAverageTime);
781 fhClusterMaxCellDiffWeightTime = new TH2F ("hClusterMaxCellDiffWeightTime","t_{cell max weighted}-t_{average weighted} per cluster", nptbins,ptmin,ptmax, 250,-500,500);
782 fhClusterMaxCellDiffWeightTime->SetXTitle("E (GeV)");
783 fhClusterMaxCellDiffWeightTime->SetYTitle("#Delta t_{cell max - average weighted} (ns)");
784 outputContainer->Add(fhClusterMaxCellDiffWeightTime);
786 fhClusterMaxCellDiffAverageNoMaxTime = new TH2F ("hClusterMaxCellDiffAverageNoMaxTime","t_{cell max}-t_{average} per cluster", nptbins,ptmin,ptmax, 250,-500,500);
787 fhClusterMaxCellDiffAverageNoMaxTime->SetXTitle("E (GeV)");
788 fhClusterMaxCellDiffAverageNoMaxTime->SetYTitle("#Delta t_{cell max - average} (ns)");
789 outputContainer->Add(fhClusterMaxCellDiffAverageNoMaxTime);
791 fhClusterMaxCellDiffWeightNoMaxTime = new TH2F ("hClusterMaxCellDiffWeightNoMaxTime","t_{cell max weighted}-t_{average weighted} per cluster", nptbins,ptmin,ptmax, 250,-500,500);
792 fhClusterMaxCellDiffWeightNoMaxTime->SetXTitle("E (GeV)");
793 fhClusterMaxCellDiffWeightNoMaxTime->SetYTitle("#Delta t_{cell max - average weighted} (ns)");
794 outputContainer->Add(fhClusterMaxCellDiffWeightNoMaxTime);
796 fhClusterDiffWeightAverTime = new TH2F ("hClusterDiffWeightAverTime","Average Time in cluster - Weighted time in cluster", nptbins,ptmin,ptmax, 250,-500,500);
797 fhClusterDiffWeightAverTime->SetXTitle("E (GeV)");
798 fhClusterDiffWeightAverTime->SetYTitle("#bar{t} - wt");
799 outputContainer->Add(fhClusterDiffWeightAverTime);
801 fhCellIdCellLargeTimeSpread= new TH1F ("hCellIdCellLargeTimeSpread","Cells with time 100 ns larger than cell max in cluster ",
802 fNMaxCols*fNMaxRows*fNModules,0,fNMaxCols*fNMaxRows*fNModules);
803 fhCellIdCellLargeTimeSpread->SetXTitle("Absolute Cell Id");
804 outputContainer->Add(fhCellIdCellLargeTimeSpread);
806 fhTime = new TH1F ("hTime","Cell Time",ntimebins,timemin,timemax);
807 fhTime->SetXTitle("Cell Time (ns)");
808 outputContainer->Add(fhTime);
810 fhTimeId = new TH2F ("hTimeId","Cell Time vs Absolute Id",
811 ntimebins,timemin,timemax,fNMaxRows*fNMaxCols*fNModules,0,fNMaxRows*fNMaxCols*fNModules);
812 fhTimeId->SetXTitle("Cell Time (ns)");
813 fhTimeId->SetYTitle("Cell Absolute Id");
814 outputContainer->Add(fhTimeId);
816 fhTimeAmp = new TH2F ("hTimeAmp","Cell Time vs Cell Energy",nptbins*2,ptmin,ptmax,ntimebins,timemin,timemax);
817 fhTimeAmp->SetYTitle("Cell Time (ns)");
818 fhTimeAmp->SetXTitle("Cell Energy (GeV)");
819 outputContainer->Add(fhTimeAmp);
824 fhClusterNoMaxCellWeight = new TH2F ("hClusterNoMaxCellWeight"," #S weight_{no max} / #S weight_{tot} per cluster", nptbins,ptmin,ptmax, 100,0,1);
825 fhClusterNoMaxCellWeight->SetXTitle("E (GeV)");
826 fhClusterNoMaxCellWeight->SetYTitle("#S weight_{no max} / #S weight_{tot}");
827 outputContainer->Add(fhClusterNoMaxCellWeight);
831 fhCaloCorrNClusters = new TH2F ("hCaloCorrNClusters","# clusters in EMCAL vs PHOS", nclbins,nclmin,nclmax,nclbins,nclmin,nclmax);
832 fhCaloCorrNClusters->SetXTitle("number of clusters in EMCAL");
833 fhCaloCorrNClusters->SetYTitle("number of clusters in PHOS");
834 outputContainer->Add(fhCaloCorrNClusters);
836 fhCaloCorrEClusters = new TH2F ("hCaloCorrEClusters","summed energy of clusters in EMCAL vs PHOS", nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
837 fhCaloCorrEClusters->SetXTitle("#Sigma E of clusters in EMCAL (GeV)");
838 fhCaloCorrEClusters->SetYTitle("#Sigma E of clusters in PHOS (GeV)");
839 outputContainer->Add(fhCaloCorrEClusters);
841 fhCaloCorrNCells = new TH2F ("hCaloCorrNCells","# Cells in EMCAL vs PHOS", ncebins,ncemin,ncemax, ncebins,ncemin,ncemax);
842 fhCaloCorrNCells->SetXTitle("number of Cells in EMCAL");
843 fhCaloCorrNCells->SetYTitle("number of Cells in PHOS");
844 outputContainer->Add(fhCaloCorrNCells);
846 fhCaloCorrECells = new TH2F ("hCaloCorrECells","summed energy of Cells in EMCAL vs PHOS", nptbins*2,ptmin,ptmax*2,nptbins*2,ptmin,ptmax*2);
847 fhCaloCorrECells->SetXTitle("#Sigma E of Cells in EMCAL (GeV)");
848 fhCaloCorrECells->SetYTitle("#Sigma E of Cells in PHOS (GeV)");
849 outputContainer->Add(fhCaloCorrECells);
851 //Calorimeter VS V0 signal
852 fhCaloV0SCorrNClusters = new TH2F ("hCaloV0SNClusters",Form("# clusters in %s vs V0 signal",fCalorimeter.Data()), nv0sbins,nv0smin,nv0smax,nclbins,nclmin,nclmax);
853 fhCaloV0SCorrNClusters->SetXTitle("V0 signal");
854 fhCaloV0SCorrNClusters->SetYTitle(Form("number of clusters in %s",fCalorimeter.Data()));
855 outputContainer->Add(fhCaloV0SCorrNClusters);
857 fhCaloV0SCorrEClusters = new TH2F ("hCaloV0SEClusters",Form("summed energy of clusters in %s vs V0 signal",fCalorimeter.Data()), nv0sbins,nv0smin,nv0smax,nptbins,ptmin,ptmax);
858 fhCaloV0SCorrEClusters->SetXTitle("V0 signal");
859 fhCaloV0SCorrEClusters->SetYTitle(Form("#Sigma E of clusters in %s (GeV)",fCalorimeter.Data()));
860 outputContainer->Add(fhCaloV0SCorrEClusters);
862 fhCaloV0SCorrNCells = new TH2F ("hCaloV0SNCells",Form("# Cells in %s vs V0 signal",fCalorimeter.Data()), nv0sbins,nv0smin,nv0smax, ncebins,ncemin,ncemax);
863 fhCaloV0SCorrNCells->SetXTitle("V0 signal");
864 fhCaloV0SCorrNCells->SetYTitle(Form("number of Cells in %s",fCalorimeter.Data()));
865 outputContainer->Add(fhCaloV0SCorrNCells);
867 fhCaloV0SCorrECells = new TH2F ("hCaloV0SECells",Form("summed energy of Cells in %s vs V0 signal",fCalorimeter.Data()), nv0sbins,nv0smin,nv0smax,nptbins,ptmin,ptmax);
868 fhCaloV0SCorrECells->SetXTitle("V0 signal");
869 fhCaloV0SCorrECells->SetYTitle(Form("#Sigma E of Cells in %s (GeV)",fCalorimeter.Data()));
870 outputContainer->Add(fhCaloV0SCorrECells);
872 //Calorimeter VS V0 multiplicity
873 fhCaloV0MCorrNClusters = new TH2F ("hCaloV0MNClusters",Form("# clusters in %s vs V0 signal",fCalorimeter.Data()), nv0mbins,nv0mmin,nv0mmax,nclbins,nclmin,nclmax);
874 fhCaloV0MCorrNClusters->SetXTitle("V0 signal");
875 fhCaloV0MCorrNClusters->SetYTitle(Form("number of clusters in %s",fCalorimeter.Data()));
876 outputContainer->Add(fhCaloV0MCorrNClusters);
878 fhCaloV0MCorrEClusters = new TH2F ("hCaloV0MEClusters",Form("summed energy of clusters in %s vs V0 signal",fCalorimeter.Data()), nv0mbins,nv0mmin,nv0mmax,nptbins,ptmin,ptmax);
879 fhCaloV0MCorrEClusters->SetXTitle("V0 signal");
880 fhCaloV0MCorrEClusters->SetYTitle(Form("#Sigma E of clusters in %s (GeV)",fCalorimeter.Data()));
881 outputContainer->Add(fhCaloV0MCorrEClusters);
883 fhCaloV0MCorrNCells = new TH2F ("hCaloV0MNCells",Form("# Cells in %s vs V0 signal",fCalorimeter.Data()), nv0mbins,nv0mmin,nv0mmax, ncebins,ncemin,ncemax);
884 fhCaloV0MCorrNCells->SetXTitle("V0 signal");
885 fhCaloV0MCorrNCells->SetYTitle(Form("number of Cells in %s",fCalorimeter.Data()));
886 outputContainer->Add(fhCaloV0MCorrNCells);
888 fhCaloV0MCorrECells = new TH2F ("hCaloV0MECells",Form("summed energy of Cells in %s vs V0 signal",fCalorimeter.Data()), nv0mbins,nv0mmin,nv0mmax,nptbins,ptmin,ptmax);
889 fhCaloV0MCorrECells->SetXTitle("V0 signal");
890 fhCaloV0MCorrECells->SetYTitle(Form("#Sigma E of Cells in %s (GeV)",fCalorimeter.Data()));
891 outputContainer->Add(fhCaloV0MCorrECells);
893 //Calorimeter VS Track multiplicity
894 fhCaloTrackMCorrNClusters = new TH2F ("hCaloTrackMNClusters",Form("# clusters in %s vs # tracks",fCalorimeter.Data()), ntrmbins,ntrmmin,ntrmmax,nclbins,nclmin,nclmax);
895 fhCaloTrackMCorrNClusters->SetXTitle("# tracks");
896 fhCaloTrackMCorrNClusters->SetYTitle(Form("number of clusters in %s",fCalorimeter.Data()));
897 outputContainer->Add(fhCaloTrackMCorrNClusters);
899 fhCaloTrackMCorrEClusters = new TH2F ("hCaloTrackMEClusters",Form("summed energy of clusters in %s vs # tracks",fCalorimeter.Data()), ntrmbins,ntrmmin,ntrmmax,nptbins,ptmin,ptmax);
900 fhCaloTrackMCorrEClusters->SetXTitle("# tracks");
901 fhCaloTrackMCorrEClusters->SetYTitle(Form("#Sigma E of clusters in %s (GeV)",fCalorimeter.Data()));
902 outputContainer->Add(fhCaloTrackMCorrEClusters);
904 fhCaloTrackMCorrNCells = new TH2F ("hCaloTrackMNCells",Form("# Cells in %s vs # tracks",fCalorimeter.Data()), ntrmbins,ntrmmin,ntrmmax, ncebins,ncemin,ncemax);
905 fhCaloTrackMCorrNCells->SetXTitle("# tracks");
906 fhCaloTrackMCorrNCells->SetYTitle(Form("number of Cells in %s",fCalorimeter.Data()));
907 outputContainer->Add(fhCaloTrackMCorrNCells);
909 fhCaloTrackMCorrECells = new TH2F ("hCaloTrackMECells",Form("summed energy of Cells in %s vs # tracks",fCalorimeter.Data()), ntrmbins,ntrmmin,ntrmmax,nptbins,ptmin,ptmax);
910 fhCaloTrackMCorrECells->SetXTitle("# tracks");
911 fhCaloTrackMCorrECells->SetYTitle(Form("#Sigma E of Cells in %s (GeV)",fCalorimeter.Data()));
912 outputContainer->Add(fhCaloTrackMCorrECells);
915 }//correlate calorimeters
919 fhEMod = new TH2F ("hE_Mod","Cluster reconstructed Energy in each present Module",nptbins,ptmin,ptmax,fNModules,0,fNModules);
920 fhEMod->SetXTitle("E (GeV)");
921 fhEMod->SetYTitle("Module");
922 outputContainer->Add(fhEMod);
924 fhNClustersMod = new TH2F ("hNClusters_Mod","# clusters vs Module", nclbins,nclmin,nclmax,fNModules,0,fNModules);
925 fhNClustersMod->SetXTitle("number of clusters");
926 fhNClustersMod->SetYTitle("Module");
927 outputContainer->Add(fhNClustersMod);
929 fhNCellsMod = new TH2F ("hNCells_Mod","# cells vs Module", fNMaxCols*fNMaxRows,0,fNMaxCols*fNMaxRows,fNModules,0,fNModules);
930 fhNCellsMod->SetXTitle("n cells");
931 fhNCellsMod->SetYTitle("Module");
932 outputContainer->Add(fhNCellsMod);
934 Int_t colmaxs = fNMaxCols;
935 Int_t rowmaxs = fNMaxRows;
936 if(fCalorimeter=="EMCAL"){
938 rowmaxs=Int_t(fNModules/2)*fNMaxRows;
941 rowmaxs=fNModules*fNMaxRows;
944 fhGridCellsMod = new TH2F ("hGridCells",Form("Entries in grid of cells"),
945 colmaxs+2,-1.5,colmaxs+0.5, rowmaxs+2,-1.5,rowmaxs+0.5);
946 fhGridCellsMod->SetYTitle("row (phi direction)");
947 fhGridCellsMod->SetXTitle("column (eta direction)");
948 outputContainer->Add(fhGridCellsMod);
950 fhGridCellsEMod = new TH2F ("hGridCellsE","Accumulated energy in grid of cells",
951 colmaxs+2,-1.5,colmaxs+0.5, rowmaxs+2,-1.5,rowmaxs+0.5);
952 fhGridCellsEMod->SetYTitle("row (phi direction)");
953 fhGridCellsEMod->SetXTitle("column (eta direction)");
954 outputContainer->Add(fhGridCellsEMod);
956 fhGridCellsTimeMod = new TH2F ("hGridCellsTime","Accumulated time in grid of cells",
957 colmaxs+2,-1.5,colmaxs+0.5, rowmaxs+2,-1.5,rowmaxs+0.5);
958 fhGridCellsTimeMod->SetYTitle("row (phi direction)");
959 fhGridCellsTimeMod->SetXTitle("column (eta direction)");
960 outputContainer->Add(fhGridCellsTimeMod);
962 fhNCellsPerClusterMod = new TH2F*[fNModules];
963 fhNCellsPerClusterModNoCut = new TH2F*[fNModules];
964 fhTimeAmpPerRCU = new TH2F*[fNModules*fNRCU];
965 fhIMMod = new TH2F*[fNModules];
967 for(Int_t imod = 0; imod < fNModules; imod++){
969 fhNCellsPerClusterMod[imod] = new TH2F (Form("hNCellsPerCluster_Mod%d",imod),
970 Form("# cells per cluster vs cluster energy in Module %d",imod),
971 nptbins,ptmin,ptmax, nceclbins,nceclmin,nceclmax);
972 fhNCellsPerClusterMod[imod]->SetXTitle("E (GeV)");
973 fhNCellsPerClusterMod[imod]->SetYTitle("n cells");
974 outputContainer->Add(fhNCellsPerClusterMod[imod]);
976 fhNCellsPerClusterModNoCut[imod] = new TH2F (Form("hNCellsPerClusterNoCut_Mod%d",imod),
977 Form("# cells per cluster vs cluster energy in Module %d, no cut",imod),
978 nptbins,ptmin,ptmax, nceclbins,nceclmin,nceclmax);
979 fhNCellsPerClusterModNoCut[imod]->SetXTitle("E (GeV)");
980 fhNCellsPerClusterModNoCut[imod]->SetYTitle("n cells");
981 outputContainer->Add(fhNCellsPerClusterModNoCut[imod]);
983 if(GetReader()->GetDataType()==AliCaloTrackReader::kESD) {
985 for(Int_t ircu = 0; ircu < fNRCU; ircu++){
986 fhTimeAmpPerRCU[imod*fNRCU+ircu] = new TH2F (Form("hTimeAmp_Mod%d_RCU%d",imod,ircu),
987 Form("Cell Energy vs Cell Time in Module %d, RCU %d ",imod,ircu),
988 nptbins,ptmin,ptmax,ntimebins,timemin,timemax);
989 fhTimeAmpPerRCU[imod*fNRCU+ircu]->SetXTitle("E (GeV)");
990 fhTimeAmpPerRCU[imod*fNRCU+ircu]->SetYTitle("time (ns)");
991 outputContainer->Add(fhTimeAmpPerRCU[imod*fNRCU+ircu]);
996 if(fFillAllPi0Histo){
997 fhIMMod[imod] = new TH2F (Form("hIM_Mod%d",imod),
998 Form("Cluster pairs Invariant mass vs reconstructed pair energy in Module %d, n cell > 1",imod),
999 nptbins,ptmin,ptmax,nmassbins,massmin,massmax);
1000 fhIMMod[imod]->SetXTitle("p_{T, cluster pairs} (GeV) ");
1001 fhIMMod[imod]->SetYTitle("M_{cluster pairs} (GeV/c^{2})");
1002 outputContainer->Add(fhIMMod[imod]);
1007 //Monte Carlo Histograms
1009 TString particleName[] = { "Photon", "Pi0", "Eta", "Electron", "NeutralHadron", "ChargedHadron" };
1012 for(Int_t iPart = 0; iPart < 6; iPart++){
1013 for(Int_t iCh = 0; iCh < 2; iCh++){
1015 fhRecoMCDeltaE[iPart][iCh] = new TH2F (Form("hRecoMCDeltaE_%s_Match%d",particleName[iPart].Data(),iCh),
1016 Form("MC - Reco E, %s, Matched %d",particleName[iPart].Data(),iCh),
1017 nptbins, ptmin, ptmax, nptbins*2,-ptmax,ptmax);
1018 fhRecoMCDeltaE[iPart][iCh]->SetXTitle("#Delta E (GeV)");
1019 outputContainer->Add(fhRecoMCDeltaE[iPart][iCh]);
1021 fhRecoMCDeltaPhi[iPart][iCh] = new TH2F (Form("hRecoMCDeltaPhi_%s_Match%d",particleName[iPart].Data(),iCh),
1022 Form("MC - Reco #phi, %s, Matched %d",particleName[iPart].Data(),iCh),
1023 nptbins, ptmin, ptmax, nphibins*2,-phimax,phimax);
1024 fhRecoMCDeltaPhi[iPart][iCh]->SetXTitle("#Delta #phi (rad)");
1025 outputContainer->Add(fhRecoMCDeltaPhi[iPart][iCh]);
1027 fhRecoMCDeltaEta[iPart][iCh] = new TH2F (Form("hRecoMCDeltaEta_%s_Match%d",particleName[iPart].Data(),iCh),
1028 Form("MC- Reco #eta, %s, Matched %d",particleName[iPart].Data(),iCh),
1029 nptbins, ptmin, ptmax,netabins*2,-etamax,etamax);
1030 fhRecoMCDeltaEta[iPart][iCh]->SetXTitle("#Delta #eta ");
1031 outputContainer->Add(fhRecoMCDeltaEta[iPart][iCh]);
1033 fhRecoMCE[iPart][iCh] = new TH2F (Form("hRecoMCE_%s_Match%d",particleName[iPart].Data(),iCh),
1034 Form("E distribution, reconstructed vs generated, %s, Matched %d",particleName[iPart].Data(),iCh),
1035 nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
1036 fhRecoMCE[iPart][iCh]->SetXTitle("E_{rec} (GeV)");
1037 fhRecoMCE[iPart][iCh]->SetYTitle("E_{gen} (GeV)");
1038 outputContainer->Add(fhRecoMCE[iPart][iCh]);
1040 fhRecoMCPhi[iPart][iCh] = new TH2F (Form("hRecoMCPhi_%s_Match%d",particleName[iPart].Data(),iCh),
1041 Form("#phi distribution, reconstructed vs generated, %s, Matched %d",particleName[iPart].Data(),iCh),
1042 nphibins,phimin,phimax, nphibins,phimin,phimax);
1043 fhRecoMCPhi[iPart][iCh]->SetXTitle("#phi_{rec} (rad)");
1044 fhRecoMCPhi[iPart][iCh]->SetYTitle("#phi_{gen} (rad)");
1045 outputContainer->Add(fhRecoMCPhi[iPart][iCh]);
1047 fhRecoMCEta[iPart][iCh] = new TH2F (Form("hRecoMCEta_%s_Match%d",particleName[iPart].Data(),iCh),
1048 Form("#eta distribution, reconstructed vs generated, %s, Matched %d",particleName[iPart].Data(),iCh),
1049 netabins,etamin,etamax,netabins,etamin,etamax);
1050 fhRecoMCEta[iPart][iCh]->SetXTitle("#eta_{rec} ");
1051 fhRecoMCEta[iPart][iCh]->SetYTitle("#eta_{gen} ");
1052 outputContainer->Add(fhRecoMCEta[iPart][iCh]);
1057 for(Int_t iPart = 0; iPart < 4; iPart++){
1058 fhGenMCE[iPart] = new TH1F(Form("hGenMCE_%s",particleName[iPart].Data()) ,
1059 Form("p_{T} of generated %s",particleName[iPart].Data()),
1060 nptbins,ptmin,ptmax);
1061 fhGenMCEtaPhi[iPart] = new TH2F(Form("hGenMCEtaPhi_%s",particleName[iPart].Data()),
1062 Form("Y vs #phi of generated %s",particleName[iPart].Data()),
1063 netabins,etamin,etamax,nphibins,phimin,phimax);
1065 fhGenMCE[iPart] ->SetXTitle("p_{T} (GeV/c)");
1066 fhGenMCEtaPhi[iPart]->SetXTitle("#eta");
1067 fhGenMCEtaPhi[iPart]->SetYTitle("#phi (rad)");
1069 outputContainer->Add(fhGenMCE[iPart]);
1070 outputContainer->Add(fhGenMCEtaPhi[iPart]);
1073 fhGenMCAccE[iPart] = new TH1F(Form("hGenMCAccE_%s",particleName[iPart].Data()) ,
1074 Form("p_{T} of generated %s",particleName[iPart].Data()),
1075 nptbins,ptmin,ptmax);
1076 fhGenMCAccEtaPhi[iPart] = new TH2F(Form("hGenMCAccEtaPhi_%s",particleName[iPart].Data()),
1077 Form("Y vs #phi of generated %s",particleName[iPart].Data()),
1078 netabins,etamin,etamax,nphibins,phimin,phimax);
1080 fhGenMCAccE[iPart] ->SetXTitle("p_{T} (GeV/c)");
1081 fhGenMCAccEtaPhi[iPart]->SetXTitle("#eta");
1082 fhGenMCAccEtaPhi[iPart]->SetYTitle("#phi (rad)");
1084 outputContainer->Add(fhGenMCAccE[iPart]);
1085 outputContainer->Add(fhGenMCAccEtaPhi[iPart]);
1089 //Vertex of generated particles
1091 fhEMVxyz = new TH2F ("hEMVxyz","Production vertex of reconstructed ElectroMagnetic particles",nvdistbins,vdistmin,vdistmax,nvdistbins,vdistmin,vdistmax);//,100,0,500);
1092 fhEMVxyz->SetXTitle("v_{x}");
1093 fhEMVxyz->SetYTitle("v_{y}");
1094 //fhEMVxyz->SetZTitle("v_{z}");
1095 outputContainer->Add(fhEMVxyz);
1097 fhHaVxyz = new TH2F ("hHaVxyz","Production vertex of reconstructed hadrons",nvdistbins,vdistmin,vdistmax,nvdistbins,vdistmin,vdistmax);//,100,0,500);
1098 fhHaVxyz->SetXTitle("v_{x}");
1099 fhHaVxyz->SetYTitle("v_{y}");
1100 //fhHaVxyz->SetZTitle("v_{z}");
1101 outputContainer->Add(fhHaVxyz);
1103 fhEMR = new TH2F ("hEMR","Distance to production vertex of reconstructed ElectroMagnetic particles vs E rec",nptbins,ptmin,ptmax,nvdistbins,vdistmin,vdistmax);
1104 fhEMR->SetXTitle("E (GeV)");
1105 fhEMR->SetYTitle("TMath::Sqrt(v_{x}^{2}+v_{y}^{2})");
1106 outputContainer->Add(fhEMR);
1108 fhHaR = new TH2F ("hHaR","Distance to production vertex of reconstructed Hadrons vs E rec",nptbins,ptmin,ptmax,nvdistbins,vdistmin,vdistmax);
1109 fhHaR->SetXTitle("E (GeV)");
1110 fhHaR->SetYTitle("TMath::Sqrt(v_{x}^{2}+v_{y}^{2})");
1111 outputContainer->Add(fhHaR);
1116 fhMCEle1pOverE = new TH2F("hMCEle1pOverE","TRACK matches p/E, MC electrons",nptbins,ptmin,ptmax, nPoverEbins,pOverEmin,pOverEmax);
1117 fhMCEle1pOverE->SetYTitle("p/E");
1118 fhMCEle1pOverE->SetXTitle("p_{T} (GeV/c)");
1119 outputContainer->Add(fhMCEle1pOverE);
1121 fhMCEle1dR = new TH1F("hMCEle1dR","TRACK matches dR, MC electrons",ndRbins,dRmin,dRmax);
1122 fhMCEle1dR->SetXTitle("#Delta R (rad)");
1123 outputContainer->Add(fhMCEle1dR) ;
1125 fhMCEle2MatchdEdx = new TH2F("hMCEle2MatchdEdx","dE/dx vs. p for all matches, MC electrons",nptbins,ptmin,ptmax,ndedxbins,dedxmin,dedxmax);
1126 fhMCEle2MatchdEdx->SetXTitle("p (GeV/c)");
1127 fhMCEle2MatchdEdx->SetYTitle("<dE/dx>");
1128 outputContainer->Add(fhMCEle2MatchdEdx);
1130 fhMCChHad1pOverE = new TH2F("hMCChHad1pOverE","TRACK matches p/E, MC charged hadrons",nptbins,ptmin,ptmax, nPoverEbins,pOverEmin,pOverEmax);
1131 fhMCChHad1pOverE->SetYTitle("p/E");
1132 fhMCChHad1pOverE->SetXTitle("p_{T} (GeV/c)");
1133 outputContainer->Add(fhMCChHad1pOverE);
1135 fhMCChHad1dR = new TH1F("hMCChHad1dR","TRACK matches dR, MC charged hadrons",ndRbins,dRmin,dRmax);
1136 fhMCChHad1dR->SetXTitle("#Delta R (rad)");
1137 outputContainer->Add(fhMCChHad1dR) ;
1139 fhMCChHad2MatchdEdx = new TH2F("hMCChHad2MatchdEdx","dE/dx vs. p for all matches, MC charged hadrons",nptbins,ptmin,ptmax,ndedxbins,dedxmin,dedxmax);
1140 fhMCChHad2MatchdEdx->SetXTitle("p (GeV/c)");
1141 fhMCChHad2MatchdEdx->SetYTitle("<dE/dx>");
1142 outputContainer->Add(fhMCChHad2MatchdEdx);
1144 fhMCNeutral1pOverE = new TH2F("hMCNeutral1pOverE","TRACK matches p/E, MC neutrals",nptbins,ptmin,ptmax, nPoverEbins,pOverEmin,pOverEmax);
1145 fhMCNeutral1pOverE->SetYTitle("p/E");
1146 fhMCNeutral1pOverE->SetXTitle("p_{T} (GeV/c)");
1147 outputContainer->Add(fhMCNeutral1pOverE);
1149 fhMCNeutral1dR = new TH1F("hMCNeutral1dR","TRACK matches dR, MC neutrals",ndRbins,dRmin,dRmax);
1150 fhMCNeutral1dR->SetXTitle("#Delta R (rad)");
1151 outputContainer->Add(fhMCNeutral1dR) ;
1153 fhMCNeutral2MatchdEdx = new TH2F("hMCNeutral2MatchdEdx","dE/dx vs. p for all matches, MC neutrals",nptbins,ptmin,ptmax,ndedxbins,dedxmin,dedxmax);
1154 fhMCNeutral2MatchdEdx->SetXTitle("p (GeV/c)");
1155 fhMCNeutral2MatchdEdx->SetYTitle("<dE/dx>");
1156 outputContainer->Add(fhMCNeutral2MatchdEdx);
1158 fhMCEle1pOverER02 = new TH2F("hMCEle1pOverER02","TRACK matches p/E, MC electrons",nptbins,ptmin,ptmax, nPoverEbins,pOverEmin,pOverEmax);
1159 fhMCEle1pOverER02->SetYTitle("p/E");
1160 fhMCEle1pOverER02->SetXTitle("p_{T} (GeV/c)");
1161 outputContainer->Add(fhMCEle1pOverER02);
1163 fhMCChHad1pOverER02 = new TH2F("hMCChHad1pOverER02","TRACK matches p/E, MC charged hadrons",nptbins,ptmin,ptmax, nPoverEbins,pOverEmin,pOverEmax);
1164 fhMCChHad1pOverER02->SetYTitle("p/E");
1165 fhMCChHad1pOverER02->SetXTitle("p_{T} (GeV/c)");
1166 outputContainer->Add(fhMCChHad1pOverER02);
1168 fhMCNeutral1pOverER02 = new TH2F("hMCNeutral1pOverER02","TRACK matches p/E, MC neutrals",nptbins,ptmin,ptmax, nPoverEbins,pOverEmin,pOverEmax);
1169 fhMCNeutral1pOverER02->SetYTitle("p/E");
1170 fhMCNeutral1pOverER02->SetXTitle("p_{T} (GeV/c)");
1171 outputContainer->Add(fhMCNeutral1pOverER02);
1174 // for(Int_t i = 0; i < outputContainer->GetEntries() ; i++)
1175 // printf("i=%d, name= %s\n",i,outputContainer->At(i)->GetName());
1177 return outputContainer;
1180 //__________________________________________________
1181 void AliAnaCalorimeterQA::Init()
1183 //Check if the data or settings are ok
1185 if(fCalorimeter != "PHOS" && fCalorimeter !="EMCAL")
1186 AliFatal(Form("Wrong calorimeter name <%s>", fCalorimeter.Data()));
1188 if(GetReader()->GetDataType()== AliCaloTrackReader::kMC)
1189 AliFatal("Analysis of reconstructed data, MC reader not aplicable");
1194 //__________________________________________________
1195 void AliAnaCalorimeterQA::InitParameters()
1197 //Initialize the parameters of the analysis.
1198 AddToHistogramsName("AnaCaloQA_");
1200 fCalorimeter = "EMCAL"; //or PHOS
1201 fNModules = 12; // set maximum to maximum number of EMCAL modules
1202 fNRCU = 2; // set maximum number of RCU in EMCAL per SM
1204 fTimeCutMax = 9999999;
1205 fEMCALCellAmpMin = 0.2;
1206 fPHOSCellAmpMin = 0.2;
1210 //__________________________________________________________________
1211 void AliAnaCalorimeterQA::Print(const Option_t * opt) const
1213 //Print some relevant parameters set for the analysis
1217 printf("**** Print %s %s ****\n", GetName(), GetTitle() ) ;
1218 AliAnaPartCorrBaseClass::Print(" ");
1220 printf("Select Calorimeter %s \n",fCalorimeter.Data());
1221 printf("Time Cut: %3.1f < TOF < %3.1f\n", fTimeCutMin, fTimeCutMax);
1222 printf("EMCAL Min Amplitude : %2.1f GeV/c\n", fEMCALCellAmpMin) ;
1223 printf("PHOS Min Amplitude : %2.1f GeV/c\n", fPHOSCellAmpMin) ;
1227 //__________________________________________________________________
1228 void AliAnaCalorimeterQA::MakeAnalysisFillHistograms()
1230 //Fill Calorimeter QA histograms
1233 TLorentzVector mom ;
1234 TLorentzVector mom2 ;
1235 TObjArray * caloClusters = NULL;
1237 Int_t *labels = 0x0;
1238 Int_t nCaloClusters = 0 ;
1239 Int_t nCaloClustersAccepted = 0 ;
1240 Int_t nCaloCellsPerCluster = 0 ;
1241 Bool_t matched = kFALSE;
1243 Int_t bc = GetReader()->GetInputEvent()->GetBunchCrossNumber();
1245 //Get vertex for photon momentum calculation and event selection
1246 Double_t v[3] = {0,0,0}; //vertex ;
1247 GetReader()->GetVertex(v);
1248 if (TMath::Abs(v[2]) > GetZvertexCut()) return ;
1250 //Play with the MC stack if available
1251 //Get the MC arrays and do some checks
1253 if(GetReader()->ReadStack()){
1256 AliFatal("Stack not available, is the MC handler called?\n");
1258 //Fill some pure MC histograms, only primaries.
1259 for(Int_t i=0 ; i<GetMCStack()->GetNprimary(); i++){//Only primary particles, for all MC transport put GetNtrack()
1260 TParticle *primary = GetMCStack()->Particle(i) ;
1262 if (primary->GetStatusCode() > 11) continue; //Working for PYTHIA and simple generators, check for HERWIG
1263 primary->Momentum(mom);
1264 MCHistograms(mom,TMath::Abs(primary->GetPdgCode()));
1267 else if(GetReader()->ReadAODMCParticles()){
1269 if(!GetReader()->GetAODMCParticles(0))
1270 AliFatal("AODMCParticles not available!");
1272 //Fill some pure MC histograms, only primaries.
1273 for(Int_t i=0 ; i < (GetReader()->GetAODMCParticles(0))->GetEntriesFast(); i++){
1274 AliAODMCParticle *aodprimary = (AliAODMCParticle*) (GetReader()->GetAODMCParticles(0))->At(i) ;
1276 if (!aodprimary->IsPrimary()) continue; //accept all which is not MC transport generated. Don't know how to avoid partons
1278 mom.SetPxPyPzE(aodprimary->Px(),aodprimary->Py(),aodprimary->Pz(),aodprimary->E());
1279 MCHistograms(mom,TMath::Abs(aodprimary->GetPdgCode()));
1285 //Get List with CaloClusters
1286 if (fCalorimeter == "PHOS") caloClusters = GetPHOSClusters();
1287 else if (fCalorimeter == "EMCAL") caloClusters = GetEMCALClusters();
1289 AliFatal(Form("AliAnaCalorimeterQA::MakeAnalysisFillHistograms() - Wrong calorimeter name <%s>, END\n", fCalorimeter.Data()));
1291 //Get List with CaloCells
1292 AliVCaloCells * cell = 0x0;
1293 if(fCalorimeter == "PHOS") cell = GetPHOSCells();
1294 else cell = GetEMCALCells();
1296 if(!caloClusters || !cell) {
1297 AliFatal(Form("AliAnaCalorimeterQA::MakeAnalysisFillHistograms() - No CaloClusters or CaloCells available\n"));
1298 return; // trick coverity
1301 //----------------------------------------------------------
1302 //Correlate Calorimeters and V0 and track Multiplicity
1303 //----------------------------------------------------------
1305 if(fCorrelate) Correlate();
1307 //----------------------------------------------------------
1309 //----------------------------------------------------------
1311 nCaloClusters = caloClusters->GetEntriesFast() ;
1312 Int_t *nClustersInModule = new Int_t[fNModules];
1313 for(Int_t imod = 0; imod < fNModules; imod++ ) nClustersInModule[imod] = 0;
1316 printf("AliAnaCalorimeterQA::MakeAnalysisFillHistograms() - In %s there are %d clusters \n", fCalorimeter.Data(), nCaloClusters);
1318 AliVTrack * track = 0x0;
1322 //Loop over CaloClusters
1323 //if(nCaloClusters > 0)printf("QA : Vertex Cut passed %f, cut %f, entries %d, %s\n",v[2], 40., nCaloClusters, fCalorimeter.Data());
1324 for(Int_t iclus = 0; iclus < nCaloClusters; iclus++){
1326 if(GetDebug() > 0) printf("AliAnaCalorimeterQA::MakeAnalysisFillHistograms() - cluster: %d/%d, data %d \n",
1327 iclus+1,nCaloClusters,GetReader()->GetDataType());
1329 AliVCluster* clus = (AliVCluster*)caloClusters->At(iclus);
1331 //Get cluster kinematics
1332 clus->GetPosition(pos);
1333 clus->GetMomentum(mom,v);
1335 //Check only certain regions
1337 if(IsFiducialCutOn()) in = GetFiducialCut()->IsInFiducialCut(mom,fCalorimeter) ;
1341 nLabel = clus->GetNLabels();
1342 labels = clus->GetLabels();
1345 nCaloCellsPerCluster = clus->GetNCells();
1346 //if(mom.E() > 10 && nCaloCellsPerCluster == 1 ) printf("%s:************** E = %f ********** ncells = %d\n",fCalorimeter.Data(), mom.E(),nCaloCellsPerCluster);
1348 // Cluster mathed with track?
1349 matched = GetCaloPID()->IsTrackMatched(clus,GetCaloUtils());
1352 //======================
1354 //======================
1356 //Get list of contributors
1357 UShort_t * indexList = clus->GetCellsAbsId() ;
1359 if(fFillAllPosHisto) FillCellPositionHistograms(nCaloCellsPerCluster,indexList,pos,mom.E());
1361 // Get the fraction of the cluster energy that carries the cell with highest energy
1362 Float_t maxCellFraction = 0.;
1363 Int_t absIdMax = GetCaloUtils()->GetMaxEnergyCell(cell, clus,maxCellFraction);
1364 Int_t smMax =0; Int_t ietaaMax=-1; Int_t iphiiMax = 0; Int_t rcuMax = 0;
1365 smMax = GetModuleNumberCellIndexes(absIdMax,fCalorimeter, ietaaMax, iphiiMax, rcuMax);
1369 if (clus->E() < 2.){
1370 fhLambda0vsClusterMaxCellDiffE0->Fill(clus->GetM02(), maxCellFraction);
1371 fhNCellsvsClusterMaxCellDiffE0 ->Fill(nCaloCellsPerCluster,maxCellFraction);
1373 else if(clus->E() < 6.){
1374 fhLambda0vsClusterMaxCellDiffE2->Fill(clus->GetM02(), maxCellFraction);
1375 fhNCellsvsClusterMaxCellDiffE2 ->Fill(nCaloCellsPerCluster,maxCellFraction);
1378 fhLambda0vsClusterMaxCellDiffE6->Fill(clus->GetM02(), maxCellFraction);
1379 fhNCellsvsClusterMaxCellDiffE6 ->Fill(nCaloCellsPerCluster,maxCellFraction);
1382 fhNCellsPerClusterNoCut ->Fill(clus->E(), nCaloCellsPerCluster);
1383 nModule = GetModuleNumber(clus);
1384 if(nModule >=0 && nModule < fNModules) fhNCellsPerClusterModNoCut[nModule]->Fill(clus->E(), nCaloCellsPerCluster);
1386 fhClusterMaxCellDiffNoCut->Fill(clus->E(),maxCellFraction);
1388 // Look inside the cluster
1389 // Calculate average time of cells in cluster and weighted average
1390 Float_t rawEnergy = 0;
1393 Double_t tmax = cell->GetCellTime(absIdMax);
1394 Double_t averTime = 0;
1395 Double_t weightedTime = 0;
1396 Double_t weight = 0;
1397 Double_t averTimeNoMax = 0;
1398 Double_t weightedTimeNoMax = 0;
1399 Double_t weightNoMax = 0;
1401 for (Int_t ipos = 0; ipos < nCaloCellsPerCluster; ipos++) {
1405 Int_t id = indexList[ipos];
1406 Int_t nModuleCell = GetModuleNumberCellIndexes(id,fCalorimeter, icol, irow, iRCU);
1408 //Recalibrate energy
1409 if(GetCaloUtils()->IsRecalibrationOn()){
1410 if(fCalorimeter == "PHOS") {
1411 factor = GetCaloUtils()->GetPHOSChannelRecalibrationFactor(nModuleCell,icol,irow);
1414 factor = GetCaloUtils()->GetEMCALChannelRecalibrationFactor(nModuleCell,icol,irow);
1420 time = cell->GetCellTime(id);
1421 if(fCalorimeter=="EMCAL" && GetCaloUtils()->GetEMCALRecoUtils()->IsTimeRecalibrationOn()){
1422 GetCaloUtils()->GetEMCALRecoUtils()->RecalibrateCellTime(id,bc,time);
1423 if(id==absIdMax) tmax = time;
1426 // printf("QA cluster : Id %d, Time org %e, Time new %e; Amp org %f, Amp new %f\n",
1427 // id, cell->GetCellTime(id),time, cell->GetCellAmplitude(id),cell->GetCellAmplitude(id)*factor);
1429 // Get the energy of the cluster without the non linearity correction
1430 rawEnergy += cell->GetCellAmplitude(id)*factor;
1432 Double_t w = GetCaloUtils()->GetEMCALRecoUtils()->GetCellWeight(cell->GetCellAmplitude(id),rawEnergy);
1433 averTime += time*1e9;
1434 weightedTime += time*1e9 * w;
1437 averTimeNoMax += time*1e9;
1438 weightedTimeNoMax += time*1e9 * w;
1444 averTime /= nCaloCellsPerCluster;
1445 if(nCaloCellsPerCluster > 1 ) averTimeNoMax /= (nCaloCellsPerCluster-1);
1447 if(weight > 0 ) tof = weightedTime / weight;
1450 tof = clus->GetTOF()*1e9;
1452 printf("AliAnaCalorimeterQA:: Null weight! Investigate: E %f GeV, ncells %d, time max cell %f ns, average time %f ns, absIdMax %d, SM %d\n",
1453 rawEnergy,nCaloCellsPerCluster, tmax, averTime,absIdMax,GetModuleNumber(clus));
1456 //printf("QA: E %f, tof org %g, tof new %g, ncells %d\n",clus->E(),clus->GetTOF()*1.e9,tof, clus->GetNCells());
1458 if(weightNoMax > 0 ) weightedTimeNoMax /= weightNoMax;
1460 //Cut on time of clusters
1461 if(tof < fTimeCutMin || tof > fTimeCutMax){
1462 printf("Remove cluster with TOF %f\n",tof);
1466 //======================
1467 //Bad clusters selection
1468 //======================
1469 //Check bad clusters if rejection was not on
1471 Bool_t badCluster = kFALSE;
1472 if(fCalorimeter=="EMCAL" && !GetCaloUtils()->GetEMCALRecoUtils()->IsRejectExoticCluster()){
1473 //Bad clusters histograms
1474 Float_t minNCells = 1;
1475 if(clus->E() > 7) minNCells = TMath::Max(1,TMath::FloorNint(1 + TMath::Log(clus->E() - 7 )*1.7 ));
1476 if(nCaloCellsPerCluster <= minNCells) {
1477 //if(clus->GetM02() < 0.05) {
1481 fhBadClusterEnergy ->Fill(clus->E());
1482 fhBadClusterTimeEnergy ->Fill(clus->E(),tof);
1483 fhBadClusterMaxCellDiff->Fill(clus->E(),maxCellFraction);
1485 if(clus->GetM02() > 0 || TMath::Abs(clus->GetM20()) > 0 || clus->GetDispersion() > 0){
1487 fhBadClusterL0 ->Fill(clus->E(),clus->GetM02());
1488 fhBadClusterL1 ->Fill(clus->E(),clus->GetM20());
1489 fhBadClusterD ->Fill(clus->E(),clus->GetDispersion());
1493 //Clusters in event time difference
1495 for(Int_t iclus2 = 0; iclus2 < nCaloClusters; iclus2++ ){
1497 AliVCluster* clus2 = (AliVCluster*)caloClusters->At(iclus2);
1499 if(clus->GetID()==clus2->GetID()) continue;
1501 if(clus->GetM02() > 0.01 && clus2->GetM02() > 0.01) {
1502 Double_t tof2 = clus2->GetTOF();
1503 // approximation in case of time recalibration
1504 if(fCalorimeter=="EMCAL" && GetCaloUtils()->GetEMCALRecoUtils()->IsTimeRecalibrationOn()){
1505 Float_t maxCellFraction2 = -1;
1506 Int_t absIdMax2 = GetCaloUtils()->GetMaxEnergyCell(cell, clus2,maxCellFraction2);
1507 tof2 = cell->GetCellTime(absIdMax2);
1508 GetCaloUtils()->GetEMCALRecoUtils()->RecalibrateCellTime(absIdMax2,bc,tof2);
1510 fhBadClusterPairDiffTimeE ->Fill(clus->E(), tof-tof2*1.e9);
1514 // Max cell compared to other cells in cluster
1515 if(GetReader()->GetDataType()==AliCaloTrackReader::kESD) {
1516 fhBadClusterMaxCellDiffAverageTime ->Fill(clus->E(),tmax-averTime);
1517 fhBadClusterMaxCellDiffWeightTime ->Fill(clus->E(),tmax-weightedTime);
1518 fhBadClusterDiffWeightAverTime ->Fill(clus->E(),weightedTime-averTime);
1519 fhBadClusterMaxCellDiffAverageNoMaxTime->Fill(clus->E(),tmax-averTimeNoMax);
1520 fhBadClusterMaxCellDiffWeightNoMaxTime ->Fill(clus->E(),tmax-weightedTimeNoMax);
1521 //printf("E %f, tmax %f, aver %f, weigh %f, averNoMax %f, weightNoMax %f\n",
1522 // clus->E(),tmax,averTime,weightedTime,averTimeNoMax,weightedTimeNoMax);
1525 if( weight > 0 ) fhBadClusterNoMaxCellWeight ->Fill(clus->E(),weightNoMax/weight);
1526 //if( weight > 0 && weightNoMax > 0.0)printf("weight %f, weightNoMax %f\n",weight,weightNoMax);
1528 for (Int_t ipos = 0; ipos < nCaloCellsPerCluster; ipos++) {
1529 Int_t absId = indexList[ipos];
1530 if(absId!=absIdMax){
1531 Float_t frac = cell->GetCellAmplitude(absId)/cell->GetCellAmplitude(absIdMax);
1532 fhBadClusterMaxCellCloseCellRatio->Fill(mom.E(),frac);
1533 fhBadClusterMaxCellCloseCellDiff ->Fill(mom.E(),cell->GetCellAmplitude(absIdMax)-cell->GetCellAmplitude(absId));
1534 time = cell->GetCellTime(absId);
1535 if(GetReader()->GetDataType()==AliCaloTrackReader::kESD) {
1537 if(fCalorimeter=="EMCAL" && GetCaloUtils()->GetEMCALRecoUtils()->IsTimeRecalibrationOn()){
1538 GetCaloUtils()->GetEMCALRecoUtils()->RecalibrateCellTime(absId,bc,time);
1540 Float_t diff = (tmax-time*1e9);
1541 fhBadCellTimeSpreadRespectToCellMax->Fill(mom.E(), diff);
1551 fhClusterMaxCellDiff->Fill(clus->E(),maxCellFraction);
1552 fhClusterTimeEnergy ->Fill(mom.E(),tof);
1554 //Clusters in event time difference
1555 for(Int_t iclus2 = 0; iclus2 < nCaloClusters; iclus2++ ){
1557 AliVCluster* clus2 = (AliVCluster*) caloClusters->At(iclus2);
1559 if(clus->GetID()==clus2->GetID()) continue;
1561 if(clus->GetM02() > 0.01 && clus2->GetM02() > 0.01) {
1563 Double_t tof2 = clus2->GetTOF();
1564 // approximation in case of time recalibration
1565 if(fCalorimeter=="EMCAL" && GetCaloUtils()->GetEMCALRecoUtils()->IsTimeRecalibrationOn()){
1566 Float_t maxCellFraction2 = -1;
1567 Int_t absIdMax2 = GetCaloUtils()->GetMaxEnergyCell(cell, clus2,maxCellFraction2);
1568 tof2 = cell->GetCellTime(absIdMax2);
1569 GetCaloUtils()->GetEMCALRecoUtils()->RecalibrateCellTime(absIdMax2,bc,tof2);
1572 fhClusterPairDiffTimeE ->Fill(clus->E(), tof-tof2*1.e9);
1576 if(nCaloCellsPerCluster > 1){
1578 // check time of cells respect to max energy cell
1580 if(GetReader()->GetDataType()==AliCaloTrackReader::kESD) {
1581 fhClusterMaxCellDiffAverageTime ->Fill(clus->E(),tmax-averTime);
1582 fhClusterMaxCellDiffWeightTime ->Fill(clus->E(),tmax-weightedTime);
1583 fhClusterDiffWeightAverTime ->Fill(clus->E(), weightedTime-averTime);
1584 fhClusterMaxCellDiffAverageNoMaxTime->Fill(clus->E(),tmax-averTimeNoMax);
1585 fhClusterMaxCellDiffWeightNoMaxTime ->Fill(clus->E(),tmax-weightedTimeNoMax);
1589 if( weight > 0 )fhClusterNoMaxCellWeight ->Fill(clus->E(),weightNoMax/weight);
1591 for (Int_t ipos = 0; ipos < nCaloCellsPerCluster; ipos++) {
1593 Int_t absId = indexList[ipos];
1594 if(absId == absIdMax) continue;
1596 Float_t frac = cell->GetCellAmplitude(absId)/cell->GetCellAmplitude(absIdMax);
1597 fhClusterMaxCellCloseCellRatio->Fill(mom.E(),frac);
1598 fhClusterMaxCellCloseCellDiff ->Fill(mom.E(),cell->GetCellAmplitude(absIdMax)-cell->GetCellAmplitude(absId));
1600 if(GetReader()->GetDataType()==AliCaloTrackReader::kESD) {
1602 time = cell->GetCellTime(absId);
1603 if(fCalorimeter=="EMCAL" && GetCaloUtils()->GetEMCALRecoUtils()->IsTimeRecalibrationOn()){
1604 GetCaloUtils()->GetEMCALRecoUtils()->RecalibrateCellTime(absId,bc,time);
1607 Float_t diff = (tmax-time*1.0e9);
1608 fhCellTimeSpreadRespectToCellMax->Fill(mom.E(), diff);
1609 if(TMath::Abs(TMath::Abs(diff) > 100) && mom.E() > 1 ) fhCellIdCellLargeTimeSpread->Fill(absId);
1612 Int_t sm =0; Int_t ietaa=-1; Int_t iphii = 0; Int_t rcu = 0;
1613 sm = GetModuleNumberCellIndexes(absId,fCalorimeter, ietaa, iphii, rcu);
1614 if(dIphi < TMath::Abs(iphii-iphiiMax)) dIphi = TMath::Abs(iphii-iphiiMax);
1616 if(dIeta < TMath::Abs(ietaa-ietaaMax)) dIeta = TMath::Abs(ietaa-ietaaMax);
1619 Int_t ietaaShift = ietaa;
1620 Int_t ietaaMaxShift = ietaaMax;
1621 if (ietaa > ietaaMax) ietaaMaxShift+=48;
1622 else ietaaShift +=48;
1623 if(dIeta < TMath::Abs(ietaaShift-ietaaMaxShift)) dIeta = TMath::Abs(ietaaShift-ietaaMaxShift);
1626 //if(TMath::Abs(clus->GetM20()) < 0.0001 && clus->GetNCells() > 3){
1627 // 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; Max icol %d, irow %d \n",
1628 // clus->E(), clus->GetNCells(),clus->GetM02(), clus->GetM20(), clus->GetDispersion(),tmax, tof,sm,ietaa,iphii, ietaaMax, iphiiMax);
1632 }// fill cell-cluster histogram loop
1634 if(nCaloCellsPerCluster > 3){
1636 if (mom.E() < 2 ) fhDeltaIEtaDeltaIPhiE0[matched]->Fill(dIeta,dIphi);
1637 else if(mom.E() < 6 ) fhDeltaIEtaDeltaIPhiE2[matched]->Fill(dIeta,dIphi);
1638 else fhDeltaIEtaDeltaIPhiE6[matched]->Fill(dIeta,dIphi);
1640 Float_t dIA = 1.*(dIphi-dIeta)/(dIeta+dIphi);
1641 fhDeltaIA[matched]->Fill(mom.E(),dIA);
1645 fhDeltaIAL0[matched]->Fill(clus->GetM02(),dIA);
1646 fhDeltaIAL1[matched]->Fill(clus->GetM20(),dIA);
1647 fhDeltaIANCells[matched]->Fill(nCaloCellsPerCluster,dIA);
1652 Int_t tag = GetMCAnalysisUtils()->CheckOrigin(labels,nLabel, GetReader(),0);
1653 if( (GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCPhoton) ||
1654 GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCPi0) ||
1655 GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCEta) ) &&
1656 !GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCConversion) ){
1657 fhDeltaIAMC[0]->Fill(mom.E(),dIA);//Pure Photon
1659 else if ( GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCElectron) &&
1660 !GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCConversion) ){
1661 fhDeltaIAMC[1]->Fill(mom.E(),dIA);//Pure electron
1663 else if ( GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCConversion) ){
1664 fhDeltaIAMC[2]->Fill(mom.E(),dIA);//Converted cluster
1667 fhDeltaIAMC[3]->Fill(mom.E(),dIA);//Hadrons
1671 }// 4 cells at least in cluster for size study
1673 }//check time and energy of cells respect to max energy cell if cluster of more than 1 cell
1676 //Get module of cluster
1678 nCaloClustersAccepted++;
1679 if(nModule >=0 && nModule < fNModules) {
1680 if (fCalorimeter=="EMCAL" && mom.E() > 2*fEMCALCellAmpMin) nClustersInModule[nModule]++;
1681 else if(fCalorimeter=="PHOS" && mom.E() > 2*fPHOSCellAmpMin ) nClustersInModule[nModule]++;
1684 //-----------------------------------------------------------
1685 //Fill histograms related to single cluster or track matching
1686 //-----------------------------------------------------------
1688 // Get Track matched
1691 if(!strcmp("AliESDCaloCluster",Form("%s",clus->ClassName())))
1693 track = dynamic_cast<AliVTrack*> (GetReader()->GetInputEvent()->GetTrack(clus->GetTrackMatchedIndex()));
1697 track = dynamic_cast<AliVTrack*> (clus->GetTrackMatched(0));
1701 ClusterHistograms(mom, pos, nCaloCellsPerCluster, nModule, matched,track, labels, nLabel);
1704 //-----------------------------------------------------------
1706 //-----------------------------------------------------------
1707 if(fFillAllPi0Histo){
1708 if(GetDebug()>1) printf("Invariant mass \n");
1710 //do not do for bad vertex
1711 // Float_t fZvtxCut = 40. ;
1712 if(v[2]<-GetZvertexCut() || v[2]> GetZvertexCut()) continue ; //Event can not be used (vertex, centrality,... cuts not fulfilled)
1714 Int_t nModule2 = -1;
1715 if (nCaloClusters > 1 && nCaloCellsPerCluster > 1) {
1716 for(Int_t jclus = iclus + 1 ; jclus < nCaloClusters ; jclus++) {
1717 AliVCluster* clus2 = (AliVCluster*)caloClusters->At(jclus);
1718 if( clus2->GetNCells() > 1 ){
1720 //Get cluster kinematics
1721 clus2->GetMomentum(mom2,v);
1723 //Check only certain regions
1725 if(IsFiducialCutOn()) in2 = GetFiducialCut()->IsInFiducialCut(mom2,fCalorimeter) ;
1728 //Get module of cluster
1729 nModule2 = GetModuleNumber(clus2);
1731 //Fill invariant mass histograms
1734 fhIM ->Fill((mom+mom2).Pt(),(mom+mom2).M());
1737 if(nModule == nModule2 && nModule >=0 && nModule < fNModules)
1738 fhIMMod[nModule]->Fill((mom+mom2).Pt(),(mom+mom2).M());
1740 //Asymetry histograms
1741 fhAsym->Fill((mom+mom2).Pt(),TMath::Abs((mom.E()-mom2.E())/(mom.E()+mom2.E())));
1743 }// 2nd cluster loop
1744 } // At least one cell in cluster and one cluster in the array
1751 //Number of clusters histograms
1752 if(nCaloClustersAccepted > 0) fhNClusters->Fill(nCaloClustersAccepted);
1753 // Number of clusters per module
1754 for(Int_t imod = 0; imod < fNModules; imod++ ){
1756 printf("AliAnaCalorimeterQA::MakeAnalysisFillHistograms() - module %d calo %s clusters %d\n", imod, fCalorimeter.Data(), nClustersInModule[imod]);
1757 fhNClustersMod->Fill(nClustersInModule[imod],imod);
1759 delete [] nClustersInModule;
1760 //delete caloClusters;
1762 //----------------------------------------------------------
1764 //----------------------------------------------------------
1766 Int_t ncells = cell->GetNumberOfCells();
1769 printf("AliAnaCalorimeterQA::MakeAnalysisFillHistograms() - %s cell entries %d\n", fCalorimeter.Data(), ncells );
1771 //Init arrays and used variables
1772 Int_t *nCellsInModule = new Int_t[fNModules];
1773 for(Int_t imod = 0; imod < fNModules; imod++ ) nCellsInModule[imod] = 0;
1780 Float_t recalF = 1.;
1782 for (Int_t iCell = 0; iCell < cell->GetNumberOfCells(); iCell++) {
1784 printf("AliAnaCalorimeterQA::MakeAnalysisFillHistograms() - Cell : amp %f, absId %d \n", cell->GetAmplitude(iCell), cell->GetCellNumber(iCell));
1785 nModule = GetModuleNumberCellIndexes(cell->GetCellNumber(iCell),fCalorimeter, icol, irow, iRCU);
1787 printf("\t module %d, column %d, row %d \n", nModule,icol,irow);
1789 if(nModule < fNModules) {
1791 //Check if the cell is a bad channel
1792 if(GetCaloUtils()->IsBadChannelsRemovalSwitchedOn()){
1793 if(fCalorimeter=="EMCAL"){
1794 if(GetCaloUtils()->GetEMCALChannelStatus(nModule,icol,irow)) continue;
1797 if(GetCaloUtils()->GetPHOSChannelStatus(nModule,icol,irow)) {
1798 printf("PHOS bad channel\n");
1802 } // use bad channel map
1804 amp = cell->GetAmplitude(iCell)*recalF;
1805 time = cell->GetTime(iCell);
1806 id = cell->GetCellNumber(iCell);
1808 //Get Recalibration factor if set
1809 if (GetCaloUtils()->IsRecalibrationOn()) {
1810 if(fCalorimeter == "PHOS") {
1811 amp *= GetCaloUtils()->GetPHOSChannelRecalibrationFactor(nModule,icol,irow);
1814 amp *= GetCaloUtils()->GetEMCALChannelRecalibrationFactor(nModule,icol,irow);
1815 GetCaloUtils()->GetEMCALRecoUtils()->RecalibrateCellTime(id,bc,time);
1816 // printf("QA : Id %d, Time org %e, Time new %e; Amp org %f, Amp new %f\n",
1817 // id, cell->GetTime(iCell),time, cell->GetAmplitude(iCell),amp);
1819 //if(fCalorimeter == "PHOS")printf("Recalibration factor (sm,row,col)=(%d,%d,%d) - %f\n",nModule,icol,irow,recalF);
1822 //Transform time to ns
1825 //Remove noisy channels, only possible in ESDs
1826 if(GetReader()->GetDataType() == AliCaloTrackReader::kESD){
1827 if(time < fTimeCutMin || time > fTimeCutMax){
1828 printf("Remove cluster with TOF %f\n",time);
1832 //if(amp > 3 && fCalorimeter=="EMCAL") printf("Amp = %f, time = %f, (mod, col, row)= (%d,%d,%d)\n",
1833 // amp,time,nModule,icol,irow);
1835 fhAmplitude->Fill(amp);
1836 fhAmpId ->Fill(amp,id);
1838 if ((fCalorimeter=="EMCAL" && amp > fEMCALCellAmpMin) ||
1839 (fCalorimeter=="PHOS" && amp > fPHOSCellAmpMin)) {
1841 nCellsInModule[nModule]++ ;
1845 if(fCalorimeter=="EMCAL"){
1846 icols = (nModule % 2) ? icol + fNMaxCols : icol;
1847 irows = irow + fNMaxRows * Int_t(nModule / 2);
1850 irows = irow + fNMaxRows * fNModules;
1853 fhGridCellsMod ->Fill(icols,irows);
1854 fhGridCellsEMod->Fill(icols,irows,amp);
1856 if(GetReader()->GetDataType() == AliCaloTrackReader::kESD){
1857 //printf("%s: time %g\n",fCalorimeter.Data(), time);
1858 fhTime ->Fill(time);
1859 fhTimeId ->Fill(time,id);
1860 fhTimeAmp ->Fill(amp,time);
1861 fhGridCellsTimeMod->Fill(icols,irows,time);
1863 fhTimeAmpPerRCU [nModule*fNRCU+iRCU]->Fill(amp, time);
1868 //Get Eta-Phi position of Cell
1869 if(fFillAllPosHisto)
1871 if(fCalorimeter=="EMCAL" && GetCaloUtils()->IsEMCALGeoMatrixSet()){
1872 Float_t celleta = 0.;
1873 Float_t cellphi = 0.;
1874 GetEMCALGeometry()->EtaPhiFromIndex(id, celleta, cellphi);
1876 fhEtaPhiAmp->Fill(celleta,cellphi,amp);
1877 Double_t cellpos[] = {0, 0, 0};
1878 GetEMCALGeometry()->GetGlobal(id, cellpos);
1879 fhXCellE->Fill(cellpos[0],amp) ;
1880 fhYCellE->Fill(cellpos[1],amp) ;
1881 fhZCellE->Fill(cellpos[2],amp) ;
1882 Float_t rcell = TMath::Sqrt(cellpos[0]*cellpos[0]+cellpos[1]*cellpos[1]);//+cellpos[2]*cellpos[2]);
1883 fhRCellE->Fill(rcell,amp) ;
1884 fhXYZCell->Fill(cellpos[0],cellpos[1],cellpos[2]) ;
1886 else if(fCalorimeter=="PHOS" && GetCaloUtils()->IsPHOSGeoMatrixSet()){
1888 Int_t relId[4], module;
1889 Float_t xCell, zCell;
1891 GetPHOSGeometry()->AbsToRelNumbering(id,relId);
1893 GetPHOSGeometry()->RelPosInModule(relId,xCell,zCell);
1894 GetPHOSGeometry()->Local2Global(module,xCell,zCell,xyz);
1895 Float_t rcell = TMath::Sqrt(xyz.X()*xyz.X()+xyz.Y()*xyz.Y());
1896 fhXCellE ->Fill(xyz.X(),amp) ;
1897 fhYCellE ->Fill(xyz.Y(),amp) ;
1898 fhZCellE ->Fill(xyz.Z(),amp) ;
1899 fhRCellE ->Fill(rcell ,amp) ;
1900 fhXYZCell->Fill(xyz.X(),xyz.Y(),xyz.Z()) ;
1902 }//fill cell position histograms
1904 if (fCalorimeter=="EMCAL" && amp > fEMCALCellAmpMin) ncells ++ ;
1905 else if(fCalorimeter=="PHOS" && amp > fPHOSCellAmpMin) ncells ++ ;
1907 // printf("AliAnaCalorimeterQA::MakeAnalysisFillHistograms() - no %s CELLS passed the analysis cut\n",fCalorimeter.Data());
1911 if(ncells > 0 )fhNCells->Fill(ncells) ; //fill the cells after the cut
1913 //Number of cells per module
1914 for(Int_t imod = 0; imod < fNModules; imod++ ) {
1916 printf("AliAnaCalorimeterQA::MakeAnalysisFillHistograms() - module %d calo %s cells %d\n", imod, fCalorimeter.Data(), nCellsInModule[imod]);
1917 fhNCellsMod->Fill(nCellsInModule[imod],imod) ;
1919 delete [] nCellsInModule;
1922 printf("AliAnaCalorimeterQA::MakeAnalysisFillHistograms() - End \n");
1926 //_____________________________________________________________________________________________
1927 void AliAnaCalorimeterQA::ClusterHistograms(const TLorentzVector mom,
1928 Float_t *pos, const Int_t nCaloCellsPerCluster,const Int_t nModule,
1929 const Bool_t matched, const AliVTrack * track,
1930 const Int_t * labels, const Int_t nLabels){
1931 //Fill CaloCluster related histograms
1933 AliAODMCParticle * aodprimary = 0x0;
1934 TParticle * primary = 0x0;
1937 Float_t e = mom.E();
1938 Float_t pt = mom.Pt();
1939 Float_t eta = mom.Eta();
1940 Float_t phi = mom.Phi();
1941 if(phi < 0) phi +=TMath::TwoPi();
1942 if(GetDebug() > 0) {
1943 printf("AliAnaCalorimeterQA::ClusterHistograms() - cluster: E %2.3f, pT %2.3f, eta %2.3f, phi %2.3f \n",e,pt,eta,phi*TMath::RadToDeg());
1945 //printf("\t Primaries: nlabels %d, labels pointer %p\n",nLabels,labels);
1946 printf("\t Primaries: nlabels %d\n",nLabels);
1947 if(!nLabels || !labels) printf("\t Strange, no labels!!!\n");
1952 if(nModule >=0 && nModule < fNModules) fhEMod->Fill(e,nModule);
1959 fhEtaPhiE->Fill(eta,phi,e);
1962 fhNCellsPerCluster ->Fill(e, nCaloCellsPerCluster);
1963 if((fCalorimeter=="EMCAL" && GetReader()->GetEMCALPtMin() < 0.3) ||
1964 (fCalorimeter=="PHOS" && GetReader()->GetPHOSPtMin() < 0.3)) fhNCellsPerClusterMIP->Fill(e, nCaloCellsPerCluster);
1967 if(fFillAllPosHisto2){
1968 fhXE ->Fill(pos[0],e);
1969 fhYE ->Fill(pos[1],e);
1970 fhZE ->Fill(pos[2],e);
1972 fhXYZ ->Fill(pos[0], pos[1],pos[2]);
1974 fhXNCells->Fill(pos[0],nCaloCellsPerCluster);
1975 fhYNCells->Fill(pos[1],nCaloCellsPerCluster);
1976 fhZNCells->Fill(pos[2],nCaloCellsPerCluster);
1977 Float_t rxyz = TMath::Sqrt(pos[0]*pos[0]+pos[1]*pos[1]);//+pos[2]*pos[2]);
1978 fhRE ->Fill(rxyz,e);
1979 fhRNCells->Fill(rxyz ,nCaloCellsPerCluster);
1982 if(nModule >=0 && nModule < fNModules) fhNCellsPerClusterMod[nModule]->Fill(e, nCaloCellsPerCluster);
1984 //Fill histograms only possible when simulation
1985 if(IsDataMC() && nLabels > 0 && labels){
1987 //Play with the MC stack if available
1988 Int_t label = labels[0];
1991 if(GetDebug() >= 0) printf("AliAnaCalorimeterQA::ClusterHistograms() *** bad label ***: label %d \n", label);
1995 Int_t pdg =-1; Int_t pdg0 =-1;Int_t status = -1; Int_t iMother = -1; Int_t iParent = -1;
1996 Float_t vxMC= 0; Float_t vyMC = 0;
1997 Float_t eMC = 0; Float_t ptMC= 0; Float_t phiMC =0; Float_t etaMC = 0;
2001 tag = GetMCAnalysisUtils()->CheckOrigin(labels,nLabels, GetReader(),0);
2003 if(GetReader()->ReadStack() && !GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCUnknown)){ //it MC stack and known tag
2005 if( label >= GetMCStack()->GetNtrack()) {
2006 if(GetDebug() >= 0) printf("AliAnaCalorimeterQA::ClusterHistograms() *** large label ***: label %d, n tracks %d \n", label, GetMCStack()->GetNtrack());
2010 primary = GetMCStack()->Particle(label);
2012 pdg0 = TMath::Abs(primary->GetPdgCode());
2014 status = primary->GetStatusCode();
2015 vxMC = primary->Vx();
2016 vyMC = primary->Vy();
2017 iParent = primary->GetFirstMother();
2019 if(GetDebug() > 1 ) {
2020 printf("AliAnaCalorimeterQA::ClusterHistograms() - Cluster most contributing mother: \n");
2021 printf("\t Mother label %d, pdg %d, %s, status %d, parent %d \n",iMother, pdg0, primary->GetName(),status, iParent);
2024 //Get final particle, no conversion products
2025 if(GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCConversion)){
2027 primary = GetMCStack()->Particle(iParent);
2028 pdg = TMath::Abs(primary->GetPdgCode());
2029 if(GetDebug() > 1 ) printf("AliAnaCalorimeterQA::ClusterHistograms() - Converted cluster!. Find before conversion: \n");
2030 while((pdg == 22 || pdg == 11) && status != 1){
2032 primary = GetMCStack()->Particle(iMother);
2033 status = primary->GetStatusCode();
2034 iParent = primary->GetFirstMother();
2035 pdg = TMath::Abs(primary->GetPdgCode());
2036 if(GetDebug() > 1 )printf("\t pdg %d, index %d, %s, status %d \n",pdg, iMother, primary->GetName(),status);
2039 if(GetDebug() > 1 ) {
2040 printf("AliAnaCalorimeterQA::ClusterHistograms() - Converted Cluster mother before conversion: \n");
2041 printf("\t Mother label %d, pdg %d, %s, status %d, parent %d \n",iMother, pdg, primary->GetName(), status, iParent);
2046 //Overlapped pi0 (or eta, there will be very few), get the meson
2047 if(GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCPi0) ||
2048 GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCEta)){
2049 if(GetDebug() > 1 ) printf("AliAnaCalorimeterQA::ClusterHistograms() - Overlapped Meson decay!, Find it: \n");
2050 while(pdg != 111 && pdg != 221){
2052 primary = GetMCStack()->Particle(iMother);
2053 status = primary->GetStatusCode();
2054 iParent = primary->GetFirstMother();
2055 pdg = TMath::Abs(primary->GetPdgCode());
2056 if(GetDebug() > 1 ) printf("\t pdg %d, %s, index %d\n",pdg, primary->GetName(),iMother);
2058 printf("AliAnaCalorimeterQA::ClusterHistograms() - Tagged as Overlapped photon but meson not found, why?\n");
2063 if(GetDebug() > 2 ) printf("AliAnaCalorimeterQA::ClusterHistograms() - Overlapped %s decay, label %d \n",
2064 primary->GetName(),iMother);
2067 eMC = primary->Energy();
2068 ptMC = primary->Pt();
2069 phiMC = primary->Phi();
2070 etaMC = primary->Eta();
2071 pdg = TMath::Abs(primary->GetPdgCode());
2072 charge = (Int_t) TDatabasePDG::Instance()->GetParticle(pdg)->Charge();
2075 else if(GetReader()->ReadAODMCParticles() && !GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCUnknown)){//it MC AOD and known tag
2076 //Get the list of MC particles
2077 if(!GetReader()->GetAODMCParticles(0))
2078 AliFatal("MCParticles not available!");
2080 aodprimary = (AliAODMCParticle*) (GetReader()->GetAODMCParticles(0))->At(label);
2082 pdg0 = TMath::Abs(aodprimary->GetPdgCode());
2084 status = aodprimary->IsPrimary();
2085 vxMC = aodprimary->Xv();
2086 vyMC = aodprimary->Yv();
2087 iParent = aodprimary->GetMother();
2089 if(GetDebug() > 1 ) {
2090 printf("AliAnaCalorimeterQA::ClusterHistograms() - Cluster most contributing mother: \n");
2091 printf("\t Mother label %d, pdg %d, Primary? %d, Physical Primary? %d, parent %d \n",
2092 iMother, pdg0, aodprimary->IsPrimary(), aodprimary->IsPhysicalPrimary(), iParent);
2095 //Get final particle, no conversion products
2096 if(GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCConversion)){
2098 printf("AliAnaCalorimeterQA::ClusterHistograms() - Converted cluster!. Find before conversion: \n");
2100 aodprimary = (AliAODMCParticle*)(GetReader()->GetAODMCParticles(0))->At(iParent);
2101 pdg = TMath::Abs(aodprimary->GetPdgCode());
2102 while ((pdg == 22 || pdg == 11) && !aodprimary->IsPhysicalPrimary()) {
2104 aodprimary = (AliAODMCParticle*)(GetReader()->GetAODMCParticles(0))->At(iMother);
2105 status = aodprimary->IsPrimary();
2106 iParent = aodprimary->GetMother();
2107 pdg = TMath::Abs(aodprimary->GetPdgCode());
2109 printf("\t pdg %d, index %d, Primary? %d, Physical Primary? %d \n",
2110 pdg, iMother, aodprimary->IsPrimary(), aodprimary->IsPhysicalPrimary());
2113 if(GetDebug() > 1 ) {
2114 printf("AliAnaCalorimeterQA::ClusterHistograms() - Converted Cluster mother before conversion: \n");
2115 printf("\t Mother label %d, pdg %d, parent %d, Primary? %d, Physical Primary? %d \n",
2116 iMother, pdg, iParent, aodprimary->IsPrimary(), aodprimary->IsPhysicalPrimary());
2121 //Overlapped pi0 (or eta, there will be very few), get the meson
2122 if(GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCPi0) ||
2123 GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCEta)){
2124 if(GetDebug() > 1 ) printf("AliAnaCalorimeterQA::ClusterHistograms() - Overlapped Meson decay!, Find it: PDG %d, mom %d \n",pdg, iMother);
2125 while(pdg != 111 && pdg != 221){
2128 aodprimary = (AliAODMCParticle*)(GetReader()->GetAODMCParticles(0))->At(iMother);
2129 status = aodprimary->IsPrimary();
2130 iParent = aodprimary->GetMother();
2131 pdg = TMath::Abs(aodprimary->GetPdgCode());
2133 if(GetDebug() > 1 ) printf("\t pdg %d, index %d\n",pdg, iMother);
2136 printf("AliAnaCalorimeterQA::ClusterHistograms() - Tagged as Overlapped photon but meson not found, why?\n");
2141 if(GetDebug() > 2 ) printf("AliAnaCalorimeterQA::ClusterHistograms() - Overlapped %s decay, label %d \n",
2142 aodprimary->GetName(),iMother);
2145 status = aodprimary->IsPrimary();
2146 eMC = aodprimary->E();
2147 ptMC = aodprimary->Pt();
2148 phiMC = aodprimary->Phi();
2149 etaMC = aodprimary->Eta();
2150 pdg = TMath::Abs(aodprimary->GetPdgCode());
2151 charge = aodprimary->Charge();
2155 //Float_t vz = primary->Vz();
2156 Float_t rVMC = TMath::Sqrt(vxMC*vxMC + vyMC*vyMC);
2157 if((pdg == 22 || TMath::Abs(pdg)==11) && status!=1) {
2158 fhEMVxyz ->Fill(vxMC,vyMC);//,vz);
2159 fhEMR ->Fill(e,rVMC);
2162 //printf("reco e %f, pt %f, phi %f, eta %f \n", e, pt, phi, eta);
2163 //printf("prim e %f, pt %f, phi %f, eta %f \n", eMC,ptMC,phiMC ,etaMC );
2164 //printf("vertex: vx %f, vy %f, vz %f, r %f \n", vxMC, vyMC, vz, r);
2167 //Overlapped pi0 (or eta, there will be very few)
2168 if(GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCPi0)){
2169 fhRecoMCE [mcPi0][matched] ->Fill(e,eMC);
2170 if(e > 0.5 && eMC > 0.5) fhRecoMCEta[mcPi0][(matched)]->Fill(eta,etaMC);
2171 if(e > 0.5 && eMC > 0.5) fhRecoMCPhi[mcPi0][(matched)]->Fill(phi,phiMC);
2172 fhRecoMCDeltaE [mcPi0][(matched)]->Fill(e,eMC-e);
2173 fhRecoMCDeltaPhi[mcPi0][(matched)]->Fill(e,phiMC-phi);
2174 fhRecoMCDeltaEta[mcPi0][(matched)]->Fill(e,etaMC-eta);
2175 }//Overlapped pizero decay
2176 else if(GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCEta)){
2177 fhRecoMCE [mcEta][(matched)] ->Fill(e,eMC);
2178 if(e > 0.5 && eMC > 0.5) fhRecoMCEta[mcEta][(matched)]->Fill(eta,etaMC);
2179 if(e > 0.5 && eMC > 0.5) fhRecoMCPhi[mcEta][(matched)]->Fill(phi,phiMC);
2180 fhRecoMCDeltaE [mcEta][(matched)]->Fill(e,eMC-e);
2181 fhRecoMCDeltaPhi[mcEta][(matched)]->Fill(e,phiMC-phi);
2182 fhRecoMCDeltaEta[mcEta][(matched)]->Fill(e,etaMC-eta);
2183 }//Overlapped eta decay
2184 else if(GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCPhoton)){
2185 fhRecoMCE [mcPhoton][(matched)] ->Fill(e,eMC);
2186 if(e > 0.5 && eMC > 0.5) fhRecoMCEta[mcPhoton][(matched)]->Fill(eta,etaMC);
2187 if(e > 0.5 && eMC > 0.5) fhRecoMCPhi[mcPhoton][(matched)]->Fill(phi,phiMC);
2188 fhRecoMCDeltaE [mcPhoton][(matched)]->Fill(e,eMC-e);
2189 fhRecoMCDeltaPhi[mcPhoton][(matched)]->Fill(e,phiMC-phi);
2190 fhRecoMCDeltaEta[mcPhoton][(matched)]->Fill(e,etaMC-eta);
2192 else if(GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCElectron)){
2193 fhRecoMCE [mcElectron][(matched)] ->Fill(e,eMC);
2194 if(e > 0.5 && eMC > 0.5) fhRecoMCEta[mcElectron][(matched)]->Fill(eta,etaMC);
2195 if(e > 0.5 && eMC > 0.5) fhRecoMCPhi[mcElectron][(matched)]->Fill(phi,phiMC);
2196 fhRecoMCDeltaE [mcElectron][(matched)]->Fill(e,eMC-e);
2197 fhRecoMCDeltaPhi[mcElectron][(matched)]->Fill(e,phiMC-phi);
2198 fhRecoMCDeltaEta[mcElectron][(matched)]->Fill(e,etaMC-eta);
2199 fhEMVxyz ->Fill(vxMC,vyMC);//,vz);
2200 fhEMR ->Fill(e,rVMC);
2202 else if(charge == 0){
2203 fhRecoMCE [mcNeHadron][(matched)] ->Fill(e,eMC);
2204 if(e > 0.5 && eMC > 0.5) fhRecoMCEta[mcNeHadron][(matched)]->Fill(eta,etaMC);
2205 if(e > 0.5 && eMC > 0.5) fhRecoMCPhi[mcNeHadron][(matched)]->Fill(phi,phiMC);
2206 fhRecoMCDeltaE [mcNeHadron][(matched)]->Fill(e,eMC-e);
2207 fhRecoMCDeltaPhi[mcNeHadron][(matched)]->Fill(e,phiMC-phi);
2208 fhRecoMCDeltaEta[mcNeHadron][(matched)]->Fill(e,etaMC-eta);
2209 fhHaVxyz ->Fill(vxMC,vyMC);//,vz);
2210 fhHaR ->Fill(e,rVMC);
2213 fhRecoMCE [mcChHadron][(matched)] ->Fill(e,eMC);
2214 if(e > 0.5 && eMC > 0.5) fhRecoMCEta[mcChHadron][(matched)]->Fill(eta,etaMC);
2215 if(e > 0.5 && eMC > 0.5) fhRecoMCPhi[mcChHadron][(matched)]->Fill(phi,phiMC);
2216 fhRecoMCDeltaE [mcChHadron][(matched)]->Fill(e,eMC-e);
2217 fhRecoMCDeltaPhi[mcChHadron][(matched)]->Fill(e,phiMC-phi);
2218 fhRecoMCDeltaEta[mcChHadron][(matched)]->Fill(e,etaMC-eta);
2219 fhHaVxyz ->Fill(vxMC,vyMC);//,vz);
2220 fhHaR ->Fill(e,rVMC);
2224 //Match tracks and clusters
2225 //To be Modified in case of AODs
2227 if( matched && fFillAllTMHisto){
2228 if(fFillAllTH12 && fFillAllTMHisto){
2229 fhECharged ->Fill(e);
2230 fhPtCharged ->Fill(pt);
2231 fhPhiCharged ->Fill(phi);
2232 fhEtaCharged ->Fill(eta);
2235 if(fFillAllTMHisto){
2236 if(fFillAllTH3)fhEtaPhiECharged->Fill(eta,phi,e);
2237 if((fCalorimeter=="EMCAL" && GetReader()->GetEMCALPtMin() < 0.3) ||
2238 (fCalorimeter=="PHOS" && GetReader()->GetPHOSPtMin() < 0.3)) fhNCellsPerClusterMIPCharged->Fill(e, nCaloCellsPerCluster);
2240 //printf("track index %d ntracks %d\n", esd->GetNumberOfTracks());
2242 //Study the track and matched cluster if track exists.
2244 Double_t emcpos[3] = {0.,0.,0.};
2245 Double_t emcmom[3] = {0.,0.,0.};
2246 Double_t radius = 441.0; //[cm] EMCAL radius +13cm
2247 Double_t bfield = 0.;
2253 Double_t tpcSignal = 0;
2254 Bool_t okpos = kFALSE;
2255 Bool_t okmom = kFALSE;
2256 Bool_t okout = kFALSE;
2260 //In case of ESDs get the parameters in this way
2261 if(GetReader()->GetDataType()==AliCaloTrackReader::kESD) {
2262 if (track->GetOuterParam() ) {
2265 bfield = GetReader()->GetInputEvent()->GetMagneticField();
2266 okpos = track->GetOuterParam()->GetXYZAt(radius,bfield,emcpos);
2267 okmom = track->GetOuterParam()->GetPxPyPzAt(radius,bfield,emcmom);
2268 if(!(okpos && okmom)) return;
2270 TVector3 position(emcpos[0],emcpos[1],emcpos[2]);
2271 TVector3 momentum(emcmom[0],emcmom[1],emcmom[2]);
2272 tphi = position.Phi();
2273 teta = position.Eta();
2274 tmom = momentum.Mag();
2278 tpcSignal = track->GetTPCsignal();
2280 nITS = track->GetNcls(0);
2281 nTPC = track->GetNcls(1);
2282 }//Outer param available
2284 else if(GetReader()->GetDataType()==AliCaloTrackReader::kAOD) {
2285 AliAODPid* pid = (AliAODPid*) ((AliAODTrack *) track)->GetDetPid();
2288 pid->GetEMCALPosition(emcpos);
2289 pid->GetEMCALMomentum(emcmom);
2291 TVector3 position(emcpos[0],emcpos[1],emcpos[2]);
2292 TVector3 momentum(emcmom[0],emcmom[1],emcmom[2]);
2293 tphi = position.Phi();
2294 teta = position.Eta();
2295 tmom = momentum.Mag();
2299 tpcSignal = pid->GetTPCsignal();
2305 //printf("okout\n");
2306 Double_t deta = teta - eta;
2307 Double_t dphi = tphi - phi;
2308 if(dphi > TMath::Pi()) dphi -= 2*TMath::Pi();
2309 if(dphi < -TMath::Pi()) dphi += 2*TMath::Pi();
2310 Double_t dR = sqrt(dphi*dphi + deta*deta);
2312 Double_t pOverE = tmom/e;
2314 fh1pOverE->Fill(tpt, pOverE);
2315 if(dR < 0.02) fh1pOverER02->Fill(tpt,pOverE);
2318 fh2MatchdEdx->Fill(tmom2,tpcSignal);
2320 if(IsDataMC() && primary){
2321 Int_t pdg = primary->GetPdgCode();
2322 Double_t charge = TDatabasePDG::Instance()->GetParticle(pdg)->Charge();
2324 if(TMath::Abs(pdg) == 11){
2325 fhMCEle1pOverE->Fill(tpt,pOverE);
2326 fhMCEle1dR->Fill(dR);
2327 fhMCEle2MatchdEdx->Fill(tmom2,tpcSignal);
2328 if(dR < 0.02) fhMCEle1pOverER02->Fill(tpt,pOverE);
2331 fhMCChHad1pOverE->Fill(tpt,pOverE);
2332 fhMCChHad1dR->Fill(dR);
2333 fhMCChHad2MatchdEdx->Fill(tmom2,tpcSignal);
2334 if(dR < 0.02) fhMCChHad1pOverER02->Fill(tpt,pOverE);
2336 else if(charge == 0){
2337 fhMCNeutral1pOverE->Fill(tpt,pOverE);
2338 fhMCNeutral1dR->Fill(dR);
2339 fhMCNeutral2MatchdEdx->Fill(tmom2,tpcSignal);
2340 if(dR < 0.02) fhMCNeutral1pOverER02->Fill(tpt,pOverE);
2344 if(dR < 0.02 && pOverE > 0.5 && pOverE < 1.5
2345 && nCaloCellsPerCluster > 1 && nITS > 3 && nTPC > 20) {
2346 fh2EledEdx->Fill(tmom2,tpcSignal);
2349 else{//no ESD external param or AODPid
2351 if(GetDebug() >= 0) printf("No ESD external param or AliAODPid \n");
2354 }//matched clusters with tracks
2359 //__________________________________
2360 void AliAnaCalorimeterQA::Correlate(){
2361 // Correlate information from PHOS and EMCAL and with V0 and track multiplicity
2364 TObjArray * caloClustersEMCAL = GetEMCALClusters();
2365 TObjArray * caloClustersPHOS = GetPHOSClusters();
2367 Int_t nclEMCAL = caloClustersEMCAL->GetEntriesFast();
2368 Int_t nclPHOS = caloClustersPHOS ->GetEntriesFast();
2370 Float_t sumClusterEnergyEMCAL = 0;
2371 Float_t sumClusterEnergyPHOS = 0;
2373 for(iclus = 0 ; iclus < caloClustersEMCAL->GetEntriesFast() ; iclus++)
2374 sumClusterEnergyEMCAL += ((AliVCluster*)caloClustersEMCAL->At(iclus))->E();
2375 for(iclus = 0 ; iclus < caloClustersPHOS->GetEntriesFast(); iclus++)
2376 sumClusterEnergyPHOS += ((AliVCluster*)caloClustersPHOS->At(iclus))->E();
2381 AliVCaloCells * cellsEMCAL = GetEMCALCells();
2382 AliVCaloCells * cellsPHOS = GetPHOSCells();
2384 Int_t ncellsEMCAL = cellsEMCAL->GetNumberOfCells();
2385 Int_t ncellsPHOS = cellsPHOS ->GetNumberOfCells();
2387 Float_t sumCellEnergyEMCAL = 0;
2388 Float_t sumCellEnergyPHOS = 0;
2390 for(icell = 0 ; icell < cellsEMCAL->GetNumberOfCells() ; icell++)
2391 sumCellEnergyEMCAL += cellsEMCAL->GetAmplitude(icell);
2392 for(icell = 0 ; icell < cellsPHOS->GetNumberOfCells(); icell++)
2393 sumCellEnergyPHOS += cellsPHOS->GetAmplitude(icell);
2397 fhCaloCorrNClusters->Fill(nclEMCAL,nclPHOS);
2398 fhCaloCorrEClusters->Fill(sumClusterEnergyEMCAL,sumClusterEnergyPHOS);
2399 fhCaloCorrNCells ->Fill(ncellsEMCAL,ncellsPHOS);
2400 fhCaloCorrECells ->Fill(sumCellEnergyEMCAL,sumCellEnergyPHOS);
2402 Int_t v0S = GetV0Signal(0)+GetV0Signal(1);
2403 Int_t v0M = GetV0Multiplicity(0)+GetV0Multiplicity(1);
2404 Int_t trM = GetTrackMultiplicity();
2405 if(fCalorimeter=="PHOS"){
2406 fhCaloV0MCorrNClusters ->Fill(v0M,nclPHOS);
2407 fhCaloV0MCorrEClusters ->Fill(v0M,sumClusterEnergyPHOS);
2408 fhCaloV0MCorrNCells ->Fill(v0M,ncellsPHOS);
2409 fhCaloV0MCorrECells ->Fill(v0M,sumCellEnergyPHOS);
2411 fhCaloV0SCorrNClusters ->Fill(v0S,nclPHOS);
2412 fhCaloV0SCorrEClusters ->Fill(v0S,sumClusterEnergyPHOS);
2413 fhCaloV0SCorrNCells ->Fill(v0S,ncellsPHOS);
2414 fhCaloV0SCorrECells ->Fill(v0S,sumCellEnergyPHOS);
2416 fhCaloTrackMCorrNClusters->Fill(trM,nclPHOS);
2417 fhCaloTrackMCorrEClusters->Fill(trM,sumClusterEnergyPHOS);
2418 fhCaloTrackMCorrNCells ->Fill(trM,ncellsPHOS);
2419 fhCaloTrackMCorrECells ->Fill(trM,sumCellEnergyPHOS);
2422 fhCaloV0MCorrNClusters ->Fill(v0M,nclEMCAL);
2423 fhCaloV0MCorrEClusters ->Fill(v0M,sumClusterEnergyEMCAL);
2424 fhCaloV0MCorrNCells ->Fill(v0M,ncellsEMCAL);
2425 fhCaloV0MCorrECells ->Fill(v0M,sumCellEnergyEMCAL);
2427 fhCaloV0SCorrNClusters ->Fill(v0S,nclEMCAL);
2428 fhCaloV0SCorrEClusters ->Fill(v0S,sumClusterEnergyEMCAL);
2429 fhCaloV0SCorrNCells ->Fill(v0S,ncellsEMCAL);
2430 fhCaloV0SCorrECells ->Fill(v0S,sumCellEnergyEMCAL);
2432 fhCaloTrackMCorrNClusters->Fill(trM,nclEMCAL);
2433 fhCaloTrackMCorrEClusters->Fill(trM,sumClusterEnergyEMCAL);
2434 fhCaloTrackMCorrNCells ->Fill(trM,ncellsEMCAL);
2435 fhCaloTrackMCorrECells ->Fill(trM,sumCellEnergyEMCAL);
2440 printf("AliAnaCalorimeterQA::Correlate(): \n");
2441 printf("\t EMCAL: N cells %d, N clusters %d, summed E cells %f, summed E clusters %f \n",
2442 ncellsEMCAL,nclEMCAL, sumCellEnergyEMCAL,sumClusterEnergyEMCAL);
2443 printf("\t PHOS : N cells %d, N clusters %d, summed E cells %f, summed E clusters %f \n",
2444 ncellsPHOS,nclPHOS,sumCellEnergyPHOS,sumClusterEnergyPHOS);
2445 printf("\t V0 : Signal %d, Multiplicity %d, Track Multiplicity %d \n", v0S,v0M,trM);
2449 //______________________________________________________________________________
2450 void AliAnaCalorimeterQA::FillCellPositionHistograms(const Int_t nCaloCellsPerCluster,const UShort_t * indexList,
2451 const Float_t pos[3], const Float_t clEnergy){
2452 // Fill histograms releated to cell position
2454 //Loop on cluster cells
2455 for (Int_t ipos = 0; ipos < nCaloCellsPerCluster; ipos++) {
2457 // printf("Index %d\n",ipos);
2458 Int_t absId = indexList[ipos];
2460 //Get position of cell compare to cluster
2462 if(fCalorimeter=="EMCAL" && GetCaloUtils()->IsEMCALGeoMatrixSet()){
2464 Double_t cellpos[] = {0, 0, 0};
2465 GetEMCALGeometry()->GetGlobal(absId, cellpos);
2467 fhDeltaCellClusterXNCells->Fill(pos[0]-cellpos[0],nCaloCellsPerCluster) ;
2468 fhDeltaCellClusterYNCells->Fill(pos[1]-cellpos[1],nCaloCellsPerCluster) ;
2469 fhDeltaCellClusterZNCells->Fill(pos[2]-cellpos[2],nCaloCellsPerCluster) ;
2471 fhDeltaCellClusterXE->Fill(pos[0]-cellpos[0],clEnergy) ;
2472 fhDeltaCellClusterYE->Fill(pos[1]-cellpos[1],clEnergy) ;
2473 fhDeltaCellClusterZE->Fill(pos[2]-cellpos[2],clEnergy) ;
2475 Float_t r = TMath::Sqrt(pos[0] *pos[0] + pos[1] * pos[1] );
2476 Float_t rcell = TMath::Sqrt(cellpos[0]*cellpos[0] + cellpos[1]* cellpos[1]);
2478 fhDeltaCellClusterRNCells->Fill(r-rcell, nCaloCellsPerCluster) ;
2479 fhDeltaCellClusterRE ->Fill(r-rcell, clEnergy) ;
2481 }//EMCAL and its matrices are available
2482 else if(fCalorimeter=="PHOS" && GetCaloUtils()->IsPHOSGeoMatrixSet()){
2484 Int_t relId[4], module;
2485 Float_t xCell, zCell;
2487 GetPHOSGeometry()->AbsToRelNumbering(absId,relId);
2489 GetPHOSGeometry()->RelPosInModule(relId,xCell,zCell);
2490 GetPHOSGeometry()->Local2Global(module,xCell,zCell,xyz);
2492 fhDeltaCellClusterXNCells->Fill(pos[0]-xyz.X(),nCaloCellsPerCluster) ;
2493 fhDeltaCellClusterYNCells->Fill(pos[1]-xyz.Y(),nCaloCellsPerCluster) ;
2494 fhDeltaCellClusterZNCells->Fill(pos[2]-xyz.Z(),nCaloCellsPerCluster) ;
2496 fhDeltaCellClusterXE->Fill(pos[0]-xyz.X(),clEnergy) ;
2497 fhDeltaCellClusterYE->Fill(pos[1]-xyz.Y(),clEnergy) ;
2498 fhDeltaCellClusterZE->Fill(pos[2]-xyz.Z(),clEnergy) ;
2500 Float_t r = TMath::Sqrt(pos[0] * pos[0] + pos[1] * pos[1] );
2501 Float_t rcell = TMath::Sqrt(xyz.X() * xyz.X() + xyz.Y() * xyz.Y());
2503 fhDeltaCellClusterRNCells->Fill(r-rcell, nCaloCellsPerCluster) ;
2504 fhDeltaCellClusterRE ->Fill(r-rcell, clEnergy) ;
2506 }//PHOS and its matrices are available
2507 }// cluster cell loop
2508 }//Fill all position histograms
2511 //______________________________________________________________________________
2512 void AliAnaCalorimeterQA::MCHistograms(const TLorentzVector mom, const Int_t pdg){
2513 //Fill pure monte carlo related histograms
2515 Float_t eMC = mom.E();
2516 Float_t phiMC = mom.Phi();
2518 phiMC += TMath::TwoPi();
2519 Float_t etaMC = mom.Eta();
2521 if (TMath::Abs(etaMC) > 1) return;
2525 //Rough stimate of acceptance for pi0, Eta and electrons
2526 if(fCalorimeter == "PHOS"){
2528 if(GetFiducialCut()->IsInFiducialCut(mom,fCalorimeter))
2530 if(GetDebug() > 2) printf("AliAnaCalorimeterQA::MCHistograms() - In %s fiducial cut acceptance? %d\n",fCalorimeter.Data(),in);
2533 else if(fCalorimeter == "EMCAL" && GetCaloUtils()->IsEMCALGeoMatrixSet()){
2534 if(GetEMCALGeometry()){
2537 GetEMCALGeometry()->GetAbsCellIdFromEtaPhi(mom.Eta(),mom.Phi(),absID);
2542 if(GetDebug() > 2) printf("AliAnaCalorimeterQA::MCHistograms() - In %s Real acceptance? %d\n",fCalorimeter.Data(),in);
2545 if(GetFiducialCut()->IsInFiducialCut(mom,fCalorimeter))
2547 if(GetDebug() > 2) printf("AliAnaCalorimeterQA::MCHistograms() - In %s fiducial cut acceptance? %d\n",fCalorimeter.Data(),in);
2552 fhGenMCE[mcPhoton] ->Fill(eMC);
2553 if(eMC > 0.5) fhGenMCEtaPhi[mcPhoton]->Fill(etaMC,phiMC);
2555 fhGenMCAccE[mcPhoton] ->Fill(eMC);
2556 if(eMC > 0.5) fhGenMCAccEtaPhi[mcPhoton]->Fill(etaMC,phiMC);
2559 else if (pdg==111) {
2560 fhGenMCE[mcPi0] ->Fill(eMC);
2561 if(eMC > 0.5) fhGenMCEtaPhi[mcPi0]->Fill(etaMC,phiMC);
2563 fhGenMCAccE[mcPi0] ->Fill(eMC);
2564 if(eMC > 0.5) fhGenMCAccEtaPhi[mcPi0]->Fill(etaMC,phiMC);
2567 else if (pdg==221) {
2568 fhGenMCE[mcEta] ->Fill(eMC);
2569 if(eMC > 0.5) fhGenMCEtaPhi[mcEta]->Fill(etaMC,phiMC);
2571 fhGenMCAccE[mcEta] ->Fill(eMC);
2572 if(eMC > 0.5) fhGenMCAccEtaPhi[mcEta]->Fill(etaMC,phiMC);
2575 else if (TMath::Abs(pdg)==11) {
2576 fhGenMCE[mcElectron] ->Fill(eMC);
2577 if(eMC > 0.5) fhGenMCEtaPhi[mcElectron]->Fill(etaMC,phiMC);
2579 fhGenMCAccE[mcElectron] ->Fill(eMC);
2580 if(eMC > 0.5) fhGenMCAccEtaPhi[mcElectron]->Fill(etaMC,phiMC);