AliPhysics  1c9c77b (1c9c77b)
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
AliAnaCalorimeterQA.cxx
Go to the documentation of this file.
1 /**************************************************************************
2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3  * *
4  * Author: The ALICE Off-line Project. *
5  * Contributors are mentioned in the code where appropriate. *
6  * *
7  * Permission to use, copy, modify and distribute this software and its *
8  * documentation strictly for non-commercial purposes is hereby granted *
9  * without fee, provided that the above copyright notice appears in all *
10  * copies and that both the copyright notice and this permission notice *
11  * appear in the supporting documentation. The authors make no claims *
12  * about the suitability of this software for any purpose. It is *
13  * provided "as is" without express or implied warranty. *
14  **************************************************************************/
15 
16 // --- ROOT system ---
17 #include <TObjArray.h>
18 #include <TParticle.h>
19 #include <TDatabasePDG.h>
20 #include <TH3F.h>
21 #include <TObjString.h>
22 
23 //---- AliRoot system ----
24 #include "AliAnaCalorimeterQA.h"
25 #include "AliCaloTrackReader.h"
26 #include "AliStack.h"
27 #include "AliVCaloCells.h"
28 #include "AliFiducialCut.h"
29 #include "AliVCluster.h"
30 #include "AliVTrack.h"
31 #include "AliVEvent.h"
32 #include "AliVEventHandler.h"
33 #include "AliAODMCParticle.h"
34 #include "AliMCAnalysisUtils.h"
35 
36 // --- Detectors ---
37 #include "AliPHOSGeoUtils.h"
38 #include "AliEMCALGeometry.h"
39 
41 ClassImp(AliAnaCalorimeterQA) ;
43 
44 //__________________________________________
46 //__________________________________________
47 AliAnaCalorimeterQA::AliAnaCalorimeterQA() :
48 AliAnaCaloTrackCorrBaseClass(),
49 
50 // Switches
51 fFillAllCellTimeHisto(kTRUE),
52 fFillAllPosHisto(kFALSE), fFillAllPosHisto2(kFALSE),
53 fFillAllTH3(kFALSE),
54 fFillAllTMHisto(kTRUE),
55 fFillAllPi0Histo(kTRUE), fFillInvMassOpenAngle(kFALSE),
56 fFillPi0PairDiffTime(kFALSE), fFillInvMassInEMCALWithPHOSDCalAcc(kFALSE),
57 fFillEBinAcceptanceHisto(kFALSE),
58 fCorrelate(kTRUE), fStudyBadClusters(kFALSE),
59 fStudyClustersAsymmetry(kFALSE), fStudyExotic(kFALSE),
60 fStudyWeight(kFALSE),
61 
62 // Parameters and cuts
63 fNModules(12), fNRCU(2),
64 fNMaxCols(48), fNMaxRows(24),
65 fTimeCutMin(-10000), fTimeCutMax(10000),
66 fCellAmpMin(0),
67 fEMCALCellAmpMin(0), fPHOSCellAmpMin(0),
68 fEMCALClusterNCellMin(0), fPHOSClusterNCellMin(0),
69 fNEBinCuts(0),
70 
71 // Invariant mass
72 fInvMassMinECut(0), fInvMassMaxECut(0),
73 fInvMassMinM02Cut(0), fInvMassMaxM02Cut(0),
74 fInvMassMaxOpenAngle(0), fInvMassMaxTimeDifference(0),
75 
76 // Exotic
77 fExoNECrossCuts(0), fExoECrossCuts(),
78 fExoNDTimeCuts(0), fExoDTimeCuts(),
79 
80 fClusterMomentum(), fClusterMomentum2(),
81 fPrimaryMomentum(),
82 fConstantTimeShift(0),
83 // Histograms
84 fhE(0), fhPt(0),
85 fhPhi(0), fhEta(0),
86 fhEtaPhi(0), fhEtaPhiE(0),
87 fhECharged(0), fhPtCharged(0),
88 fhPhiCharged(0), fhEtaCharged(0),
89 fhEtaPhiCharged(0), fhEtaPhiECharged(0),
90 
91 // Invariant mass
92 fhIM(0), fhIMSame(0), fhIMDiff(0),
93 fhIMDCAL(0), fhIMDCALSame(0), fhIMDCALDiff(0),
94 fhIMDCALPHOS(0), fhIMDCALPHOSSame(0),
95 fhIMEMCALPHOS(0), fhIMEMCALPHOSSame(0),
96 fhAsym(0),
97 fhOpAngle(0), fhIMvsOpAngle(0),
98 fhNCellsPerCluster(0), fhNCellsPerClusterNoCut(0), fhNClusters(0),
99 
100 // Timing
101 fhClusterTimeEnergy(0), fhCellTimeSpreadRespectToCellMax(0),
102 fhCellIdCellLargeTimeSpread(0), fhClusterPairDiffTimeE(0), fhClusterPairDiffTimeESameMod(0),
103 fhClusterMaxCellCloseCellRatio(0), fhClusterMaxCellCloseCellDiff(0),
104 fhClusterMaxCellDiff(0), fhClusterMaxCellDiffNoCut(0),
105 //fhClusterMaxCellDiffAverageTime(0), fhClusterMaxCellDiffWeightedTime(0),
106 fhClusterMaxCellECross(0),
107 fhLambda0(0), fhLambda1(0), // fhDispersion(0),
108 
109 // bad clusters
110 fhBadClusterEnergy(0), fhBadClusterTimeEnergy(0), fhBadClusterEtaPhi(0),
111 fhBadClusterPairDiffTimeE(0), fhBadCellTimeSpreadRespectToCellMax(0),
112 fhBadClusterMaxCellCloseCellRatio(0), fhBadClusterMaxCellCloseCellDiff(0), fhBadClusterMaxCellDiff(0),
113 fhBadClusterMaxCellDiffAverageTime(0), fhBadClusterMaxCellDiffWeightedTime(0),
114 fhBadClusterMaxCellECross(0),
115 fhBadClusterLambda0(0), fhBadClusterLambda1(0),
116 fhBadClusterDeltaIEtaDeltaIPhiE0(0), fhBadClusterDeltaIEtaDeltaIPhiE2(0),
117 fhBadClusterDeltaIEtaDeltaIPhiE6(0), fhBadClusterDeltaIA(0),
118 
119 // Position
120 fhRNCells(0), fhXNCells(0),
121 fhYNCells(0), fhZNCells(0),
122 fhRE(0), fhXE(0),
123 fhYE(0), fhZE(0),
124 fhXYZ(0),
125 fhRCellE(0), fhXCellE(0),
126 fhYCellE(0), fhZCellE(0),
127 fhXYZCell(0),
128 fhDeltaCellClusterRNCells(0), fhDeltaCellClusterXNCells(0),
129 fhDeltaCellClusterYNCells(0), fhDeltaCellClusterZNCells(0),
130 fhDeltaCellClusterRE(0), fhDeltaCellClusterXE(0),
131 fhDeltaCellClusterYE(0), fhDeltaCellClusterZE(0),
132 
133 // Cells
134 fhNCells(0), fhNCellsCutAmpMin(0),
135 fhAmplitude(0), fhAmpId(0),
136 fhEtaPhiAmpCell(0), fhEtaPhiCell(0),
137 fhTime(0), //fhTimeVz(0),
138 fhTimeId(0), fhTimeAmp(0),
139 fhAmpIdLowGain(0), fhTimeIdLowGain(0), fhTimeAmpLowGain(0),
140 
141 fhCellECross(0),
142 
143 fhEMCALPHOSCorrNClusters(0), fhEMCALPHOSCorrEClusters(0),
144 fhEMCALPHOSCorrNCells(0), fhEMCALPHOSCorrECells(0),
145 fhEMCALDCALCorrNClusters(0), fhEMCALDCALCorrEClusters(0),
146 fhEMCALDCALCorrNCells(0), fhEMCALDCALCorrECells(0),
147 fhDCALPHOSCorrNClusters(0), fhDCALPHOSCorrEClusters(0),
148 fhDCALPHOSCorrNCells(0), fhDCALPHOSCorrECells(0),
149 fhCaloV0SCorrNClusters(0), fhCaloV0SCorrEClusters(0),
150 fhCaloV0SCorrNCells(0), fhCaloV0SCorrECells(0),
151 fhCaloV0MCorrNClusters(0), fhCaloV0MCorrEClusters(0),
152 fhCaloV0MCorrNCells(0), fhCaloV0MCorrECells(0),
153 fhCaloTrackMCorrNClusters(0), fhCaloTrackMCorrEClusters(0),
154 fhCaloTrackMCorrNCells(0), fhCaloTrackMCorrECells(0),
155 fhCaloCenNClusters(0), fhCaloCenEClusters(0),
156 fhCaloCenNCells(0), fhCaloCenECells(0),
157 fhCaloEvPNClusters(0), fhCaloEvPEClusters(0),
158 fhCaloEvPNCells(0), fhCaloEvPECells(0),
159 
160 // Super-Module dependent histograms
161 fhEMod(0), fhAmpMod(0),
162 fhEWeirdMod(0), fhAmpWeirdMod(0),
163 fhTimeMod(0),
164 fhNClustersMod(0), fhNCellsMod(0),
165 fhNCellsSumAmpPerMod(0), fhNClustersSumEnergyPerMod(0),
166 fhNCellsPerClusterMod(0), fhNCellsPerClusterModNoCut(0),
167 fhNCellsPerClusterWeirdMod(0), fhNCellsPerClusterWeirdModNoCut(0),
168 
169 fhGridCells(0), fhGridCellsE(0), fhGridCellsTime(0),
170 fhGridCellsLowGain(0), fhGridCellsELowGain(0), fhGridCellsTimeLowGain(0),
171 fhTimeAmpPerRCU(0), fhIMMod(0),
172 
173 // Weight studies
174 fhECellClusterRatio(0), fhECellClusterLogRatio(0),
175 fhEMaxCellClusterRatio(0), fhEMaxCellClusterLogRatio(0),
176 fhECellTotalRatio(0), fhECellTotalLogRatio(0),
177 fhECellTotalRatioMod(0), fhECellTotalLogRatioMod(0),
178 
179 fhExoL0ECross(0), fhExoL1ECross(0),
180 
181 // MC and reco
182 fhRecoMCE(), fhRecoMCPhi(), fhRecoMCEta(),
183 fhRecoMCDeltaE(), fhRecoMCRatioE(),
184 fhRecoMCDeltaPhi(), fhRecoMCDeltaEta(),
185 
186 // MC only
187 fhGenMCE(), fhGenMCPt(), fhGenMCEtaPhi(),
188 fhGenMCAccE(), fhGenMCAccPt(), fhGenMCAccEtaPhi(),
189 
190 // Matched MC
191 fhEMVxyz(0), fhEMR(0),
192 fhHaVxyz(0), fhHaR(0),
193 fh1EOverP(0), fh2dR(0),
194 fh2EledEdx(0), fh2MatchdEdx(0),
195 fh1EOverPR02(0), fh1EleEOverP(0),
196 fhMCEle1EOverP(0), fhMCEle1dR(0), fhMCEle2MatchdEdx(0),
197 fhMCChHad1EOverP(0), fhMCChHad1dR(0), fhMCChHad2MatchdEdx(0),
198 fhMCNeutral1EOverP(0), fhMCNeutral1dR(0), fhMCNeutral2MatchdEdx(0),
199 fhMCEle1EOverPR02(0), fhMCChHad1EOverPR02(0), fhMCNeutral1EOverPR02(0),
200 fhMCEle1EleEOverP(0), fhMCChHad1EleEOverP(0), fhMCNeutral1EleEOverP(0),
201 fhTrackMatchedDEtaNeg(0), fhTrackMatchedDPhiNeg(0), fhTrackMatchedDEtaDPhiNeg(0),
202 fhTrackMatchedDEtaPos(0), fhTrackMatchedDPhiPos(0), fhTrackMatchedDEtaDPhiPos(0),
203 fhTrackMatchedDEtaNegMod(0), fhTrackMatchedDPhiNegMod(0),
204 fhTrackMatchedDEtaPosMod(0), fhTrackMatchedDPhiPosMod(0)
205 {
206  // Weight studies
207  for(Int_t i =0; i < 12; i++)
208  {
209  for(Int_t j = 0; j < 4; j++)
210  {
211  for(Int_t k = 0; k < 3; k++)
212  {
213  fhLambda0ForW0AndCellCuts [i][j][k] = 0;
214 // fhLambda1ForW0AndCellCuts [i][j][k] = 0;
215  fhLambda0ForW0AndCellCutsEta0[i][j][k] = 0;
216  }
217  }
218  for(Int_t j = 0; j < 5; j++)
219  {
220  fhLambda0ForW0MC[i][j] = 0;
221 // fhLambda1ForW0MC[i][j] = 0;
222  }
223  }
224 
225  // Cluster size
226  fhDeltaIEtaDeltaIPhiE0[0] = 0 ; fhDeltaIEtaDeltaIPhiE2[0] = 0; fhDeltaIEtaDeltaIPhiE6[0] = 0;
227  fhDeltaIEtaDeltaIPhiE0[1] = 0 ; fhDeltaIEtaDeltaIPhiE2[1] = 0; fhDeltaIEtaDeltaIPhiE6[1] = 0;
228  fhDeltaIA[0] = 0 ; fhDeltaIAL0[0] = 0; fhDeltaIAL1[0] = 0;
229  fhDeltaIA[1] = 0 ; fhDeltaIAL0[1] = 0; fhDeltaIAL1[1] = 0;
230  fhDeltaIANCells[0] = 0 ; fhDeltaIANCells[1] = 0;
231  fhDeltaIAMC[0] = 0 ; fhDeltaIAMC[1] = 0;
232  fhDeltaIAMC[2] = 0 ; fhDeltaIAMC[3] = 0;
233 
234  // Exotic
235  for (Int_t ie = 0; ie < 10 ; ie++)
236  {
237  fhExoDTime[ie] = 0;
238  for (Int_t idt = 0; idt < 5 ; idt++)
239  {
240  fhExoNCell [ie][idt] = 0;
241  fhExoL0 [ie][idt] = 0;
242  fhExoL1 [ie][idt] = 0;
243  fhExoECross [ie][idt] = 0;
244  fhExoTime [ie][idt] = 0;
245  fhExoL0NCell [ie][idt] = 0;
246  fhExoL1NCell [ie][idt] = 0;
247  }
248  }
249 
250  // MC
251 
252  for(Int_t i = 0; i < 7; i++)
253  {
254  fhRecoMCE[i][0] = 0; fhRecoMCE[i][1] = 0;
255  fhRecoMCPhi[i][0] = 0; fhRecoMCPhi[i][1] = 0;
256  fhRecoMCEta[i][0] = 0; fhRecoMCEta[i][1] = 0;
257  fhRecoMCDeltaE[i][0] = 0; fhRecoMCDeltaE[i][1] = 0;
258  fhRecoMCRatioE[i][0] = 0; fhRecoMCRatioE[i][1] = 0;
259  fhRecoMCDeltaPhi[i][0] = 0; fhRecoMCDeltaPhi[i][1] = 0;
260  fhRecoMCDeltaEta[i][0] = 0; fhRecoMCDeltaEta[i][1] = 0;
261  }
262 
263  for(Int_t i = 0; i < 4; i++)
264  {
265  fhGenMCE[i] = 0;
266  fhGenMCPt[i] = 0;
267  fhGenMCEtaPhi[i] = 0;
268  fhGenMCAccE[i] = 0;
269  fhGenMCAccPt[i] = 0;
270  fhGenMCAccEtaPhi[i] = 0;
271  }
272 
273  for(Int_t i = 0; i < 14; i++)
274  {
275  fhEBinClusterEtaPhi[i] = 0 ;
276  fhEBinClusterColRow[i] = 0 ;
277  fhEBinCellColRow [i] = 0 ;
278  }
279 
280  for(Int_t bc = 0; bc < 4; bc++) fhTimePerSMPerBC[bc] = 0 ;
281 
282  InitParameters();
283 }
284 
285 //______________________________________________________________________________________________________________________
295 //______________________________________________________________________________________________________________________
296 void AliAnaCalorimeterQA::BadClusterHistograms(AliVCluster* clus, const TObjArray *caloClusters, AliVCaloCells * cells,
297  Int_t absIdMax, Double_t maxCellFraction, Float_t eCrossFrac,
298  Double_t tmax)
299 {
300  // printf("AliAnaCalorimeterQA::BadClusterHistograms() - Event %d - Calorimeter %s \n \t E %f, n cells %d, max cell absId %d, maxCellFrac %f\n",
301  // GetReader()->GetEventNumber(), GetCalorimeterString().Data(),
302  // clus->E(),clus->GetNCells(),absIdMax,maxCellFraction);
303 
304  Double_t tof = clus->GetTOF()*1.e9;
305  if(tof > 400) tof-=fConstantTimeShift;
306 
307  Float_t energy = clus->E();
308 
309  fhBadClusterEnergy ->Fill(energy, GetEventWeight());
310  fhBadClusterTimeEnergy ->Fill(energy, tof , GetEventWeight());
311  fhBadClusterMaxCellDiff ->Fill(energy, maxCellFraction, GetEventWeight());
312  fhBadClusterMaxCellECross->Fill(energy, eCrossFrac , GetEventWeight());
313 
314  Float_t phi = fClusterMomentum.Phi();
315  if(phi < 0) phi += TMath::TwoPi();
316 
317  if(energy > 0.5) fhBadClusterEtaPhi->Fill(fClusterMomentum.Eta(), phi, GetEventWeight());
318 
319  fhBadClusterLambda0->Fill(energy, clus->GetM02());
320  fhBadClusterLambda1->Fill(energy, clus->GetM20());
321 
322  if(fStudyClustersAsymmetry) ClusterAsymmetryHistograms(clus,absIdMax,kFALSE);
323 
324  // Clusters in event time difference bad minus good
325 
326  for(Int_t iclus2 = 0; iclus2 < caloClusters->GetEntriesFast(); iclus2++ )
327  {
328  AliVCluster* clus2 = (AliVCluster*) caloClusters->At(iclus2);
329 
330  if(clus->GetID() == clus2->GetID()) continue;
331 
332  Float_t maxCellFraction2 = 0.;
333  Int_t absIdMax2 = GetCaloUtils()->GetMaxEnergyCell(cells, clus2,maxCellFraction2);
334 
335  if(IsGoodCluster(absIdMax2, clus->GetM02(), clus->GetNCells(), cells) && clus2->GetM02() > 0.1 )
336  {
337  Double_t tof2 = clus2->GetTOF()*1.e9;
338  if(tof2>400) tof2-=fConstantTimeShift;
339 
340  fhBadClusterPairDiffTimeE ->Fill(clus->E(), (tof-tof2), GetEventWeight());
341  }
342  } // loop
343 
344  // Max cell compared to other cells in cluster
345  if(fFillAllCellTimeHisto)
346  {
347  // Get some time averages
348  Double_t timeAverages[2] = {0.,0.};
349  CalculateAverageTime(clus, cells, timeAverages);
350 
351  fhBadClusterMaxCellDiffAverageTime ->Fill(clus->E(), tmax-timeAverages[0], GetEventWeight());
352  fhBadClusterMaxCellDiffWeightedTime->Fill(clus->E(), tmax-timeAverages[1], GetEventWeight());
353  }
354 
355  for (Int_t ipos = 0; ipos < clus->GetNCells(); ipos++)
356  {
357  Int_t absId = clus->GetCellsAbsId()[ipos];
358  if(absId!=absIdMax && cells->GetCellAmplitude(absIdMax) > 0.01)
359  {
360  Float_t frac = cells->GetCellAmplitude(absId)/cells->GetCellAmplitude(absIdMax);
361 
362  fhBadClusterMaxCellCloseCellRatio->Fill(clus->E(), frac, GetEventWeight());
363  fhBadClusterMaxCellCloseCellDiff ->Fill(clus->E(), cells->GetCellAmplitude(absIdMax)-cells->GetCellAmplitude(absId), GetEventWeight());
364 
365  if(fFillAllCellTimeHisto)
366  {
367  Double_t time = cells->GetCellTime(absId);
368  GetCaloUtils()->RecalibrateCellTime(time, GetCalorimeter(), absId,GetReader()->GetInputEvent()->GetBunchCrossNumber());
369 
370  Float_t diff = (tmax-(time*1e9-fConstantTimeShift));
371  fhBadCellTimeSpreadRespectToCellMax->Fill(clus->E(), diff, GetEventWeight());
372 
373  }
374  } // Not max
375  } // loop
376 }
377 
378 //______________________________________________________________________
383 //______________________________________________________________________
384 void AliAnaCalorimeterQA::CalculateAverageTime(AliVCluster *clus,
385  AliVCaloCells* cells,
386  Double_t timeAverages[2])
387 {
388  // First recalculate energy in case non linearity was applied
389  Float_t energy = 0;
390  Float_t ampMax = 0, amp = 0;
391 //Int_t absIdMax = -1;
392 
393  for (Int_t ipos = 0; ipos < clus->GetNCells(); ipos++)
394  {
395  Int_t id = clus->GetCellsAbsId()[ipos];
396 
397  // Recalibrate cell energy if needed
398  amp = cells->GetCellAmplitude(id);
399  GetCaloUtils()->RecalibrateCellAmplitude(amp,GetCalorimeter(), id);
400 
401  energy += amp;
402 
403  if(amp> ampMax)
404  {
405  ampMax = amp;
406 // absIdMax = id;
407  }
408  } // energy loop
409 
410  // Calculate average time of cells in cluster and weighted average
411  Double_t aTime = 0;
412  Double_t wTime = 0;
413  Float_t wTot = 0;
414  Double_t time = 0;
415  Int_t id =-1;
416  Double_t w = 0;
417  Int_t ncells = clus->GetNCells();
418 
419  for (Int_t ipos = 0; ipos < ncells; ipos++)
420  {
421  id = clus ->GetCellsAbsId()[ipos];
422  amp = cells->GetCellAmplitude(id);
423  time = cells->GetCellTime(id);
424 
425  // Recalibrate energy and time
426  GetCaloUtils()->RecalibrateCellAmplitude(amp , GetCalorimeter(), id);
427  GetCaloUtils()->RecalibrateCellTime (time, GetCalorimeter(), id, GetReader()->GetInputEvent()->GetBunchCrossNumber());
428 
429  w = GetCaloUtils()->GetEMCALRecoUtils()->GetCellWeight(cells->GetCellAmplitude(id),energy);
430  aTime += time*1e9;
431  wTime += time*1e9 * w;
432  wTot += w;
433  }
434 
435  if(ncells > 0) aTime /= ncells;
436  else aTime = 0;
437 
438  if(wTot > 0) wTime /= wTot;
439  else wTime = 0;
440 
441  timeAverages[0] = aTime;
442  timeAverages[1] = wTime;
443 }
444 
445 //____________________________________________________________
447 //____________________________________________________________
448 void AliAnaCalorimeterQA::CellHistograms(AliVCaloCells *cells)
449 {
450  Int_t ncells = cells->GetNumberOfCells();
451  if( ncells > 0 ) fhNCells->Fill(ncells, GetEventWeight()) ; // Not ok for PHOS with CPV
452 
453  Int_t ncellsCut = 0;
454  Float_t ecellsCut = 0;
455 
456  AliDebug(1,Form("%s cell entries %d", GetCalorimeterString().Data(), ncells));
457 
458  // Init arrays and used variables
459  Int_t *nCellsInModule = new Int_t [fNModules];
460  Float_t *eCellsInModule = new Float_t[fNModules];
461 
462  for(Int_t imod = 0; imod < fNModules; imod++ )
463  {
464  nCellsInModule[imod] = 0 ;
465  eCellsInModule[imod] = 0.;
466  }
467 
468  Int_t icol = -1, icolAbs = -1;
469  Int_t irow = -1, irowAbs = -1;
470  Int_t iRCU = -1;
471  Float_t amp = 0.;
472  Double_t time = 0.;
473  Int_t id = -1;
474  Bool_t highG = kFALSE;
475  Float_t recalF = 1.;
476  Int_t bc = GetReader()->GetInputEvent()->GetBunchCrossNumber();
477 
478  for (Int_t iCell = 0; iCell < cells->GetNumberOfCells(); iCell++)
479  {
480  if ( cells->GetCellNumber(iCell) < 0 ) continue; // CPV
481 
482  AliDebug(2,Form("Cell : amp %f, absId %d", cells->GetAmplitude(iCell), cells->GetCellNumber(iCell)));
483 
484  Int_t nModule = GetModuleNumberCellIndexesAbsCaloMap(cells->GetCellNumber(iCell),GetCalorimeter(),
485  icol , irow, iRCU,
486  icolAbs, irowAbs );
487 
488  AliDebug(2,Form("\t module %d, column %d (%d), row %d (%d)", nModule,icolAbs,icol,irowAbs,irow));
489 
490  if(nModule < fNModules)
491  {
492  //Check if the cell is a bad channel
493  if(GetCaloUtils()->IsBadChannelsRemovalSwitchedOn())
494  {
495  if(GetCalorimeter()==kEMCAL)
496  {
497  if(GetCaloUtils()->GetEMCALChannelStatus(nModule,icol,irow)) continue;
498  }
499  else
500  {
501  if(GetCaloUtils()->GetPHOSChannelStatus(nModule,icol,irow) ) continue;
502  }
503  } // use bad channel map
504 
505  amp = cells->GetAmplitude(iCell)*recalF;
506  time = cells->GetTime(iCell);
507  id = cells->GetCellNumber(iCell);
508  highG = cells->GetCellHighGain(id);
509  if(IsDataMC()) highG = kTRUE; // MC does not distinguish High and Low, put them all in high
510 
511  // Amplitude recalibration if set
512  GetCaloUtils()->RecalibrateCellAmplitude(amp, GetCalorimeter(), id);
513 
514  // Time recalibration if set
515  GetCaloUtils()->RecalibrateCellTime (time, GetCalorimeter(), id, GetReader()->GetInputEvent()->GetBunchCrossNumber());
516 
517  // Transform time to ns
518  time *= 1.0e9;
519  time-=fConstantTimeShift;
520 
521  if(time < fTimeCutMin || time > fTimeCutMax)
522  {
523  AliDebug(1,Form("Remove cell with Time %f",time));
524  continue;
525  }
526 
527  // Remove exotic cells, defined only for EMCAL
528  if(GetCalorimeter()==kEMCAL &&
529  GetCaloUtils()->GetEMCALRecoUtils()->IsExoticCell(id, cells, bc)) continue;
530 
531  fhAmplitude->Fill(amp, GetEventWeight());
532  fhAmpId ->Fill(amp, id , GetEventWeight());
533  fhAmpMod ->Fill(amp, nModule, GetEventWeight());
534  fhAmpWeirdMod->Fill(amp, nModule, GetEventWeight());
535  if(!highG) fhAmpIdLowGain->Fill(amp, id, GetEventWeight());
536 
537  if(fFillEBinAcceptanceHisto)
538  {
539  for(Int_t ie = 0; ie < fNEBinCuts; ie++)
540  {
541  if( amp >= fEBinCuts[ie] && amp < fEBinCuts[ie+1] )
542  {
543  fhEBinCellColRow[ie]->Fill(icolAbs,irowAbs,GetEventWeight()) ;
544  }
545  }
546  }
547 
548  // E cross for exotic cells
549  if(amp > 0.05)
550  {
551  fhCellECross->Fill(amp, 1-GetECross(id,cells)/amp, GetEventWeight());
552  ecellsCut+=amp ;
553  if(fStudyWeight) eCellsInModule[nModule]+=amp ;
554  }
555 
556  if ( amp > fCellAmpMin )
557  {
558  ncellsCut++ ;
559  nCellsInModule[nModule]++ ;
560 
561  if(!fStudyWeight) eCellsInModule[nModule]+=amp ;
562 
563  fhGridCells ->Fill(icolAbs, irowAbs, GetEventWeight());
564  fhGridCellsE->Fill(icolAbs, irowAbs, amp );
565 
566  if(!highG)
567  {
568  fhGridCellsLowGain ->Fill(icolAbs, irowAbs, GetEventWeight());
569  fhGridCellsELowGain->Fill(icolAbs, irowAbs, amp );
570  }
571 
572  if(fFillAllCellTimeHisto)
573  {
574  //printf("%s: time %g\n",GetCalorimeterString().Data(), time);
575 
576 // Double_t v[3] = {0,0,0}; //vertex ;
577 // GetReader()->GetVertex(v);
578 // if(amp > 0.5) fhTimeVz ->Fill(TMath::Abs(v[2]), time, GetEventWeight());
579 
580  fhTime ->Fill(time, GetEventWeight());
581  fhTimeId ->Fill(time, id , GetEventWeight());
582  fhTimeAmp->Fill(amp , time, GetEventWeight());
583 
584  Int_t bc = (GetReader()->GetInputEvent()->GetBunchCrossNumber())%4;
585  fhTimePerSMPerBC[bc]->Fill(time, nModule, GetEventWeight());
586 
587  fhGridCellsTime->Fill(icolAbs, irowAbs, time);
588  if(!highG) fhGridCellsTimeLowGain->Fill(icolAbs, irowAbs, time);
589 
590  fhTimeMod->Fill(time, nModule, GetEventWeight());
591  fhTimeAmpPerRCU[nModule*fNRCU+iRCU]->Fill(amp, time, GetEventWeight());
592 
593  if(!highG)
594  {
595  fhTimeIdLowGain ->Fill(time, id , GetEventWeight());
596  fhTimeAmpLowGain->Fill(amp , time, GetEventWeight());
597  }
598  }
599  }
600 
601  // Get Eta-Phi position of Cell
602  if(fFillAllPosHisto)
603  {
604  if ( GetCalorimeter() == kEMCAL && GetCaloUtils()->IsEMCALGeoMatrixSet() )
605  {
606  Float_t celleta = 0.;
607  Float_t cellphi = 0.;
608  GetEMCALGeometry()->EtaPhiFromIndex(id, celleta, cellphi);
609 
610  if ( cellphi < 0 ) cellphi+=TMath::TwoPi();
611 
612  if(fFillAllTH3)
613  fhEtaPhiAmpCell->Fill(celleta, cellphi, amp, GetEventWeight());
614  else
615  fhEtaPhiCell ->Fill(celleta, cellphi, GetEventWeight());
616 
617  Double_t cellpos[] = {0, 0, 0};
618  GetEMCALGeometry()->GetGlobal(id, cellpos);
619 
620  fhXCellE->Fill(cellpos[0], amp, GetEventWeight()) ;
621  fhYCellE->Fill(cellpos[1], amp, GetEventWeight()) ;
622  fhZCellE->Fill(cellpos[2], amp, GetEventWeight()) ;
623 
624  Float_t rcell = TMath::Sqrt(cellpos[0]*cellpos[0]+cellpos[1]*cellpos[1]);//+cellpos[2]*cellpos[2]);
625  fhRCellE ->Fill(rcell, amp, GetEventWeight()) ;
626 
627  fhXYZCell->Fill(cellpos[0], cellpos[1], cellpos[2], GetEventWeight()) ;
628  } // EMCAL Cells
629  else if ( GetCalorimeter() == kPHOS && GetCaloUtils()->IsPHOSGeoMatrixSet() )
630  {
631  TVector3 xyz;
632  Int_t relId[4], module;
633  Float_t xCell, zCell;
634 
635  GetPHOSGeometry()->AbsToRelNumbering(id,relId);
636  module = relId[0];
637  GetPHOSGeometry()->RelPosInModule(relId,xCell,zCell);
638  GetPHOSGeometry()->Local2Global(module,xCell,zCell,xyz);
639 
640  Float_t rcell = TMath::Sqrt(xyz.X()*xyz.X()+xyz.Y()*xyz.Y());
641 
642  fhXCellE ->Fill(xyz.X(), amp, GetEventWeight()) ;
643  fhYCellE ->Fill(xyz.Y(), amp, GetEventWeight()) ;
644  fhZCellE ->Fill(xyz.Z(), amp, GetEventWeight()) ;
645  fhRCellE ->Fill(rcell , amp, GetEventWeight()) ;
646 
647  fhXYZCell->Fill(xyz.X(), xyz.Y(), xyz.Z(), GetEventWeight()) ;
648  } // PHOS cells
649  } // Fill cell position histograms
650 
651  } // N modules
652  } // Cell loop
653 
654  // Fill the cells after the cut on min amplitude and bad/exotic channels
655  if( ncellsCut > 0 ) fhNCellsCutAmpMin->Fill(ncellsCut, GetEventWeight()) ;
656 
657  // Number of cells per module
658  for(Int_t imod = 0; imod < fNModules; imod++ )
659  {
660  AliDebug(1,Form("Module %d, calo %s, N cells %d, sum Amp %f", imod, GetCalorimeterString().Data(), nCellsInModule[imod], eCellsInModule[imod]));
661 
662  fhNCellsMod ->Fill(nCellsInModule[imod], imod, GetEventWeight()) ;
663  fhSumCellsAmpMod->Fill(eCellsInModule[imod], imod, GetEventWeight()) ;
664 
665  fhNCellsSumAmpPerMod[imod]->Fill(eCellsInModule[imod], nCellsInModule[imod], GetEventWeight());
666  }
667 
668  // Check energy distribution in calorimeter for selected cells
669  if(fStudyWeight)
670  {
671  for (Int_t iCell = 0; iCell < cells->GetNumberOfCells(); iCell++)
672  {
673  if ( cells->GetCellNumber(iCell) < 0 ) continue; // CPV
674 
675  AliDebug(2,Form("Cell : amp %f, absId %d", cells->GetAmplitude(iCell), cells->GetCellNumber(iCell)));
676 
677  Int_t nModule = GetModuleNumberCellIndexes(cells->GetCellNumber(iCell),GetCalorimeter(), icol, irow, iRCU);
678 
679  AliDebug(2,Form("\t module %d, column %d, row %d", nModule,icol,irow));
680 
681  if(nModule < fNModules)
682  {
683  //Check if the cell is a bad channel
684  if(GetCaloUtils()->IsBadChannelsRemovalSwitchedOn())
685  {
686  if(GetCalorimeter()==kEMCAL)
687  {
688  if(GetCaloUtils()->GetEMCALChannelStatus(nModule, icol, irow)) continue;
689  }
690  else
691  {
692  if(GetCaloUtils()->GetPHOSChannelStatus (nModule, icol, irow) ) continue;
693  }
694  } // use bad channel map
695 
696  amp = cells->GetAmplitude(iCell)*recalF;
697  time = cells->GetTime(iCell);
698  id = cells->GetCellNumber(iCell);
699 
700  // Amplitude recalibration if set
701  GetCaloUtils()->RecalibrateCellAmplitude(amp, GetCalorimeter(), id);
702 
703  // Time recalibration if set
704  GetCaloUtils()->RecalibrateCellTime (time, GetCalorimeter(), id,
705  GetReader()->GetInputEvent()->GetBunchCrossNumber());
706 
707  // Transform time to ns
708  time *= 1.0e9;
709  time -= fConstantTimeShift;
710 
711  if(time < fTimeCutMin || time > fTimeCutMax)
712  {
713  AliDebug(1,Form("Remove cell with Time %f",time));
714  continue;
715  }
716 
717  // Remove exotic cells, defined only for EMCAL
718  if(GetCalorimeter()==kEMCAL &&
719  GetCaloUtils()->GetEMCALRecoUtils()->IsExoticCell(id, cells, bc)) continue;
720 
721  // E cross for exotic cells
722  if(amp > 0.05)
723  {
724  if(ecellsCut > 0)
725  {
726  Float_t ratio = amp/ecellsCut;
727  fhECellTotalRatio ->Fill(ecellsCut, ratio , GetEventWeight());
728  fhECellTotalLogRatio ->Fill(ecellsCut,TMath::Log(ratio), GetEventWeight());
729  }
730 
731  if(eCellsInModule[nModule] > 0)
732  {
733  Float_t ratioMod = amp/eCellsInModule[nModule];
734  fhECellTotalRatioMod [nModule]->Fill(eCellsInModule[nModule], ratioMod , GetEventWeight());
735  fhECellTotalLogRatioMod[nModule]->Fill(eCellsInModule[nModule],TMath::Log(ratioMod), GetEventWeight());
736  }
737  } // amp > 0.5
738  } // nMod > 0 < Max
739  } // cell loop
740  } // weight studies
741 
742  delete [] nCellsInModule;
743  delete [] eCellsInModule;
744 }
745 
746 //__________________________________________________________________________
748 //__________________________________________________________________________
749 void AliAnaCalorimeterQA::CellInClusterPositionHistograms(AliVCluster* clus)
750 {
751  Int_t nCaloCellsPerCluster = clus->GetNCells();
752 
753  UShort_t * indexList = clus->GetCellsAbsId();
754 
755  Float_t pos[3];
756  clus->GetPosition(pos);
757 
758  Float_t clEnergy = clus->E();
759 
760  // Loop on cluster cells
761  for (Int_t ipos = 0; ipos < nCaloCellsPerCluster; ipos++)
762  {
763  // printf("Index %d\n",ipos);
764  Int_t absId = indexList[ipos];
765 
766  //Get position of cell compare to cluster
767 
768  if ( GetCalorimeter() == kEMCAL && GetCaloUtils()->IsEMCALGeoMatrixSet() )
769  {
770  Double_t cellpos[] = {0, 0, 0};
771  GetEMCALGeometry()->GetGlobal(absId, cellpos);
772 
773  fhDeltaCellClusterXNCells->Fill(pos[0]-cellpos[0], nCaloCellsPerCluster, GetEventWeight()) ;
774  fhDeltaCellClusterYNCells->Fill(pos[1]-cellpos[1], nCaloCellsPerCluster, GetEventWeight()) ;
775  fhDeltaCellClusterZNCells->Fill(pos[2]-cellpos[2], nCaloCellsPerCluster, GetEventWeight()) ;
776 
777  fhDeltaCellClusterXE->Fill(pos[0]-cellpos[0], clEnergy, GetEventWeight()) ;
778  fhDeltaCellClusterYE->Fill(pos[1]-cellpos[1], clEnergy, GetEventWeight()) ;
779  fhDeltaCellClusterZE->Fill(pos[2]-cellpos[2], clEnergy, GetEventWeight()) ;
780 
781  Float_t r = TMath::Sqrt(pos[0] *pos[0] + pos[1] * pos[1] );
782  Float_t rcell = TMath::Sqrt(cellpos[0]*cellpos[0] + cellpos[1]* cellpos[1]);
783 
784  fhDeltaCellClusterRNCells->Fill(r-rcell, nCaloCellsPerCluster, GetEventWeight()) ;
785  fhDeltaCellClusterRE ->Fill(r-rcell, clEnergy , GetEventWeight()) ;
786  } // EMCAL and its matrices are available
787  else if ( GetCalorimeter() == kPHOS && GetCaloUtils()->IsPHOSGeoMatrixSet() )
788  {
789  TVector3 xyz;
790  Int_t relId[4], module;
791  Float_t xCell, zCell;
792 
793  GetPHOSGeometry()->AbsToRelNumbering(absId,relId);
794  module = relId[0];
795  GetPHOSGeometry()->RelPosInModule(relId,xCell,zCell);
796  GetPHOSGeometry()->Local2Global(module,xCell,zCell,xyz);
797 
798  fhDeltaCellClusterXNCells->Fill(pos[0]-xyz.X(), nCaloCellsPerCluster, GetEventWeight()) ;
799  fhDeltaCellClusterYNCells->Fill(pos[1]-xyz.Y(), nCaloCellsPerCluster, GetEventWeight()) ;
800  fhDeltaCellClusterZNCells->Fill(pos[2]-xyz.Z(), nCaloCellsPerCluster, GetEventWeight()) ;
801 
802  fhDeltaCellClusterXE->Fill(pos[0]-xyz.X(), clEnergy, GetEventWeight()) ;
803  fhDeltaCellClusterYE->Fill(pos[1]-xyz.Y(), clEnergy, GetEventWeight()) ;
804  fhDeltaCellClusterZE->Fill(pos[2]-xyz.Z(), clEnergy, GetEventWeight()) ;
805 
806  Float_t r = TMath::Sqrt(pos[0] * pos[0] + pos[1] * pos[1] );
807  Float_t rcell = TMath::Sqrt(xyz.X() * xyz.X() + xyz.Y() * xyz.Y());
808 
809  fhDeltaCellClusterRNCells->Fill(r-rcell, nCaloCellsPerCluster, GetEventWeight()) ;
810  fhDeltaCellClusterRE ->Fill(r-rcell, clEnergy , GetEventWeight()) ;
811  } // PHOS and its matrices are available
812  } // cluster cell loop
813 }
814 
815 //_____________________________________________________________________________________
816 // Study the shape of the cluster in cell units terms.
817 //_____________________________________________________________________________________
818 void AliAnaCalorimeterQA::ClusterAsymmetryHistograms(AliVCluster* clus, Int_t absIdMax,
819  Bool_t goodCluster)
820 {
821  // No use to study clusters with less than 4 cells
822  if( clus->GetNCells() <= 3 ) return;
823 
824  Int_t dIeta = 0;
825  Int_t dIphi = 0;
826 
827  Int_t ietaMax=-1; Int_t iphiMax = 0; Int_t rcuMax = 0;
828  Int_t smMax = GetModuleNumberCellIndexes(absIdMax,GetCalorimeter(), ietaMax, iphiMax, rcuMax);
829 
830  for (Int_t ipos = 0; ipos < clus->GetNCells(); ipos++)
831  {
832  Int_t absId = clus->GetCellsAbsId()[ipos];
833 
834  Int_t ieta=-1; Int_t iphi = 0; Int_t rcu = 0;
835  Int_t sm = GetModuleNumberCellIndexes(absId,GetCalorimeter(), ieta, iphi, rcu);
836 
837  if(dIphi < TMath::Abs(iphi-iphiMax)) dIphi = TMath::Abs(iphi-iphiMax);
838 
839  if(smMax==sm)
840  {
841  if(dIeta < TMath::Abs(ieta-ietaMax)) dIeta = TMath::Abs(ieta-ietaMax);
842  }
843  else
844  {
845  Int_t ietaShift = ieta;
846  Int_t ietaMaxShift = ietaMax;
847  if (ieta > ietaMax) ietaMaxShift+=48;
848  else ietaShift +=48;
849  if(dIeta < TMath::Abs(ietaShift-ietaMaxShift)) dIeta = TMath::Abs(ietaShift-ietaMaxShift);
850  }
851  }// Fill cell-cluster histogram loop
852 
853  Float_t dIA = 1.*(dIphi-dIeta)/(dIeta+dIphi);
854 
855  if(goodCluster)
856  {
857  // Was cluster matched?
858  Bool_t matched = GetCaloPID()->IsTrackMatched(clus,GetCaloUtils(),GetReader()->GetInputEvent());
859 
860  if (clus->E() < 2 ) fhDeltaIEtaDeltaIPhiE0[matched]->Fill(dIeta, dIphi, GetEventWeight());
861  else if(clus->E() < 6 ) fhDeltaIEtaDeltaIPhiE2[matched]->Fill(dIeta, dIphi, GetEventWeight());
862  else fhDeltaIEtaDeltaIPhiE6[matched]->Fill(dIeta, dIphi, GetEventWeight());
863 
864  fhDeltaIA[matched]->Fill(clus->E(), dIA, GetEventWeight());
865 
866  if(clus->E() > 0.5)
867  {
868  fhDeltaIAL0 [matched]->Fill(clus->GetM02() , dIA, GetEventWeight());
869  fhDeltaIAL1 [matched]->Fill(clus->GetM20() , dIA, GetEventWeight());
870  fhDeltaIANCells[matched]->Fill(clus->GetNCells(), dIA, GetEventWeight());
871  }
872 
873  // Origin of clusters
874  Int_t nLabel = clus->GetNLabels();
875  Int_t* labels = clus->GetLabels();
876 
877  if(IsDataMC())
878  {
879  Int_t tag = GetMCAnalysisUtils()->CheckOrigin(labels,nLabel, GetReader(),GetCalorimeter());
880  if( GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCPhoton) &&
881  !GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCPi0) &&
882  !GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCEta) &&
883  !GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCConversion) ){
884  fhDeltaIAMC[0]->Fill(clus->E(), dIA, GetEventWeight()); // Pure Photon
885  }
886  else if ( GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCElectron) &&
887  !GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCConversion) ){
888  fhDeltaIAMC[1]->Fill(clus->E(), dIA, GetEventWeight()); // Pure electron
889  }
890  else if ( GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCPhoton) &&
891  GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCConversion) ){
892  fhDeltaIAMC[2]->Fill(clus->E(), dIA, GetEventWeight()); // Converted cluster
893  }
894  else if(!GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCPhoton)){
895  fhDeltaIAMC[3]->Fill(clus->E(), dIA, GetEventWeight()); // Hadrons
896  }
897 
898  } // MC
899  } // good cluster
900  else
901  {
902  if (clus->E() < 2 ) fhBadClusterDeltaIEtaDeltaIPhiE0->Fill(dIeta, dIphi, GetEventWeight());
903  else if(clus->E() < 6 ) fhBadClusterDeltaIEtaDeltaIPhiE2->Fill(dIeta, dIphi, GetEventWeight());
904  else fhBadClusterDeltaIEtaDeltaIPhiE6->Fill(dIeta, dIphi, GetEventWeight());
905 
906  fhBadClusterDeltaIA->Fill(clus->E(), dIA, GetEventWeight());
907  }
908 }
909 
910 //__________________________________________________________________________________________________________________
919 //__________________________________________________________________________________________________________________
920 void AliAnaCalorimeterQA::ClusterHistograms(AliVCluster* clus, const TObjArray *caloClusters, AliVCaloCells * cells,
921  Int_t absIdMax, Double_t maxCellFraction, Float_t eCrossFrac,
922  Double_t tmax)
923 {
924  Double_t tof = clus->GetTOF()*1.e9;
925  if(tof>400) tof-=fConstantTimeShift;
926 
927  fhLambda0 ->Fill(clus->E(), clus->GetM02() , GetEventWeight());
928  fhLambda1 ->Fill(clus->E(), clus->GetM20() , GetEventWeight());
929  // fhDispersion ->Fill(clus->E(), clus->GetDispersion(), GetEventWeight());
930 
931  fhClusterMaxCellDiff ->Fill(clus->E(), maxCellFraction, GetEventWeight());
932  fhClusterMaxCellECross->Fill(clus->E(), eCrossFrac , GetEventWeight());
933  fhClusterTimeEnergy ->Fill(clus->E(), tof , GetEventWeight());
934 
935  if(fStudyClustersAsymmetry) ClusterAsymmetryHistograms(clus,absIdMax,kTRUE);
936 
937  Int_t nModule = GetModuleNumber(clus);
938  Int_t nCaloCellsPerCluster = clus->GetNCells();
939 
940  // Clusters in event time difference
941  for(Int_t iclus2 = 0; iclus2 < caloClusters->GetEntriesFast(); iclus2++ )
942  {
943  AliVCluster* clus2 = (AliVCluster*) caloClusters->At(iclus2);
944 
945  if( clus->GetID() == clus2->GetID() ) continue;
946 
947  if( clus->GetM02() > 0.01 && clus2->GetM02() > 0.01 )
948  {
949  Int_t nModule2 = GetModuleNumber(clus2);
950 
951  Double_t tof2 = clus2->GetTOF()*1.e9;
952  if(tof2>400) tof2-=fConstantTimeShift;
953 
954  fhClusterPairDiffTimeE ->Fill(clus->E(), tof-tof2, GetEventWeight());
955 
956  if ( nModule2 == nModule )
957  fhClusterPairDiffTimeESameMod->Fill(clus->E(), tof-tof2, GetEventWeight());
958  }
959  }
960 
961  if(nCaloCellsPerCluster > 1)
962  {
963  // Check time of cells respect to max energy cell
964 
965 // if(fFillAllCellTimeHisto)
966 // {
967 // // Get some time averages
968 // Double_t timeAverages[2] = {0.,0.};
969 // CalculateAverageTime(clus, cells, timeAverages);
970 //
971 // fhClusterMaxCellDiffAverageTime ->Fill(clus->E(), tmax-timeAverages[0], GetEventWeight());
972 // fhClusterMaxCellDiffWeightedTime ->Fill(clus->E(), tmax-timeAverages[1], GetEventWeight());
973 // }
974 
975  for (Int_t ipos = 0; ipos < nCaloCellsPerCluster; ipos++)
976  {
977  Int_t absId = clus->GetCellsAbsId()[ipos];
978  if( absId == absIdMax || cells->GetCellAmplitude(absIdMax) < 0.01 ) continue;
979 
980  Float_t frac = cells->GetCellAmplitude(absId)/cells->GetCellAmplitude(absIdMax);
981  fhClusterMaxCellCloseCellRatio->Fill(clus->E(), frac, GetEventWeight());
982  fhClusterMaxCellCloseCellDiff ->Fill(clus->E(), cells->GetCellAmplitude(absIdMax)-cells->GetCellAmplitude(absId), GetEventWeight());
983 
984  if(fFillAllCellTimeHisto)
985  {
986  Double_t time = cells->GetCellTime(absId);
987  GetCaloUtils()->RecalibrateCellTime(time, GetCalorimeter(), absId,GetReader()->GetInputEvent()->GetBunchCrossNumber());
988 
989  Float_t diff = (tmax-(time*1.0e9-fConstantTimeShift));
990  fhCellTimeSpreadRespectToCellMax->Fill(clus->E(), diff, GetEventWeight());
991  if(TMath::Abs(TMath::Abs(diff) > 100) && clus->E() > 1 ) fhCellIdCellLargeTimeSpread->Fill(absId, GetEventWeight());
992  }
993 
994  } // Fill cell-cluster histogram loop
995 
996  } // Check time and energy of cells respect to max energy cell if cluster of more than 1 cell
997 
998  Float_t e = fClusterMomentum.E();
999  Float_t pt = fClusterMomentum.Pt();
1000  Float_t eta = fClusterMomentum.Eta();
1001  Float_t phi = GetPhi(fClusterMomentum.Phi());
1002 
1003  AliDebug(1,Form("cluster: E %2.3f, pT %2.3f, eta %2.3f, phi %2.3f",e,pt,eta,phi*TMath::RadToDeg()));
1004 
1005  fhE ->Fill(e, GetEventWeight());
1006  if(nModule >=0 && nModule < fNModules)
1007  {
1008  fhEMod ->Fill(e, nModule, GetEventWeight());
1009  fhEWeirdMod->Fill(e, nModule, GetEventWeight()); // different binning
1010  }
1011 
1012  fhPt ->Fill(pt , GetEventWeight());
1013  fhPhi ->Fill(phi, GetEventWeight());
1014  fhEta ->Fill(eta, GetEventWeight());
1015 
1016  if(fFillEBinAcceptanceHisto)
1017  {
1018  Float_t maxCellFraction = 0;
1019  Int_t absIdMax = GetCaloUtils()->GetMaxEnergyCell(cells,clus,maxCellFraction);
1020 
1021  Int_t icol = -1, irow = -1, iRCU = -1, icolAbs = -1, irowAbs = -1;
1022  GetModuleNumberCellIndexesAbsCaloMap(absIdMax,GetCalorimeter(), icol, irow, iRCU, icolAbs, irowAbs);
1023 
1024  for(Int_t ie = 0; ie < fNEBinCuts; ie++)
1025  {
1026  if( e >= fEBinCuts[ie] && e < fEBinCuts[ie+1] )
1027  {
1028  fhEBinClusterEtaPhi[ie]->Fill(eta,phi,GetEventWeight()) ;
1029 
1030  fhEBinClusterColRow[ie]->Fill(icolAbs,irowAbs,GetEventWeight()) ;
1031  }
1032  }
1033  }
1034  else if ( fFillAllTH3 ) fhEtaPhiE->Fill(eta, phi, e, GetEventWeight());
1035  else if ( e > 0.5 ) fhEtaPhi ->Fill(eta, phi, GetEventWeight());
1036 
1037  // Cells per cluster
1038  fhNCellsPerCluster->Fill(e, nCaloCellsPerCluster, GetEventWeight());
1039 
1040  if(e > 100) fhNCellsPerClusterWeirdMod->Fill(nCaloCellsPerCluster, nModule, GetEventWeight());
1041 
1042  // Position
1043  if(fFillAllPosHisto2)
1044  {
1045  Float_t pos[3] ;
1046  clus->GetPosition(pos);
1047 
1048  fhXE ->Fill(pos[0], e, GetEventWeight());
1049  fhYE ->Fill(pos[1], e, GetEventWeight());
1050  fhZE ->Fill(pos[2], e, GetEventWeight());
1051 
1052  fhXYZ ->Fill(pos[0], pos[1], pos[2], GetEventWeight());
1053 
1054  fhXNCells->Fill(pos[0], nCaloCellsPerCluster, GetEventWeight());
1055  fhYNCells->Fill(pos[1], nCaloCellsPerCluster, GetEventWeight());
1056  fhZNCells->Fill(pos[2], nCaloCellsPerCluster, GetEventWeight());
1057 
1058  Float_t rxyz = TMath::Sqrt(pos[0]*pos[0]+pos[1]*pos[1]);//+pos[2]*pos[2]);
1059 
1060  fhRE ->Fill(rxyz, e , GetEventWeight());
1061  fhRNCells->Fill(rxyz, nCaloCellsPerCluster, GetEventWeight());
1062  }
1063 
1064  if( nModule >= 0 && nModule < fNModules ) fhNCellsPerClusterMod[nModule]->Fill(e, nCaloCellsPerCluster, GetEventWeight());
1065 }
1066 
1067 //____________________________________________________________________________
1079 //____________________________________________________________________________
1080 void AliAnaCalorimeterQA::ClusterLoopHistograms(const TObjArray *caloClusters,
1081  AliVCaloCells* cells)
1082 {
1083  Int_t nLabel = 0 ;
1084  Int_t *labels = 0x0;
1085  Int_t nCaloClusters = caloClusters->GetEntriesFast() ;
1086  Int_t nCaloClustersAccepted = 0 ;
1087  Int_t nCaloCellsPerCluster = 0 ;
1088  Bool_t matched = kFALSE;
1089  Int_t nModule =-1 ;
1090 
1091  // Get vertex for photon momentum calculation and event selection
1092  Double_t v[3] = {0,0,0}; //vertex ;
1093 //GetReader()->GetVertex(v);
1094 
1095  Int_t *nClustersInModule = new Int_t [fNModules];
1096  Float_t *energyInModule = new Float_t[fNModules];
1097  for(Int_t imod = 0; imod < fNModules; imod++ )
1098  {
1099  nClustersInModule[imod] = 0;
1100  energyInModule [imod] = 0;
1101  }
1102 
1103  AliDebug(1,Form("In %s there are %d clusters", GetCalorimeterString().Data(), nCaloClusters));
1104 
1105  // Loop over CaloClusters
1106  for(Int_t iclus = 0; iclus < nCaloClusters; iclus++)
1107  {
1108  AliDebug(1,Form("Cluster: %d/%d, data %d",iclus+1,nCaloClusters,GetReader()->GetDataType()));
1109 
1110  AliVCluster* clus = (AliVCluster*) caloClusters->At(iclus);
1111 
1112  // Get the fraction of the cluster energy that carries the cell with highest energy and its absId
1113  Float_t maxCellFraction = 0.;
1114  Int_t absIdMax = GetCaloUtils()->GetMaxEnergyCell(cells, clus,maxCellFraction);
1115 
1116  // Cut on time of clusters
1117  Double_t tof = clus->GetTOF()*1.e9;
1118  if(tof>400) tof-=fConstantTimeShift;
1119 
1120  if( tof < fTimeCutMin || tof > fTimeCutMax )
1121  {
1122  AliDebug(1,Form("Remove cluster with TOF %f",tof));
1123  continue;
1124  }
1125 
1126  // Get cluster kinematics
1127  clus->GetMomentum(fClusterMomentum,v);
1128 
1129  // Check only certain regions
1130  Bool_t in = kTRUE;
1131  if(IsFiducialCutOn()) in = GetFiducialCut()->IsInFiducialCut(fClusterMomentum.Eta(),fClusterMomentum.Phi(),GetCalorimeter()) ;
1132  if(!in) continue;
1133 
1134  // MC labels
1135  nLabel = clus->GetNLabels();
1136  labels = clus->GetLabels();
1137 
1138  // SuperModule number of cluster
1139  nModule = GetModuleNumber(clus);
1140 
1141  // Cells per cluster
1142  nCaloCellsPerCluster = clus->GetNCells();
1143 
1144  // Cluster mathed with track?
1145  matched = GetCaloPID()->IsTrackMatched(clus,GetCaloUtils(), GetReader()->GetInputEvent());
1146 
1147  // Get time of max cell
1148  Double_t tmax = cells->GetCellTime(absIdMax);
1149  GetCaloUtils()->RecalibrateCellTime(tmax, GetCalorimeter(), absIdMax,GetReader()->GetInputEvent()->GetBunchCrossNumber());
1150  tmax*=1.e9;
1151  tmax-=fConstantTimeShift;
1152  //
1153  // Fill histograms related to single cluster
1154  //
1155 
1156  // Fill some histograms before applying the exotic cell / bad map cut
1157  if(fStudyBadClusters)
1158  {
1159  fhNCellsPerClusterNoCut->Fill(clus->E(), nCaloCellsPerCluster, GetEventWeight());
1160 
1161  if(nModule >=0 && nModule < fNModules)
1162  {
1163  fhNCellsPerClusterModNoCut[nModule]->Fill(clus->E(), nCaloCellsPerCluster, GetEventWeight());
1164  if(clus->E() > 100) fhNCellsPerClusterWeirdModNoCut->Fill(nCaloCellsPerCluster, nModule, GetEventWeight());
1165  }
1166 
1167  fhClusterMaxCellDiffNoCut->Fill(clus->E(), maxCellFraction, GetEventWeight());
1168  }
1169 
1170  Float_t ampMax = cells->GetCellAmplitude(absIdMax);
1171  GetCaloUtils()->RecalibrateCellAmplitude(ampMax,GetCalorimeter(), absIdMax);
1172 
1173  if(fStudyExotic) ExoticHistograms(absIdMax, ampMax, clus, cells);
1174 
1175  // Check bad clusters if requested and rejection was not on
1176  Bool_t goodCluster = IsGoodCluster(absIdMax, clus->GetM02(), nCaloCellsPerCluster, cells);
1177 
1178  Float_t eCrossFrac = 0;
1179  if(ampMax > 0.01) eCrossFrac = 1-GetECross(absIdMax,cells)/ampMax;
1180 
1181  AliDebug(1,Form("Accept cluster? %d",goodCluster));
1182 
1183  if(!goodCluster)
1184  {
1185  if ( fStudyBadClusters ) BadClusterHistograms(clus, caloClusters,
1186  cells, absIdMax,
1187  maxCellFraction,
1188  eCrossFrac, tmax);
1189  continue;
1190  }
1191 
1192 
1193  ClusterHistograms(clus, caloClusters, cells, absIdMax,
1194  maxCellFraction, eCrossFrac, tmax);
1195 
1196  nCaloClustersAccepted++;
1197  nModule = GetModuleNumber(clus);
1198  if(nModule >=0 && nModule < fNModules && fClusterMomentum.E() > 2*fCellAmpMin)
1199  {
1200  nClustersInModule[nModule]++;
1201  if(clus->E() > 0.5)
1202  energyInModule [nModule] += clus->E();
1203  }
1204 
1205  // Cluster weights
1206  if(fStudyWeight) WeightHistograms(clus, cells);
1207 
1208  // Cells in cluster position
1209  if(fFillAllPosHisto) CellInClusterPositionHistograms(clus);
1210 
1211  // Fill histograms related to single cluster, mc vs data
1212  Int_t mcOK = kFALSE;
1213  Int_t pdg = -1;
1214  if(IsDataMC() && nLabel > 0 && labels)
1215  mcOK = ClusterMCHistograms(matched, labels, nLabel, pdg);
1216 
1217  // Matched clusters with tracks, also do some MC comparison, needs input from ClusterMCHistograms
1218  if( matched && fFillAllTMHisto)
1219  ClusterMatchedWithTrackHistograms(clus,mcOK,pdg);
1220 
1221  // Invariant mass
1222  // Try to reduce background with a mild shower shape cut and no more than 1 maxima
1223  // in cluster and remove low energy clusters
1224 
1225  if ( fFillAllPi0Histo
1226  && nCaloClusters > 1
1227  && nCaloCellsPerCluster > 1
1228  && GetCaloUtils()->GetNumberOfLocalMaxima(clus,cells) == 1
1229  && clus->GetM02() < fInvMassMaxM02Cut
1230  && clus->GetM02() > fInvMassMinM02Cut
1231  && clus->E() > fInvMassMinECut
1232  && clus->E() < fInvMassMaxECut
1233  )
1234  InvariantMassHistograms(iclus, nModule, caloClusters,cells);
1235 
1236  } // Cluster loop
1237 
1238  // Number of clusters histograms
1239  if(nCaloClustersAccepted > 0) fhNClusters->Fill(nCaloClustersAccepted, GetEventWeight());
1240 
1241  // Number of clusters per module
1242  for(Int_t imod = 0; imod < fNModules; imod++ )
1243  {
1244  AliDebug(1,Form("Module %d calo %s clusters %d, sum E %f", imod, GetCalorimeterString().Data(), nClustersInModule[imod], energyInModule[imod]));
1245 
1246  fhNClustersMod ->Fill(nClustersInModule[imod], imod, GetEventWeight());
1247  fhSumClustersEnergyMod->Fill(energyInModule [imod], imod, GetEventWeight());
1248 
1249  fhNClustersSumEnergyPerMod[imod]->Fill(energyInModule[imod], nClustersInModule[imod], GetEventWeight());
1250  }
1251 
1252  delete [] nClustersInModule;
1253  delete [] energyInModule;
1254 }
1255 
1256 //__________________________________________________________________________________
1259 //__________________________________________________________________________________
1260 Bool_t AliAnaCalorimeterQA::ClusterMCHistograms(Bool_t matched,const Int_t * labels,
1261  Int_t nLabels, Int_t & pdg )
1262 {
1263  if(!labels || nLabels<=0)
1264  {
1265  AliWarning(Form("Strange, labels array %p, n labels %d", labels,nLabels));
1266  return kFALSE;
1267  }
1268 
1269  AliDebug(1,Form("Primaries: nlabels %d",nLabels));
1270 
1271  Float_t e = fClusterMomentum.E();
1272  Float_t eta = fClusterMomentum.Eta();
1273  Float_t phi = fClusterMomentum.Phi();
1274  if(phi < 0) phi +=TMath::TwoPi();
1275 
1276  AliAODMCParticle * aodprimary = 0x0;
1277  TParticle * primary = 0x0;
1278 
1279  // Play with the MC stack if available
1280  Int_t label = labels[0];
1281 
1282  if(label < 0)
1283  {
1284  AliDebug(1,Form(" *** bad label ***: label %d", label));
1285  return kFALSE;
1286  }
1287 
1288  Int_t pdg0 =-1; Int_t status = -1; Int_t iMother = -1; Int_t iParent = -1;
1289  Float_t vxMC = 0; Float_t vyMC = 0;
1290  Float_t eMC = 0; //Float_t ptMC= 0;
1291  Float_t phiMC = 0; Float_t etaMC = 0;
1292  Int_t charge = 0;
1293 
1294  // Check the origin.
1295  Int_t tag = GetMCAnalysisUtils()->CheckOrigin(labels,nLabels, GetReader(),GetCalorimeter());
1296 
1297  if ( GetReader()->ReadStack() &&
1298  !GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCUnknown))
1299  { // If MC stack and known tag
1300 
1301  if( label >= GetMCStack()->GetNtrack())
1302  {
1303  AliDebug(1,Form("*** large label ***: label %d, n tracks %d", label, GetMCStack()->GetNtrack()));
1304  return kFALSE;
1305  }
1306 
1307  primary = GetMCStack()->Particle(label);
1308  iMother = label;
1309  pdg0 = TMath::Abs(primary->GetPdgCode());
1310  pdg = pdg0;
1311  status = primary->GetStatusCode();
1312  vxMC = primary->Vx();
1313  vyMC = primary->Vy();
1314  iParent = primary->GetFirstMother();
1315 
1316  AliDebug(1,"Cluster most contributing mother:");
1317  AliDebug(1,Form("\t Mother label %d, pdg %d, %s, status %d, parent %d",iMother, pdg0, primary->GetName(),status, iParent));
1318 
1319 
1320  // Get final particle, no conversion products
1321  if(GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCConversion))
1322  {
1323  // Get the parent
1324  primary = GetMCStack()->Particle(iParent);
1325  pdg = TMath::Abs(primary->GetPdgCode());
1326 
1327  AliDebug(2,"Converted cluster!. Find before conversion:");
1328 
1329  while((pdg == 22 || pdg == 11) && status != 1)
1330  {
1331  Int_t iMotherOrg = iMother;
1332  iMother = iParent;
1333  primary = GetMCStack()->Particle(iMother);
1334  status = primary->GetStatusCode();
1335  pdg = TMath::Abs(primary->GetPdgCode());
1336  iParent = primary->GetFirstMother();
1337 
1338  // If gone too back and non stable, assign the decay photon/electron
1339  // there are other possible decays, ignore them for the moment
1340  if(pdg==111 || pdg==221)
1341  {
1342  primary = GetMCStack()->Particle(iMotherOrg);
1343  break;
1344  }
1345 
1346  if( iParent < 0 )
1347  {
1348  iParent = iMother;
1349  break;
1350  }
1351 
1352  AliDebug(2,Form("\t pdg %d, index %d, %s, status %d",pdg, iMother, primary->GetName(),status));
1353  }
1354 
1355  AliDebug(1,"Converted Cluster mother before conversion:");
1356  AliDebug(1,Form("\t Mother label %d, pdg %d, %s, status %d, parent %d",iMother, pdg, primary->GetName(), status, iParent));
1357 
1358  }
1359 
1360  // Overlapped pi0 (or eta, there will be very few), get the meson
1361  if(GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCPi0) ||
1362  GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCEta))
1363  {
1364  AliDebug(2,"Overlapped Meson decay!, Find it:");
1365 
1366  while(pdg != 111 && pdg != 221)
1367  {
1368  //printf("iMother %d, pdg %d, iParent %d, pdg %d\n",iMother,pdg,iParent,GetMCStack()->Particle(iParent)->GetPdgCode());
1369  iMother = iParent;
1370  primary = GetMCStack()->Particle(iMother);
1371  status = primary->GetStatusCode();
1372  pdg = TMath::Abs(primary->GetPdgCode());
1373  iParent = primary->GetFirstMother();
1374 
1375  if( iParent < 0 ) break;
1376 
1377  AliDebug(2,Form("\t pdg %d, %s, index %d",pdg, primary->GetName(),iMother));
1378 
1379  if(iMother==-1)
1380  {
1381  AliWarning("Tagged as Overlapped photon but meson not found, why?");
1382  //break;
1383  }
1384  }
1385 
1386  AliDebug(2,Form("Overlapped %s decay, label %d",primary->GetName(),iMother));
1387  }
1388 
1389  eMC = primary->Energy();
1390  //ptMC = primary->Pt();
1391  phiMC = primary->Phi();
1392  etaMC = primary->Eta();
1393  pdg = TMath::Abs(primary->GetPdgCode());
1394  charge = (Int_t) TDatabasePDG::Instance()->GetParticle(pdg)->Charge();
1395 
1396  }
1397  else if( GetReader()->ReadAODMCParticles() &&
1398  !GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCUnknown))
1399  {
1400  // If MC AOD and known tag
1401  // Get the list of MC particles
1402  if(!GetReader()->GetAODMCParticles())
1403  AliFatal("MCParticles not available!");
1404 
1405  aodprimary = (AliAODMCParticle*) (GetReader()->GetAODMCParticles())->At(label);
1406  iMother = label;
1407  pdg0 = TMath::Abs(aodprimary->GetPdgCode());
1408  pdg = pdg0;
1409  status = aodprimary->IsPrimary();
1410  vxMC = aodprimary->Xv();
1411  vyMC = aodprimary->Yv();
1412  iParent = aodprimary->GetMother();
1413 
1414  AliDebug(1,"Cluster most contributing mother:");
1415  AliDebug(1,Form("\t Mother label %d, pdg %d, Primary? %d, Physical Primary? %d, parent %d",
1416  iMother, pdg0, aodprimary->IsPrimary(), aodprimary->IsPhysicalPrimary(), iParent));
1417 
1418  //Get final particle, no conversion products
1419  if( GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCConversion) )
1420  {
1421  AliDebug(2,"Converted cluster!. Find before conversion:");
1422 
1423  // Get the parent
1424  aodprimary = (AliAODMCParticle*)(GetReader()->GetAODMCParticles())->At(iParent);
1425  pdg = TMath::Abs(aodprimary->GetPdgCode());
1426 
1427  while ((pdg == 22 || pdg == 11) && !aodprimary->IsPhysicalPrimary())
1428  {
1429  Int_t iMotherOrg = iMother;
1430  iMother = iParent;
1431  aodprimary = (AliAODMCParticle*)(GetReader()->GetAODMCParticles())->At(iMother);
1432  status = aodprimary->IsPrimary();
1433  iParent = aodprimary->GetMother();
1434  pdg = TMath::Abs(aodprimary->GetPdgCode());
1435 
1436  // If gone too back and non stable, assign the decay photon/electron
1437  // there are other possible decays, ignore them for the moment
1438  if( pdg == 111 || pdg == 221 )
1439  {
1440  aodprimary = (AliAODMCParticle*)(GetReader()->GetAODMCParticles())->At(iMotherOrg);
1441  break;
1442  }
1443 
1444  if( iParent < 0 )
1445  {
1446  iParent = iMother;
1447  break;
1448  }
1449 
1450  AliDebug(2,Form("\t pdg %d, index %d, Primary? %d, Physical Primary? %d",
1451  pdg, iMother, aodprimary->IsPrimary(), aodprimary->IsPhysicalPrimary()));
1452  }
1453 
1454  AliDebug(1,"Converted Cluster mother before conversion:");
1455  AliDebug(1,Form("\t Mother label %d, pdg %d, parent %d, Primary? %d, Physical Primary? %d",
1456  iMother, pdg, iParent, aodprimary->IsPrimary(), aodprimary->IsPhysicalPrimary()));
1457 
1458  }
1459 
1460  //Overlapped pi0 (or eta, there will be very few), get the meson
1461  if(GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCPi0) ||
1462  GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCEta))
1463  {
1464  AliDebug(2,Form("Overlapped Meson decay!, Find it: PDG %d, mom %d",pdg, iMother));
1465 
1466  while(pdg != 111 && pdg != 221)
1467  {
1468  iMother = iParent;
1469  aodprimary = (AliAODMCParticle*)(GetReader()->GetAODMCParticles())->At(iMother);
1470  status = aodprimary->IsPrimary();
1471  iParent = aodprimary->GetMother();
1472  pdg = TMath::Abs(aodprimary->GetPdgCode());
1473 
1474  if( iParent < 0 ) break;
1475 
1476  AliDebug(2,Form("\t pdg %d, index %d",pdg, iMother));
1477 
1478  if(iMother==-1)
1479  {
1480  AliWarning("Tagged as Overlapped photon but meson not found, why?");
1481  //break;
1482  }
1483  }
1484 
1485  AliDebug(2,Form("Overlapped %s decay, label %d",aodprimary->GetName(),iMother));
1486  }
1487 
1488  status = aodprimary->IsPrimary();
1489  eMC = aodprimary->E();
1490  //ptMC = aodprimary->Pt();
1491  phiMC = aodprimary->Phi();
1492  etaMC = aodprimary->Eta();
1493  pdg = TMath::Abs(aodprimary->GetPdgCode());
1494  charge = aodprimary->Charge();
1495  }
1496 
1497  //Float_t vz = primary->Vz();
1498  Float_t rVMC = TMath::Sqrt(vxMC*vxMC + vyMC*vyMC);
1499  if( ( pdg == 22 || TMath::Abs(pdg) == 11 ) && status != 1 )
1500  {
1501  fhEMVxyz ->Fill(vxMC, vyMC, GetEventWeight());//,vz);
1502  fhEMR ->Fill(e , rVMC, GetEventWeight());
1503  }
1504 
1505  //printf("reco e %f, pt %f, phi %f, eta %f \n", e, pt, phi, eta);
1506  //printf("prim e %f, pt %f, phi %f, eta %f \n", eMC,ptMC,phiMC ,etaMC );
1507  //printf("vertex: vx %f, vy %f, vz %f, r %f \n", vxMC, vyMC, vz, r);
1508 
1509  // Overlapped pi0 (or eta, there will be very few)
1510  Int_t mcIndex = -1;
1511  if ( GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCPi0 ) )
1512  {
1513  mcIndex = kmcPi0;
1514  } // Overlapped pizero decay
1515  else if( GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCEta ) )
1516  {
1517  mcIndex = kmcEta;
1518  } // Overlapped eta decay
1519  else if( GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCPhoton ) )
1520  {
1521  if( GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCConversion))
1522  mcIndex = kmcPhotonConv ;
1523  else
1524  mcIndex = kmcPhoton ;
1525  } // photon
1526  else if( GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCElectron) )
1527  {
1528  mcIndex = kmcElectron;
1529  fhEMVxyz->Fill(vxMC, vyMC, GetEventWeight());//,vz);
1530  fhEMR ->Fill(e , rVMC, GetEventWeight());
1531  }
1532  else if(charge == 0)
1533  {
1534  mcIndex = kmcNeHadron;
1535  fhHaVxyz->Fill(vxMC, vyMC, GetEventWeight());//,vz);
1536  fhHaR ->Fill(e , rVMC, GetEventWeight());
1537  }
1538  else if(charge!=0)
1539  {
1540  mcIndex = kmcChHadron;
1541  fhHaVxyz->Fill(vxMC, vyMC, GetEventWeight());//,vz);
1542  fhHaR ->Fill(e , rVMC, GetEventWeight());
1543  }
1544 
1545  //printf("mc index %d\n",mcIndex);
1546 
1547  if( mcIndex >= 0 && mcIndex < 7 && e > 0.5 && eMC > 0.5)
1548  {
1549  fhRecoMCE [mcIndex][(matched)]->Fill(e , eMC , GetEventWeight());
1550  fhRecoMCEta [mcIndex][(matched)]->Fill(eta, etaMC , GetEventWeight());
1551  fhRecoMCPhi [mcIndex][(matched)]->Fill(phi, phiMC , GetEventWeight());
1552  fhRecoMCRatioE [mcIndex][(matched)]->Fill(e , e/eMC , GetEventWeight());
1553  fhRecoMCDeltaE [mcIndex][(matched)]->Fill(e , eMC-e , GetEventWeight());
1554  fhRecoMCDeltaPhi[mcIndex][(matched)]->Fill(e , phiMC-phi, GetEventWeight());
1555  fhRecoMCDeltaEta[mcIndex][(matched)]->Fill(e , etaMC-eta, GetEventWeight());
1556  }
1557 
1558  if( primary || aodprimary ) return kTRUE ;
1559  else return kFALSE;
1560 }
1561 
1562 //_________________________________________________________________________________________________________
1567 //_________________________________________________________________________________________________________
1568 void AliAnaCalorimeterQA::ClusterMatchedWithTrackHistograms(AliVCluster *clus, Bool_t okPrimary, Int_t pdg)
1569 {
1570  Float_t e = fClusterMomentum.E();
1571  Float_t pt = fClusterMomentum.Pt();
1572  Float_t eta = fClusterMomentum.Eta();
1573  Float_t phi = fClusterMomentum.Phi();
1574  if(phi < 0) phi +=TMath::TwoPi();
1575 
1576  fhECharged ->Fill(e , GetEventWeight());
1577  fhPtCharged ->Fill(pt , GetEventWeight());
1578  fhPhiCharged ->Fill(phi, GetEventWeight());
1579  fhEtaCharged ->Fill(eta, GetEventWeight());
1580 
1581  if ( fFillAllTH3 ) fhEtaPhiECharged->Fill(eta, phi, e, GetEventWeight());
1582  else if ( e > 0.5 ) fhEtaPhiCharged ->Fill(eta, phi, GetEventWeight());
1583 
1584  // Study the track and matched cluster if track exists.
1585 
1586  AliVTrack *track = GetCaloUtils()->GetMatchedTrack(clus, GetReader()->GetInputEvent());
1587 
1588  if(!track) return ;
1589 
1590  Double_t tpt = track->Pt();
1591  Double_t tmom = track->P();
1592  Double_t dedx = track->GetTPCsignal();
1593  Int_t nITS = track->GetNcls(0);
1594  Int_t nTPC = track->GetNcls(1);
1595  Bool_t positive = kFALSE;
1596  if(track) positive = (track->Charge()>0);
1597 
1598  // Residuals
1599  Float_t deta = clus->GetTrackDz();
1600  Float_t dphi = clus->GetTrackDx();
1601  Double_t dR = TMath::Sqrt(dphi*dphi + deta*deta);
1602  Int_t nModule = GetModuleNumber(clus);
1603 
1604  if( TMath::Abs(dphi) < 999 )
1605  {
1606  if(positive)
1607  {
1608  fhTrackMatchedDEtaPos->Fill(e, deta, GetEventWeight());
1609  fhTrackMatchedDPhiPos->Fill(e, dphi, GetEventWeight());
1610 
1611  if(e > 0.5)
1612  {
1613  fhTrackMatchedDEtaPosMod ->Fill(deta, nModule, GetEventWeight());
1614  fhTrackMatchedDPhiPosMod ->Fill(dphi, nModule, GetEventWeight());
1615  fhTrackMatchedDEtaDPhiPos->Fill(deta, dphi , GetEventWeight());
1616  }
1617  }
1618  else
1619  {
1620  fhTrackMatchedDEtaNeg->Fill(e, deta, GetEventWeight());
1621  fhTrackMatchedDPhiNeg->Fill(e, dphi, GetEventWeight());
1622 
1623  if(e > 0.5)
1624  {
1625  fhTrackMatchedDEtaNegMod ->Fill(deta, nModule, GetEventWeight());
1626  fhTrackMatchedDPhiNegMod ->Fill(dphi, nModule, GetEventWeight());
1627  fhTrackMatchedDEtaDPhiNeg->Fill(deta, dphi , GetEventWeight());
1628  }
1629  }
1630  }
1631 
1632  Double_t eOverP = e/tmom;
1633  fh1EOverP->Fill(tpt, eOverP, GetEventWeight());
1634  if(e > 0.5 && tpt > 0.5) fh1EOverPMod->Fill(eOverP, nModule, GetEventWeight());
1635 
1636  if(dR < 0.02)
1637  {
1638  fh1EOverPR02->Fill(tpt, eOverP, GetEventWeight());
1639  if(e > 0.5 && tpt > 0.5) fh1EOverPR02Mod->Fill(eOverP, nModule, GetEventWeight());
1640 
1641  if(dedx > 60 && dedx < 100)
1642  {
1643  fh1EleEOverP->Fill(tpt, eOverP, GetEventWeight());
1644  if(e > 0.5 && tpt > 0.5) fh1EleEOverPMod->Fill(eOverP, nModule, GetEventWeight());
1645  }
1646  }
1647 
1648  fh2dR->Fill(e, dR, GetEventWeight());
1649  fh2MatchdEdx->Fill(tmom, dedx, GetEventWeight());
1650 
1651  if(e > 0.5 && tmom > 0.5)
1652  {
1653  fh2dRMod->Fill(dR, nModule, GetEventWeight());
1654  fh2MatchdEdxMod->Fill(dedx, nModule, GetEventWeight());
1655  }
1656 
1657  if(dR < 0.02 && eOverP > 0.6 && eOverP < 1.2
1658  && clus->GetNCells() > 1 && nITS > 3 && nTPC > 20)
1659  {
1660  fh2EledEdx->Fill(tmom, dedx, GetEventWeight());
1661  if(e > 0.5 && tmom > 0.5) fh2EledEdxMod->Fill(dedx, nModule, GetEventWeight());
1662  }
1663 
1664  if(IsDataMC() && okPrimary)
1665  {
1666  Double_t charge = TDatabasePDG::Instance()->GetParticle(pdg)->Charge();
1667 
1668  if(TMath::Abs(pdg) == 11)
1669  {
1670  fhMCEle1EOverP ->Fill(tpt , eOverP, GetEventWeight());
1671  fhMCEle1dR ->Fill(dR , GetEventWeight());
1672  fhMCEle2MatchdEdx->Fill(tmom, dedx , GetEventWeight());
1673 
1674  if(dR < 0.02)
1675  {
1676  fhMCEle1EOverPR02->Fill(tpt, eOverP, GetEventWeight());
1677  if(dedx > 60 && dedx < 100) fhMCEle1EleEOverP->Fill(tpt, eOverP, GetEventWeight());
1678  }
1679  }
1680  else if(charge!=0)
1681  {
1682  fhMCChHad1EOverP ->Fill(tpt , eOverP, GetEventWeight());
1683  fhMCChHad1dR ->Fill(dR , GetEventWeight());
1684  fhMCChHad2MatchdEdx->Fill(tmom, dedx , GetEventWeight());
1685 
1686  if(dR < 0.02)
1687  {
1688  fhMCChHad1EOverPR02->Fill(tpt, eOverP, GetEventWeight());
1689  if(dedx > 60 && dedx < 100) fhMCChHad1EleEOverP->Fill(tpt, eOverP, GetEventWeight());
1690  }
1691  }
1692  else if(charge == 0)
1693  {
1694  fhMCNeutral1EOverP ->Fill(tpt , eOverP, GetEventWeight());
1695  fhMCNeutral1dR ->Fill(dR , GetEventWeight());
1696  fhMCNeutral2MatchdEdx->Fill(tmom, dedx , GetEventWeight());
1697 
1698  if(dR < 0.02)
1699  {
1700  fhMCNeutral1EOverPR02->Fill(tpt, eOverP, GetEventWeight());
1701  if(dedx > 60 && dedx < 100) fhMCNeutral1EleEOverP->Fill(tpt, eOverP, GetEventWeight());
1702  }
1703  }
1704  } // DataMC
1705 }
1706 
1707 //___________________________________
1710 //___________________________________
1711 void AliAnaCalorimeterQA::Correlate()
1712 {
1713  // Clusters arrays
1714  TObjArray * caloClustersEMCAL = GetEMCALClusters();
1715  TObjArray * caloClustersPHOS = GetPHOSClusters();
1716 
1717  if(!caloClustersEMCAL || !caloClustersPHOS)
1718  {
1719  AliDebug(1,Form("PHOS (%p) or EMCAL (%p) clusters array not available, do not correlate",caloClustersPHOS,caloClustersEMCAL));
1720  return ;
1721  }
1722 
1723  // Cells arrays
1724  AliVCaloCells * cellsEMCAL = GetEMCALCells();
1725  AliVCaloCells * cellsPHOS = GetPHOSCells();
1726 
1727  if(!cellsEMCAL || !cellsPHOS)
1728  {
1729  AliDebug(1,Form("PHOS (%p) or EMCAL (%p) cells array ot available, do not correlate",cellsPHOS,cellsEMCAL));
1730  return ;
1731  }
1732 
1733  // Clusters parameters
1734  Int_t nclEMCAL = 0;
1735  Int_t nclDCAL = 0;
1736  Int_t nclPHOS = 0;
1737 
1738  Float_t sumClusterEnergyEMCAL = 0;
1739  Float_t sumClusterEnergyDCAL = 0;
1740  Float_t sumClusterEnergyPHOS = 0;
1741 
1742  Int_t iclus = 0;
1743  Float_t energy = 0;
1744  AliVCluster* cluster = 0;
1745  for(iclus = 0 ; iclus < caloClustersEMCAL->GetEntriesFast() ; iclus++)
1746  {
1747  cluster = (AliVCluster*)caloClustersEMCAL->At(iclus);
1748  Float_t energy = cluster->E();
1749 
1750  if( energy < 0.5 ) continue;
1751 
1752  if(cluster->GetCellsAbsId()[0] < 12288)
1753  {
1754  nclEMCAL++;
1755  sumClusterEnergyEMCAL += energy;
1756  }
1757  else
1758  {
1759  nclDCAL++;
1760  sumClusterEnergyDCAL += energy;
1761  }
1762  }
1763 
1764  for(iclus = 0 ; iclus < caloClustersPHOS ->GetEntriesFast(); iclus++)
1765  {
1766  cluster = (AliVCluster*) caloClustersPHOS->At(iclus);
1767 
1768  energy = cluster->E();
1769 
1770  if( energy < 0.5 ) continue;
1771 
1772  nclPHOS++;
1773  sumClusterEnergyPHOS += energy;
1774  }
1775 
1776  // Cells parameters
1777  Int_t ncellsEMCAL = 0 ;
1778  Int_t ncellsDCAL = 0 ;
1779  Int_t ncellsPHOS = 0;
1780 
1781  Float_t sumCellEnergyEMCAL = 0;
1782  Float_t sumCellEnergyDCAL = 0;
1783  Float_t sumCellEnergyPHOS = 0;
1784  Int_t icell = 0;
1785  for(icell = 0 ; icell < cellsEMCAL->GetNumberOfCells() ; icell++)
1786  {
1787  Float_t amp = cellsEMCAL->GetAmplitude(icell);
1788  Int_t cellId = cellsEMCAL->GetCellNumber(icell);
1789 
1790  if (amp < fEMCALCellAmpMin) continue;
1791 
1792  if( cellId < 12288)
1793  {
1794  ncellsEMCAL++;
1795  sumCellEnergyEMCAL += amp;
1796  }
1797  else
1798  {
1799  ncellsDCAL++;
1800  sumCellEnergyDCAL += amp;
1801  }
1802  }
1803 
1804  for(icell = 0 ; icell < cellsPHOS->GetNumberOfCells(); icell++)
1805  {
1806  Float_t amp = cellsPHOS->GetAmplitude(icell);
1807  Int_t cellId = cellsPHOS->GetCellNumber(icell);
1808 
1809  if ( cellId < 0 ) continue ; // CPV
1810 
1811  if ( amp < fPHOSCellAmpMin ) continue;
1812 
1813  ncellsPHOS++;
1814  sumCellEnergyPHOS += amp;
1815  }
1816 
1817  // Fill Histograms
1818 
1819  fhEMCALPHOSCorrNClusters->Fill(nclEMCAL , nclPHOS , GetEventWeight());
1820  fhEMCALPHOSCorrEClusters->Fill(sumClusterEnergyEMCAL, sumClusterEnergyPHOS, GetEventWeight());
1821  fhEMCALPHOSCorrNCells ->Fill(ncellsEMCAL , ncellsPHOS , GetEventWeight());
1822  fhEMCALPHOSCorrECells ->Fill(sumCellEnergyEMCAL , sumCellEnergyPHOS , GetEventWeight());
1823 
1824  fhEMCALDCALCorrNClusters->Fill(nclEMCAL , nclDCAL , GetEventWeight());
1825  fhEMCALDCALCorrEClusters->Fill(sumClusterEnergyEMCAL, sumClusterEnergyDCAL, GetEventWeight());
1826  fhEMCALDCALCorrNCells ->Fill(ncellsEMCAL , ncellsDCAL , GetEventWeight());
1827  fhEMCALDCALCorrECells ->Fill(sumCellEnergyEMCAL , sumCellEnergyDCAL , GetEventWeight());
1828 
1829  fhDCALPHOSCorrNClusters ->Fill(nclDCAL , nclPHOS , GetEventWeight());
1830  fhDCALPHOSCorrEClusters ->Fill(sumClusterEnergyDCAL , sumClusterEnergyPHOS, GetEventWeight());
1831  fhDCALPHOSCorrNCells ->Fill(ncellsDCAL , ncellsPHOS , GetEventWeight());
1832  fhDCALPHOSCorrECells ->Fill(sumCellEnergyDCAL , sumCellEnergyPHOS , GetEventWeight());
1833 
1834  Int_t v0S = GetV0Signal(0) + GetV0Signal(1);
1835  Int_t v0M = GetV0Multiplicity(0) + GetV0Multiplicity(1);
1836  Int_t trM = GetTrackMultiplicity();
1837  Float_t cen = GetEventCentrality();
1838  Float_t ep = GetEventPlaneAngle();
1839 
1840  Int_t ncl = nclPHOS;
1841  Float_t sumClusterEnergy = sumClusterEnergyPHOS;
1842  Int_t ncells = ncellsPHOS;
1843  Float_t sumCellEnergy = sumCellEnergyPHOS;
1844 
1845  if ( GetCalorimeter() == kEMCAL )
1846  {
1847  ncl = nclEMCAL + nclDCAL;
1848  sumClusterEnergy = sumClusterEnergyEMCAL + sumClusterEnergyDCAL;
1849  ncells = ncellsEMCAL + ncellsDCAL;
1850  sumCellEnergy = sumCellEnergyEMCAL + sumCellEnergyDCAL;
1851  }
1852 
1853  fhCaloV0MCorrNClusters ->Fill(v0M, ncl , GetEventWeight());
1854  fhCaloV0MCorrEClusters ->Fill(v0M, sumClusterEnergy, GetEventWeight());
1855  fhCaloV0MCorrNCells ->Fill(v0M, ncells , GetEventWeight());
1856  fhCaloV0MCorrECells ->Fill(v0M, sumCellEnergy , GetEventWeight());
1857 
1858  fhCaloV0SCorrNClusters ->Fill(v0S, ncl , GetEventWeight());
1859  fhCaloV0SCorrEClusters ->Fill(v0S, sumClusterEnergy, GetEventWeight());
1860  fhCaloV0SCorrNCells ->Fill(v0S, ncells , GetEventWeight());
1861  fhCaloV0SCorrECells ->Fill(v0S, sumCellEnergy , GetEventWeight());
1862 
1863  fhCaloTrackMCorrNClusters->Fill(trM, ncl , GetEventWeight());
1864  fhCaloTrackMCorrEClusters->Fill(trM, sumClusterEnergy, GetEventWeight());
1865  fhCaloTrackMCorrNCells ->Fill(trM, ncells , GetEventWeight());
1866  fhCaloTrackMCorrECells ->Fill(trM, sumCellEnergy , GetEventWeight());
1867 
1868  fhCaloCenNClusters ->Fill(cen, ncl , GetEventWeight());
1869  fhCaloCenEClusters ->Fill(cen, sumClusterEnergy, GetEventWeight());
1870  fhCaloCenNCells ->Fill(cen, ncells , GetEventWeight());
1871  fhCaloCenECells ->Fill(cen, sumCellEnergy , GetEventWeight());
1872 
1873  fhCaloEvPNClusters ->Fill(ep , ncl , GetEventWeight());
1874  fhCaloEvPEClusters ->Fill(ep , sumClusterEnergy, GetEventWeight());
1875  fhCaloEvPNCells ->Fill(ep , ncells , GetEventWeight());
1876  fhCaloEvPECells ->Fill(ep , sumCellEnergy , GetEventWeight());
1877 
1878  AliDebug(1,"Correlate():");
1879  AliDebug(1,Form("\t EMCAL: N cells %d, N clusters %d, summed E cells %f, summed E clusters %f",
1880  ncellsEMCAL,nclEMCAL, sumCellEnergyEMCAL,sumClusterEnergyEMCAL));
1881  AliDebug(1,Form("\t DCAL : N cells %d, N clusters %d, summed E cells %f, summed E clusters %f",
1882  ncellsDCAL,nclDCAL, sumCellEnergyDCAL,sumClusterEnergyDCAL));
1883  AliDebug(1,Form("\t PHOS : N cells %d, N clusters %d, summed E cells %f, summed E clusters %f",
1884  ncellsPHOS,nclPHOS,sumCellEnergyPHOS,sumClusterEnergyPHOS));
1885  AliDebug(1,Form("\t V0 : Signal %d, Multiplicity %d, Track Multiplicity %d", v0S,v0M,trM));
1886  AliDebug(1,Form("\t centrality : %f, Event plane angle %f", cen,ep));
1887 }
1888 
1889 //_________________________________________________
1891 //_________________________________________________
1892 TObjString * AliAnaCalorimeterQA::GetAnalysisCuts()
1893 {
1894  TString parList ; //this will be list of parameters used for this analysis.
1895  const Int_t buffersize = 255;
1896  char onePar[buffersize] ;
1897 
1898  snprintf(onePar,buffersize,"--- AliAnaCalorimeterQA ---:") ;
1899  parList+=onePar ;
1900  snprintf(onePar,buffersize,"Calorimeter: %s;",GetCalorimeterString().Data()) ;
1901  parList+=onePar ;
1902  snprintf(onePar,buffersize,"Time Cut : %2.2f < T < %2.2f ns;",fTimeCutMin, fTimeCutMax) ;
1903  parList+=onePar ;
1904  snprintf(onePar,buffersize,"Cell Amplitude: PHOS > %2.2f GeV, EMCAL > %2.2f GeV;",fPHOSCellAmpMin, fEMCALCellAmpMin) ;
1905  parList+=onePar ;
1906  snprintf(onePar,buffersize,"N cells per cluster: PHOS >= %d, EMCAL >= %d;",fPHOSClusterNCellMin, fEMCALClusterNCellMin) ;
1907  parList+=onePar ;
1908  snprintf(onePar,buffersize,"Inv. Mass %2.1f < E_cl < %2.1f GeV;",fInvMassMinECut, fInvMassMaxECut) ;
1909  parList+=onePar ;
1910  snprintf(onePar,buffersize,"Inv. Mass %2.1f < M02 < %2.1f GeV;",fInvMassMinM02Cut, fInvMassMaxM02Cut) ;
1911  parList+=onePar ;
1912  snprintf(onePar,buffersize,"Cluster pair opening angle < %2.1f rad;",fInvMassMaxOpenAngle) ;
1913  parList+=onePar ;
1914  snprintf(onePar,buffersize,"Cluster pair time difference < %2.1f rad;",fInvMassMaxTimeDifference) ;
1915  parList+=onePar ;
1916 
1917  //Get parameters set in base class.
1918  //parList += GetBaseParametersList() ;
1919 
1920  //Get parameters set in FiducialCut class (not available yet)
1921  //parlist += GetFidCut()->GetFidCutParametersList()
1922 
1923  return new TObjString(parList) ;
1924 }
1925 
1926 //_________________________________________________________________________________
1928 //_________________________________________________________________________________
1929 void AliAnaCalorimeterQA::ExoticHistograms(Int_t absIdMax, Float_t ampMax,
1930  AliVCluster *clus, AliVCaloCells* cells)
1931 {
1932  if(ampMax < 0.01)
1933  {
1934  AliDebug(1,Form("Low amplitude energy %f",ampMax));
1935  return;
1936  }
1937 
1938  Float_t l0 = clus->GetM02();
1939  Float_t l1 = clus->GetM20();
1940  Float_t en = clus->E();
1941  Int_t nc = clus->GetNCells();
1942  Double_t tmax = clus->GetTOF()*1.e9-fConstantTimeShift; // recalibrated elsewhere
1943 
1944  Float_t eCrossFrac = 1-GetECross(absIdMax,cells, 10000000)/ampMax;
1945 
1946  if(en > 5)
1947  {
1948  fhExoL0ECross->Fill(eCrossFrac, l0, GetEventWeight());
1949  fhExoL1ECross->Fill(eCrossFrac, l1, GetEventWeight());
1950  }
1951 
1952  for(Int_t ie = 0; ie < fExoNECrossCuts; ie++)
1953  {
1954  for(Int_t idt = 0; idt < fExoNDTimeCuts; idt++)
1955  {
1956  eCrossFrac = 1-GetECross(absIdMax,cells, fExoDTimeCuts[idt])/ampMax;
1957 
1958  if(eCrossFrac > fExoECrossCuts[ie])
1959  {
1960  // Exotic
1961  fhExoL0 [ie][idt]->Fill(en,l0 , GetEventWeight());
1962  fhExoL1 [ie][idt]->Fill(en,l1 , GetEventWeight());
1963  fhExoTime [ie][idt]->Fill(en,tmax, GetEventWeight());
1964 
1965  if(en > 5)
1966  {
1967  fhExoL0NCell[ie][idt]->Fill(nc, l0, GetEventWeight());
1968  fhExoL1NCell[ie][idt]->Fill(nc, l1, GetEventWeight());
1969  }
1970 
1971  // Diff time, do for one cut in time
1972  if(idt == 0)
1973  {
1974  for (Int_t icell = 0; icell < clus->GetNCells(); icell++)
1975  {
1976  Int_t absId = clus->GetCellsAbsId()[icell];
1977  Double_t time = cells->GetCellTime(absId);
1978  GetCaloUtils()->RecalibrateCellTime(time, GetCalorimeter(), absId,GetReader()->GetInputEvent()->GetBunchCrossNumber());
1979  time *= 1e9;
1980  time -= fConstantTimeShift;
1981 
1982  Float_t diff = tmax-time;
1983  fhExoDTime[ie]->Fill(en, diff, GetEventWeight());
1984  }
1985  }
1986  }
1987  else
1988  {
1989  fhExoECross[ie][idt]->Fill(en, eCrossFrac, GetEventWeight());
1990  fhExoNCell [ie][idt]->Fill(en, nc , GetEventWeight());
1991  }
1992  } // D time cut loop
1993  } // e cross cut loop
1994 }
1995 
1996 //___________________________________________________
1999 //___________________________________________________
2000 TList * AliAnaCalorimeterQA::GetCreateOutputObjects()
2001 {
2002  TList * outputContainer = new TList() ;
2003  outputContainer->SetName("QAHistos") ;
2004 
2005  // Init the number of modules, set in the class AliCalorimeterUtils
2006 
2007  fNModules = GetCaloUtils()->GetNumberOfSuperModulesUsed();
2008  if(GetCalorimeter()==kPHOS && fNModules > 4) fNModules = 4;
2009 
2010  // Histogram binning and ranges
2011 
2012  Int_t nptbins = GetHistogramRanges()->GetHistoPtBins(); Float_t ptmax = GetHistogramRanges()->GetHistoPtMax(); Float_t ptmin = GetHistogramRanges()->GetHistoPtMin();
2013  Int_t nfineptbins = GetHistogramRanges()->GetHistoFinePtBins(); Float_t ptfinemax = GetHistogramRanges()->GetHistoFinePtMax(); Float_t ptfinemin = GetHistogramRanges()->GetHistoFinePtMin();
2014  Int_t nphibins = GetHistogramRanges()->GetHistoPhiBins(); Float_t phimax = GetHistogramRanges()->GetHistoPhiMax(); Float_t phimin = GetHistogramRanges()->GetHistoPhiMin();
2015  Int_t netabins = GetHistogramRanges()->GetHistoEtaBins(); Float_t etamax = GetHistogramRanges()->GetHistoEtaMax(); Float_t etamin = GetHistogramRanges()->GetHistoEtaMin();
2016  Int_t nmassbins = GetHistogramRanges()->GetHistoMassBins(); Float_t massmax = GetHistogramRanges()->GetHistoMassMax(); Float_t massmin = GetHistogramRanges()->GetHistoMassMin();
2017  Int_t nasymbins = GetHistogramRanges()->GetHistoAsymmetryBins(); Float_t asymmax = GetHistogramRanges()->GetHistoAsymmetryMax(); Float_t asymmin = GetHistogramRanges()->GetHistoAsymmetryMin();
2018  Int_t nPoverEbins = GetHistogramRanges()->GetHistoPOverEBins(); Float_t eOverPmax = GetHistogramRanges()->GetHistoPOverEMax(); Float_t eOverPmin = GetHistogramRanges()->GetHistoPOverEMin();
2019  Int_t ndedxbins = GetHistogramRanges()->GetHistodEdxBins(); Float_t dedxmax = GetHistogramRanges()->GetHistodEdxMax(); Float_t dedxmin = GetHistogramRanges()->GetHistodEdxMin();
2020  Int_t ndRbins = GetHistogramRanges()->GetHistodRBins(); Float_t dRmax = GetHistogramRanges()->GetHistodRMax(); Float_t dRmin = GetHistogramRanges()->GetHistodRMin();
2021  Int_t ntimebins = GetHistogramRanges()->GetHistoTimeBins(); Float_t timemax = GetHistogramRanges()->GetHistoTimeMax(); Float_t timemin = GetHistogramRanges()->GetHistoTimeMin();
2022  Int_t nclbins = GetHistogramRanges()->GetHistoNClustersBins(); Int_t nclmax = GetHistogramRanges()->GetHistoNClustersMax(); Int_t nclmin = GetHistogramRanges()->GetHistoNClustersMin();
2023  Int_t ncebins = GetHistogramRanges()->GetHistoNCellsBins(); Int_t ncemax = GetHistogramRanges()->GetHistoNCellsMax(); Int_t ncemin = GetHistogramRanges()->GetHistoNCellsMin();
2024  Int_t nceclbins = GetHistogramRanges()->GetHistoNClusterCellBins(); Int_t nceclmax = GetHistogramRanges()->GetHistoNClusterCellMax(); Int_t nceclmin = GetHistogramRanges()->GetHistoNClusterCellMin();
2025  Int_t nvdistbins = GetHistogramRanges()->GetHistoVertexDistBins(); Float_t vdistmax = GetHistogramRanges()->GetHistoVertexDistMax(); Float_t vdistmin = GetHistogramRanges()->GetHistoVertexDistMin();
2026  Int_t rbins = GetHistogramRanges()->GetHistoRBins(); Float_t rmax = GetHistogramRanges()->GetHistoRMax(); Float_t rmin = GetHistogramRanges()->GetHistoRMin();
2027  Int_t xbins = GetHistogramRanges()->GetHistoXBins(); Float_t xmax = GetHistogramRanges()->GetHistoXMax(); Float_t xmin = GetHistogramRanges()->GetHistoXMin();
2028  Int_t ybins = GetHistogramRanges()->GetHistoYBins(); Float_t ymax = GetHistogramRanges()->GetHistoYMax(); Float_t ymin = GetHistogramRanges()->GetHistoYMin();
2029  Int_t zbins = GetHistogramRanges()->GetHistoZBins(); Float_t zmax = GetHistogramRanges()->GetHistoZMax(); Float_t zmin = GetHistogramRanges()->GetHistoZMin();
2030  Int_t ssbins = GetHistogramRanges()->GetHistoShowerShapeBins(); Float_t ssmax = GetHistogramRanges()->GetHistoShowerShapeMax(); Float_t ssmin = GetHistogramRanges()->GetHistoShowerShapeMin();
2031  Int_t tdbins = GetHistogramRanges()->GetHistoDiffTimeBins() ; Float_t tdmax = GetHistogramRanges()->GetHistoDiffTimeMax(); Float_t tdmin = GetHistogramRanges()->GetHistoDiffTimeMin();
2032 
2033  Int_t nv0sbins = GetHistogramRanges()->GetHistoV0SignalBins(); Int_t nv0smax = GetHistogramRanges()->GetHistoV0SignalMax(); Int_t nv0smin = GetHistogramRanges()->GetHistoV0SignalMin();
2034  Int_t nv0mbins = GetHistogramRanges()->GetHistoV0MultiplicityBins(); Int_t nv0mmax = GetHistogramRanges()->GetHistoV0MultiplicityMax(); Int_t nv0mmin = GetHistogramRanges()->GetHistoV0MultiplicityMin();
2035  Int_t ntrmbins = GetHistogramRanges()->GetHistoTrackMultiplicityBins(); Int_t ntrmmax = GetHistogramRanges()->GetHistoTrackMultiplicityMax(); Int_t ntrmmin = GetHistogramRanges()->GetHistoTrackMultiplicityMin();
2036 
2037  // Calorimeter grid
2038 
2039  // EMCAL
2040  fNMaxCols = 48;
2041  fNMaxRows = 24;
2042  fNRCU = 2 ;
2043  // PHOS
2044  if(GetCalorimeter()==kPHOS)
2045  {
2046  fNMaxCols = 56;
2047  fNMaxRows = 64;
2048  fNRCU = 4 ;
2049  }
2050 
2051 
2052  // Init histograms
2053 
2054  fhE = new TH1F ("hE","#it{E} reconstructed clusters ", nptbins*5,ptmin,ptmax*5);
2055  fhE->SetXTitle("#it{E} (GeV)");
2056  outputContainer->Add(fhE);
2057 
2058  fhPt = new TH1F ("hPt","#it{p}_{T} reconstructed clusters", nptbins,ptmin,ptmax);
2059  fhPt->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2060  outputContainer->Add(fhPt);
2061 
2062  fhPhi = new TH1F ("hPhi","#phi reconstructed clusters ",nphibins,phimin,phimax);
2063  fhPhi->SetXTitle("#phi (rad)");
2064  outputContainer->Add(fhPhi);
2065 
2066  fhEta = new TH1F ("hEta","#eta reconstructed clusters ",netabins,etamin,etamax);
2067  fhEta->SetXTitle("#eta ");
2068  outputContainer->Add(fhEta);
2069 
2070  if(!fFillEBinAcceptanceHisto)
2071  {
2072  if(fFillAllTH3)
2073  {
2074  fhEtaPhiE = new TH3F ("hEtaPhiE","#eta vs #phi vs energy, reconstructed clusters",
2075  netabins,etamin,etamax,nphibins,phimin,phimax,nptbins,ptmin,ptmax);
2076  fhEtaPhiE->SetXTitle("#eta ");
2077  fhEtaPhiE->SetYTitle("#phi (rad)");
2078  fhEtaPhiE->SetZTitle("#it{E} (GeV) ");
2079  outputContainer->Add(fhEtaPhiE);
2080  }
2081  else
2082  {
2083  fhEtaPhi = new TH2F ("hEtaPhi","#eta vs #phi for #it{E} > 0.5 GeV, reconstructed clusters",
2084  netabins,etamin,etamax,nphibins,phimin,phimax);
2085  fhEtaPhi->SetXTitle("#eta ");
2086  fhEtaPhi->SetYTitle("#phi (rad)");
2087  outputContainer->Add(fhEtaPhi);
2088  }
2089  }
2090 
2091  fhClusterTimeEnergy = new TH2F ("hClusterTimeEnergy","energy vs TOF, reconstructed clusters",
2092  nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
2093  fhClusterTimeEnergy->SetXTitle("#it{E} (GeV) ");
2094  fhClusterTimeEnergy->SetYTitle("TOF (ns)");
2095  outputContainer->Add(fhClusterTimeEnergy);
2096 
2097  fhClusterPairDiffTimeE = new TH2F("hClusterPairDiffTimeE","cluster pair time difference vs E, only good clusters",
2098  nptbins,ptmin,ptmax, tdbins,tdmin,tdmax);
2099  fhClusterPairDiffTimeE->SetXTitle("#it{E}_{cluster} (GeV)");
2100  fhClusterPairDiffTimeE->SetYTitle("#Delta #it{t} (ns)");
2101  outputContainer->Add(fhClusterPairDiffTimeE);
2102 
2103  fhClusterPairDiffTimeESameMod = new TH2F("hClusterPairDiffTimeESameMod","cluster pair time difference vs E, only good clusters",
2104  nptbins,ptmin,ptmax, tdbins,tdmin,tdmax);
2105  fhClusterPairDiffTimeESameMod->SetXTitle("#it{E}_{cluster} (GeV)");
2106  fhClusterPairDiffTimeESameMod->SetYTitle("#Delta #it{t} (ns)");
2107  outputContainer->Add(fhClusterPairDiffTimeESameMod);
2108 
2109  fhLambda0 = new TH2F ("hLambda0","shower shape, #lambda^{2}_{0} vs E",
2110  nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
2111  fhLambda0->SetXTitle("#it{E}_{cluster}");
2112  fhLambda0->SetYTitle("#lambda^{2}_{0}");
2113  outputContainer->Add(fhLambda0);
2114 
2115  fhLambda1 = new TH2F ("hLambda1","shower shape, #lambda^{2}_{1} vs E for bad cluster ",
2116  nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
2117  fhLambda1->SetXTitle("#it{E}_{cluster}");
2118  fhLambda1->SetYTitle("#lambda^{2}_{1}");
2119  outputContainer->Add(fhLambda1);
2120 
2121 // fhDispersion = new TH2F ("hDispersion","shower shape, Dispersion^{2} vs E for bad cluster ",
2122 // nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
2123 // fhDispersion->SetXTitle("#it{E}_{cluster}");
2124 // fhDispersion->SetYTitle("Dispersion");
2125 // outputContainer->Add(fhDispersion);
2126 
2127  fhClusterMaxCellCloseCellRatio = new TH2F ("hClusterMaxCellCloseCellRatio","energy vs ratio of max cell / neighbour cell, reconstructed clusters",
2128  nptbins,ptmin,ptmax, 100,0,1.);
2129  fhClusterMaxCellCloseCellRatio->SetXTitle("#it{E}_{cluster} (GeV) ");
2130  fhClusterMaxCellCloseCellRatio->SetYTitle("#it{E}_{cell i}/#it{E}_{cell max}");
2131  outputContainer->Add(fhClusterMaxCellCloseCellRatio);
2132 
2133  fhClusterMaxCellCloseCellDiff = new TH2F ("hClusterMaxCellCloseCellDiff","energy vs ratio of max cell / neighbour cell, reconstructed clusters",
2134  nptbins,ptmin,ptmax, 500,0,100.);
2135  fhClusterMaxCellCloseCellDiff->SetXTitle("#it{E}_{cluster} (GeV) ");
2136  fhClusterMaxCellCloseCellDiff->SetYTitle("#it{E}_{cell max}-#it{E}_{cell i} (GeV)");
2137  outputContainer->Add(fhClusterMaxCellCloseCellDiff);
2138 
2139  fhClusterMaxCellDiff = new TH2F ("hClusterMaxCellDiff","energy vs difference of cluster energy - max cell energy / cluster energy, good clusters",
2140  nptbins,ptmin,ptmax, 500,0,1.);
2141  fhClusterMaxCellDiff->SetXTitle("#it{E}_{cluster} (GeV) ");
2142  fhClusterMaxCellDiff->SetYTitle("(#it{E}_{cluster} - #it{E}_{cell max})/ #it{E}_{cluster}");
2143  outputContainer->Add(fhClusterMaxCellDiff);
2144 
2145  if(fStudyBadClusters)
2146  {
2147  fhClusterMaxCellDiffNoCut = new TH2F ("hClusterMaxCellDiffNoCut","energy vs difference of cluster energy - max cell energy / cluster energy",
2148  nptbins,ptmin,ptmax, 500,0,1.);
2149  fhClusterMaxCellDiffNoCut->SetXTitle("#it{E}_{cluster} (GeV) ");
2150  fhClusterMaxCellDiffNoCut->SetYTitle("(#it{E}_{cluster} - #it{E}_{cell max})/ #it{E}_{cluster}");
2151  outputContainer->Add(fhClusterMaxCellDiffNoCut);
2152  }
2153 
2154  fhClusterMaxCellECross = new TH2F ("hClusterMaxCellECross","1 - Energy in cross around max energy cell / max energy cell vs cluster energy, good clusters",
2155  nptbins,ptmin,ptmax, 400,-1,1.);
2156  fhClusterMaxCellECross->SetXTitle("#it{E}_{cluster} (GeV) ");
2157  fhClusterMaxCellECross->SetYTitle("1- #it{E}_{cross}/#it{E}_{cell max}");
2158  outputContainer->Add(fhClusterMaxCellECross);
2159 
2160  if(fStudyBadClusters)
2161  {
2162  fhNCellsPerClusterNoCut = new TH2F ("hNCellsPerClusterNoCut","# cells per cluster vs energy, no bad clusters cut",
2163  nptbins,ptmin,ptmax, nceclbins,nceclmin,nceclmax);
2164  fhNCellsPerClusterNoCut->SetXTitle("#it{E} (GeV)");
2165  fhNCellsPerClusterNoCut->SetYTitle("#it{n}_{cells}");
2166  outputContainer->Add(fhNCellsPerClusterNoCut);
2167  }
2168 
2169  fhNCellsPerCluster = new TH2F ("hNCellsPerCluster","# cells per cluster vs energy",nptbins,ptmin,ptmax, nceclbins,nceclmin,nceclmax);
2170  fhNCellsPerCluster->SetXTitle("#it{E} (GeV)");
2171  fhNCellsPerCluster->SetYTitle("#it{n}_{cells}");
2172  outputContainer->Add(fhNCellsPerCluster);
2173 
2174  if(fStudyBadClusters)
2175  {
2176  fhNCellsPerClusterWeirdModNoCut = new TH2F ("hNCellsPerClusterWeirdNoCutMod","# cells per cluster vs energy, E > 100, no bad clusters cut, per SM",
2177  nceclbins,nceclmin,nceclmax,fNModules,0,fNModules);
2178  fhNCellsPerClusterWeirdModNoCut->SetYTitle("SM number");
2179  fhNCellsPerClusterWeirdModNoCut->SetXTitle("#it{n}_{cells}");
2180  outputContainer->Add(fhNCellsPerClusterWeirdModNoCut);
2181  }
2182 
2183  fhNCellsPerClusterWeirdMod = new TH2F ("hNCellsPerClusterWeirdMod","# cells per cluster, E > 100 GeV, per SM",
2184  nceclbins*2,nceclmin,nceclmax*2,fNModules,0,fNModules);
2185  fhNCellsPerClusterWeirdMod->SetYTitle("SM number");
2186  fhNCellsPerClusterWeirdMod->SetXTitle("#it{n}_{cells}");
2187  outputContainer->Add(fhNCellsPerClusterWeirdMod);
2188 
2189  fhNClusters = new TH1F ("hNClusters","# clusters", nclbins,nclmin,nclmax);
2190  fhNClusters->SetXTitle("#it{n}_{clusters}");
2191  outputContainer->Add(fhNClusters);
2192 
2193  if(fStudyBadClusters)
2194  {
2195  fhBadClusterEnergy = new TH1F ("hBadClusterEnergy","Bad cluster energy", nptbins,ptmin,ptmax);
2196  fhBadClusterEnergy->SetXTitle("#it{E}_{cluster} (GeV) ");
2197  outputContainer->Add(fhBadClusterEnergy);
2198 
2199  fhBadClusterEtaPhi = new TH2F ("hBadClusterEtaPhi","Bad cluster, #eta vs #phi, #it{E} > 0.5 GeV",
2200  netabins,etamin,etamax,nphibins,phimin,phimax);
2201  fhBadClusterEtaPhi->SetXTitle("#eta ");
2202  fhBadClusterEtaPhi->SetXTitle("#phi (rad) ");
2203  outputContainer->Add(fhBadClusterEtaPhi);
2204 
2205  fhBadClusterLambda0 = new TH2F ("hBadClusterLambda0","Bad cluster,shower shape, #lambda^{2}_{0} vs E",
2206  nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
2207  fhBadClusterLambda0->SetXTitle("#it{E}_{cluster}");
2208  fhBadClusterLambda0->SetYTitle("#lambda^{2}_{0}");
2209  outputContainer->Add(fhBadClusterLambda0);
2210 
2211  fhBadClusterLambda1 = new TH2F ("hBadClusterLambda1","Bad cluster,shower shape, #lambda^{2}_{1} vs E for bad cluster ",
2212  nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
2213  fhBadClusterLambda1->SetXTitle("#it{E}_{cluster}");
2214  fhBadClusterLambda1->SetYTitle("#lambda^{2}_{1}");
2215  outputContainer->Add(fhBadClusterLambda1);
2216 
2217  fhBadClusterMaxCellCloseCellRatio = new TH2F ("hBadClusterMaxCellCloseCellRatio","energy vs ratio of max cell / neighbour cell constributing cell, reconstructed bad clusters",
2218  nptbins,ptmin,ptmax, 100,0,1.);
2219  fhBadClusterMaxCellCloseCellRatio->SetXTitle("#it{E}_{cluster} (GeV) ");
2220  fhBadClusterMaxCellCloseCellRatio->SetYTitle("ratio");
2221  outputContainer->Add(fhBadClusterMaxCellCloseCellRatio);
2222 
2223  fhBadClusterMaxCellCloseCellDiff = new TH2F ("hBadClusterMaxCellCloseCellDiff","energy vs ratio of max cell - neighbour cell constributing cell, reconstructed bad clusters",
2224  nptbins,ptmin,ptmax, 500,0,100);
2225  fhBadClusterMaxCellCloseCellDiff->SetXTitle("#it{E}_{cluster} (GeV) ");
2226  fhBadClusterMaxCellCloseCellDiff->SetYTitle("#it{E}_{cell max} - #it{E}_{cell i} (GeV)");
2227  outputContainer->Add(fhBadClusterMaxCellCloseCellDiff);
2228 
2229  fhBadClusterMaxCellDiff = new TH2F ("hBadClusterMaxCellDiff","energy vs difference of cluster energy - max cell energy / cluster energy for bad clusters",
2230  nptbins,ptmin,ptmax, 500,0,1.);
2231  fhBadClusterMaxCellDiff->SetXTitle("#it{E}_{cluster} (GeV) ");
2232  fhBadClusterMaxCellDiff->SetYTitle("(#it{E}_{cluster} - #it{E}_{cell max}) / #it{E}_{cluster}");
2233  outputContainer->Add(fhBadClusterMaxCellDiff);
2234 
2235  fhBadClusterTimeEnergy = new TH2F ("hBadClusterTimeEnergy","energy vs TOF of reconstructed bad clusters",
2236  nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
2237  fhBadClusterTimeEnergy->SetXTitle("#it{E}_{cluster} (GeV) ");
2238  fhBadClusterTimeEnergy->SetYTitle("#it{t} (ns)");
2239  outputContainer->Add(fhBadClusterTimeEnergy);
2240 
2241  fhBadClusterPairDiffTimeE = new TH2F("hBadClusterPairDiffTimeE","cluster pair time difference (bad - good) vs E from bad cluster",nptbins,ptmin,ptmax, tdbins,tdmin,tdmax);
2242  fhBadClusterPairDiffTimeE->SetXTitle("#it{E}_{bad cluster} (GeV)");
2243  fhBadClusterPairDiffTimeE->SetYTitle("#Delta #it{t} (ns)");
2244  outputContainer->Add(fhBadClusterPairDiffTimeE);
2245 
2246  fhBadClusterMaxCellECross = new TH2F ("hBadClusterMaxCellECross","1 - #it{E}_{+} around max energy cell / max energy cell vs cluster energy, bad clusters",
2247  nptbins,ptmin,ptmax, 400,-1,1.);
2248  fhBadClusterMaxCellECross->SetXTitle("#it{E}_{cluster} (GeV) ");
2249  fhBadClusterMaxCellECross->SetYTitle("1- #it{E}_{cross}/#it{E}_{cell max}");
2250  outputContainer->Add(fhBadClusterMaxCellECross);
2251 
2252  if(fFillAllCellTimeHisto)
2253  {
2254  fhBadCellTimeSpreadRespectToCellMax = new TH2F ("hBadCellTimeSpreadRespectToCellMax","#it{t}_{cell max}-#it{t}_{cell i} from bad cluster", nptbins,ptmin,ptmax, tdbins,tdmin,tdmax);
2255  fhBadCellTimeSpreadRespectToCellMax->SetXTitle("#it{E} (GeV)");
2256  fhBadCellTimeSpreadRespectToCellMax->SetYTitle("#Delta #it{t}_{cell max - i} (ns)");
2257  outputContainer->Add(fhBadCellTimeSpreadRespectToCellMax);
2258 
2259  fhBadClusterMaxCellDiffAverageTime = new TH2F ("hBadClusterMaxCellDiffAverageTime","#it{t}_{cell max}-#it{t}_{average} from bad cluster", nptbins,ptmin,ptmax, tdbins,tdmin,tdmax);
2260  fhBadClusterMaxCellDiffAverageTime->SetXTitle("#it{E} (GeV)");
2261  fhBadClusterMaxCellDiffAverageTime->SetYTitle("#Delta #it{t}_{cell max - average} (ns)");
2262  outputContainer->Add(fhBadClusterMaxCellDiffAverageTime);
2263 
2264  fhBadClusterMaxCellDiffWeightedTime = new TH2F ("hBadClusterMaxCellDiffWeightedTime","#it{t}_{cell max}-#it{t}_{weighted} from bad cluster", nptbins,ptmin,ptmax, tdbins,tdmin,tdmax);
2265  fhBadClusterMaxCellDiffWeightedTime->SetXTitle("#it{E} (GeV)");
2266  fhBadClusterMaxCellDiffWeightedTime->SetYTitle("#Delta #it{t}_{cell max - weighted} (ns)");
2267  outputContainer->Add(fhBadClusterMaxCellDiffWeightedTime);
2268 
2269  }
2270 
2271  }
2272 
2273  if(fStudyExotic)
2274  {
2275  fhExoL0ECross = new TH2F("hExoL0_ECross",
2276  "#lambda^{2}_{0} vs 1-#it{E}_{+}/#it{E}_{max} for E > 5 GeV",
2277  400,0,1,ssbins,ssmin,ssmax);
2278  fhExoL0ECross ->SetXTitle("1-#it{E}_{+}/#it{E}_{cell max}");
2279  fhExoL0ECross ->SetYTitle("#lambda^{2}_{0}");
2280  outputContainer->Add(fhExoL0ECross) ;
2281 
2282  fhExoL1ECross = new TH2F("hExoL1_ECross",
2283  "#lambda^{2}_{1} vs 1-#it{E}_{+}/#it{E}_{max} for E > 5 GeV",
2284  400,0,1,ssbins,ssmin,ssmax);
2285  fhExoL1ECross ->SetXTitle("1-#it{E}_{+}/#it{E}_{cell max}");
2286  fhExoL1ECross ->SetYTitle("#lambda^{2}_{1}");
2287  outputContainer->Add(fhExoL1ECross) ;
2288 
2289  for(Int_t ie = 0; ie <fExoNECrossCuts; ie++)
2290  {
2291 
2292  fhExoDTime[ie] = new TH2F(Form("hExoDTime_ECross%d",ie),
2293  Form("#Delta time = t_{max}-t_{cells} vs #it{E}_{cluster} for exotic, 1-#it{E}_{+}/#it{E}_{max} < %2.2f",fExoECrossCuts[ie]),
2294  nptbins,ptmin,ptmax,tdbins,tdmin,tdmax);
2295  fhExoDTime[ie] ->SetYTitle("#Delta #it{t} (ns)");
2296  fhExoDTime[ie] ->SetXTitle("#it{E} (GeV)");
2297  outputContainer->Add(fhExoDTime[ie]) ;
2298 
2299  for(Int_t idt = 0; idt < fExoNDTimeCuts; idt++)
2300  {
2301  fhExoNCell[ie][idt] = new TH2F(Form("hExoNCell_ECross%d_DT%d",ie,idt),
2302  Form("N cells per cluster vs E cluster, 1-#it{E}_{+}/#it{E}_{max} < %2.2f, #Delta t < %2.0f",fExoECrossCuts[ie],fExoDTimeCuts[idt]),
2303  nptbins,ptmin,ptmax,nceclbins,nceclmin,nceclmax);
2304  fhExoNCell[ie][idt] ->SetYTitle("#it{n}_cells");
2305  fhExoNCell[ie][idt] ->SetXTitle("#it{E} (GeV)");
2306  outputContainer->Add(fhExoNCell[ie][idt]) ;
2307 
2308  fhExoL0 [ie][idt] = new TH2F(Form("hExoL0_ECross%d_DT%d",ie,idt),
2309  Form("#lambda^{2}_{0} vs E cluster for exotic, 1-#it{E}_{+}/#it{E}_{max} < %2.2f, #Delta t = %2.0f",fExoECrossCuts[ie],fExoDTimeCuts[idt]),
2310  nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
2311  fhExoL0 [ie][idt] ->SetYTitle("#lambda^{2}_{0}");
2312  fhExoL0 [ie][idt] ->SetXTitle("#it{E} (GeV)");
2313  outputContainer->Add(fhExoL0[ie][idt]) ;
2314 
2315  fhExoL1 [ie][idt] = new TH2F(Form("hExoL1_ECross%d_DT%d",ie,idt),
2316  Form("#lambda^{2}_{1} vs E cluster for exotic, 1-#it{E}_{+}/#it{E}_{max} < %2.2f, #Delta t = %2.0f",fExoECrossCuts[ie],fExoDTimeCuts[idt]),
2317  nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
2318  fhExoL1 [ie][idt] ->SetYTitle("#lambda^{2}_{1}");
2319  fhExoL1 [ie][idt] ->SetXTitle("#it{E} (GeV)");
2320  outputContainer->Add(fhExoL1[ie][idt]) ;
2321 
2322  fhExoECross[ie][idt] = new TH2F(Form("hExoECross_ECross%d_DT%d",ie,idt),
2323  Form("#it{E} cross for cells vs E cell, 1-#it{E}_{+}/#it{E}_{max} < %2.2f, #Delta t < %2.0f",fExoECrossCuts[ie],fExoDTimeCuts[idt]),
2324  nptbins,ptmin,ptmax,400,0,1);
2325  fhExoECross[ie][idt] ->SetYTitle("1-#it{E}_{+}/#it{E}_{cell max}");
2326  fhExoECross[ie][idt] ->SetXTitle("#it{E}_{cell} (GeV)");
2327  outputContainer->Add(fhExoECross[ie][idt]) ;
2328 
2329  fhExoTime [ie][idt] = new TH2F(Form("hExoTime_ECross%d_DT%d",ie,idt),
2330  Form("Time of cluster (max cell) vs E cluster for exotic, 1-#it{E}_{+}/#it{E}_{max} < %2.2f, #Delta t = %2.0f",fExoECrossCuts[ie],fExoDTimeCuts[idt]),
2331  nptbins,ptmin,ptmax,ntimebins,timemin,timemax);
2332  fhExoTime [ie][idt] ->SetYTitle("#it{t}_{max} (ns)");
2333  fhExoTime [ie][idt] ->SetXTitle("#it{E} (GeV)");
2334  outputContainer->Add(fhExoTime[ie][idt]) ;
2335 
2336  fhExoL0NCell[ie][idt] = new TH2F(Form("hExoL0_NCell%d_DT%d",ie,idt),
2337  Form("#lambda^{2}_{0} vs N cells per clusters for E > 5 GeV, 1-#it{E}_{+}/#it{E}_{max} < %2.2f, #Delta t = %2.0f",fExoECrossCuts[ie],fExoDTimeCuts[idt]),
2338  nptbins,ptmin,ptmax,ntimebins,timemin,timemax);
2339  fhExoL0NCell[ie][idt] ->SetYTitle("#it{n}_{cells}");
2340  fhExoL0NCell[ie][idt] ->SetXTitle("#lambda^{2}_{0}");
2341  outputContainer->Add(fhExoL0NCell[ie][idt]) ;
2342 
2343  fhExoL1NCell[ie][idt] = new TH2F(Form("hExoL1_NCell%d_DT%d",ie,idt),
2344  Form("#lambda^{2}_{1} vs N cells per clusters for E > 5 GeV, 1-#it{E}_{+}/#it{E}_{max} < %2.2f, #Delta t = %2.0f",fExoECrossCuts[ie],fExoDTimeCuts[idt]),
2345  nptbins,ptmin,ptmax,ntimebins,timemin,timemax);
2346  fhExoL1NCell[ie][idt] ->SetYTitle("#it{n}_{cells}");
2347  fhExoL1NCell[ie][idt] ->SetXTitle("#lambda^{2}_{1}");
2348  outputContainer->Add(fhExoL1NCell[ie][idt]) ;
2349  }
2350  }
2351  }
2352 
2353  // Cluster size in terms of cells
2354 
2355  if(fStudyClustersAsymmetry)
2356  {
2357  fhDeltaIEtaDeltaIPhiE0[0] = new TH2F ("hDeltaIEtaDeltaIPhiE0"," Cluster size in columns vs rows for E < 2 GeV, #it{n}_{cells} > 3",
2358  50,0,50,50,0,50);
2359  fhDeltaIEtaDeltaIPhiE0[0]->SetXTitle("#Delta Column");
2360  fhDeltaIEtaDeltaIPhiE0[0]->SetYTitle("#Delta Row");
2361  outputContainer->Add(fhDeltaIEtaDeltaIPhiE0[0]);
2362 
2363  fhDeltaIEtaDeltaIPhiE2[0] = new TH2F ("hDeltaIEtaDeltaIPhiE2"," Cluster size in columns vs rows for 2 <E < 6 GeV, #it{n}_{cells} > 3",
2364  50,0,50,50,0,50);
2365  fhDeltaIEtaDeltaIPhiE2[0]->SetXTitle("#Delta Column");
2366  fhDeltaIEtaDeltaIPhiE2[0]->SetYTitle("#Delta Row");
2367  outputContainer->Add(fhDeltaIEtaDeltaIPhiE2[0]);
2368 
2369  fhDeltaIEtaDeltaIPhiE6[0] = new TH2F ("hDeltaIEtaDeltaIPhiE6"," Cluster size in columns vs rows for E > 6 GeV, #it{n}_{cells} > 3",
2370  50,0,50,50,0,50);
2371  fhDeltaIEtaDeltaIPhiE6[0]->SetXTitle("#Delta Column");
2372  fhDeltaIEtaDeltaIPhiE6[0]->SetYTitle("#Delta Row");
2373  outputContainer->Add(fhDeltaIEtaDeltaIPhiE6[0]);
2374 
2375  fhDeltaIA[0] = new TH2F ("hDeltaIA"," Cluster *asymmetry* in cell units vs E",
2376  nptbins,ptmin,ptmax,21,-1.05,1.05);
2377  fhDeltaIA[0]->SetXTitle("#it{E}_{cluster}");
2378  fhDeltaIA[0]->SetYTitle("#it{A}_{cell in cluster}");
2379  outputContainer->Add(fhDeltaIA[0]);
2380 
2381  fhDeltaIAL0[0] = new TH2F ("hDeltaIAL0"," Cluster *asymmetry* in cell units vs #lambda^{2}_{0}",
2382  ssbins,ssmin,ssmax,21,-1.05,1.05);
2383  fhDeltaIAL0[0]->SetXTitle("#lambda^{2}_{0}");
2384  fhDeltaIAL0[0]->SetYTitle("#it{A}_{cell in cluster}");
2385  outputContainer->Add(fhDeltaIAL0[0]);
2386 
2387  fhDeltaIAL1[0] = new TH2F ("hDeltaIAL1"," Cluster *asymmetry* in cell units vs #lambda^{2}_{1}",
2388  ssbins,ssmin,ssmax,21,-1.05,1.05);
2389  fhDeltaIAL1[0]->SetXTitle("#lambda^{2}_{1}");
2390  fhDeltaIAL1[0]->SetYTitle("#it{A}_{cell in cluster}");
2391  outputContainer->Add(fhDeltaIAL1[0]);
2392 
2393  fhDeltaIANCells[0] = new TH2F ("hDeltaIANCells"," Cluster *asymmetry* in cell units vs N cells in cluster",
2394  nceclbins,nceclmin,nceclmax,21,-1.05,1.05);
2395  fhDeltaIANCells[0]->SetXTitle("#it{n}_{cell in cluster}");
2396  fhDeltaIANCells[0]->SetYTitle("#it{A}_{cell in cluster}");
2397  outputContainer->Add(fhDeltaIANCells[0]);
2398 
2399 
2400  fhDeltaIEtaDeltaIPhiE0[1] = new TH2F ("hDeltaIEtaDeltaIPhiE0Charged"," Cluster size in columns vs rows for E < 2 GeV, #it{n}_{cells} > 3, matched with track",
2401  50,0,50,50,0,50);
2402  fhDeltaIEtaDeltaIPhiE0[1]->SetXTitle("#Delta Column");
2403  fhDeltaIEtaDeltaIPhiE0[1]->SetYTitle("#Delta Row");
2404  outputContainer->Add(fhDeltaIEtaDeltaIPhiE0[1]);
2405 
2406  fhDeltaIEtaDeltaIPhiE2[1] = new TH2F ("hDeltaIEtaDeltaIPhiE2Charged"," Cluster size in columns vs rows for 2 <E < 6 GeV, #it{n}_{cells} > 3, matched with track",
2407  50,0,50,50,0,50);
2408  fhDeltaIEtaDeltaIPhiE2[1]->SetXTitle("#Delta Column");
2409  fhDeltaIEtaDeltaIPhiE2[1]->SetYTitle("#Delta Row");
2410  outputContainer->Add(fhDeltaIEtaDeltaIPhiE2[1]);
2411 
2412  fhDeltaIEtaDeltaIPhiE6[1] = new TH2F ("hDeltaIEtaDeltaIPhiE6Charged"," Cluster size in columns vs rows for E > 6 GeV, #it{n}_{cells} > 3, matched with track",
2413  50,0,50,50,0,50);
2414  fhDeltaIEtaDeltaIPhiE6[1]->SetXTitle("#Delta Column");
2415  fhDeltaIEtaDeltaIPhiE6[1]->SetYTitle("#Delta Row");
2416  outputContainer->Add(fhDeltaIEtaDeltaIPhiE6[1]);
2417 
2418  fhDeltaIA[1] = new TH2F ("hDeltaIACharged"," Cluster *asymmetry* in cell units vs E, matched with track",
2419  nptbins,ptmin,ptmax,21,-1.05,1.05);
2420  fhDeltaIA[1]->SetXTitle("#it{E}_{cluster}");
2421  fhDeltaIA[1]->SetYTitle("#it{A}_{cell in cluster}");
2422  outputContainer->Add(fhDeltaIA[1]);
2423 
2424  fhDeltaIAL0[1] = new TH2F ("hDeltaIAL0Charged"," Cluster *asymmetry* in cell units vs #lambda^{2}_{0}, matched with track",
2425  ssbins,ssmin,ssmax,21,-1.05,1.05);
2426  fhDeltaIAL0[1]->SetXTitle("#lambda^{2}_{0}");
2427  fhDeltaIAL0[1]->SetYTitle("#it{A}_{cell in cluster}");
2428  outputContainer->Add(fhDeltaIAL0[1]);
2429 
2430  fhDeltaIAL1[1] = new TH2F ("hDeltaIAL1Charged"," Cluster *asymmetry* in cell units vs #lambda^{2}_{1}, matched with track",
2431  ssbins,ssmin,ssmax,21,-1.05,1.05);
2432  fhDeltaIAL1[1]->SetXTitle("#lambda^{2}_{1}");
2433  fhDeltaIAL1[1]->SetYTitle("#it{A}_{cell in cluster}");
2434  outputContainer->Add(fhDeltaIAL1[1]);
2435 
2436  fhDeltaIANCells[1] = new TH2F ("hDeltaIANCellsCharged"," Cluster *asymmetry* in cell units vs N cells in cluster, matched with track",
2437  nceclbins,nceclmin,nceclmax,21,-1.05,1.05);
2438  fhDeltaIANCells[1]->SetXTitle("#it{n}_{cell in cluster}");
2439  fhDeltaIANCells[1]->SetYTitle("#it{A}_{cell in cluster}");
2440  outputContainer->Add(fhDeltaIANCells[1]);
2441 
2442  if(IsDataMC()){
2443  TString particle[]={"Photon","Electron","Conversion","Hadron"};
2444  for (Int_t iPart = 0; iPart < 4; iPart++) {
2445 
2446  fhDeltaIAMC[iPart] = new TH2F (Form("hDeltaIA_MC%s",particle[iPart].Data()),Form(" Cluster *asymmetry* in cell units vs E, from %s",particle[iPart].Data()),
2447  nptbins,ptmin,ptmax,21,-1.05,1.05);
2448  fhDeltaIAMC[iPart]->SetXTitle("#it{E}_{cluster}");
2449  fhDeltaIAMC[iPart]->SetYTitle("#it{A}_{cell in cluster}");
2450  outputContainer->Add(fhDeltaIAMC[iPart]);
2451  }
2452  }
2453 
2454  if(fStudyBadClusters)
2455  {
2456  fhBadClusterDeltaIEtaDeltaIPhiE0 = new TH2F ("hBadClusterDeltaIEtaDeltaIPhiE0"," Cluster size in columns vs rows for E < 2 GeV, #it{n}_{cells} > 3",
2457  50,0,50,50,0,50);
2458  fhBadClusterDeltaIEtaDeltaIPhiE0->SetXTitle("#Delta Column");
2459  fhBadClusterDeltaIEtaDeltaIPhiE0->SetYTitle("#Delta Row");
2460  outputContainer->Add(fhBadClusterDeltaIEtaDeltaIPhiE0);
2461 
2462  fhBadClusterDeltaIEtaDeltaIPhiE2 = new TH2F ("hBadClusterDeltaIEtaDeltaIPhiE2"," Cluster size in columns vs rows for 2 <E < 6 GeV, #it{n}_{cells} > 3",
2463  50,0,50,50,0,50);
2464  fhBadClusterDeltaIEtaDeltaIPhiE2->SetXTitle("#Delta Column");
2465  fhBadClusterDeltaIEtaDeltaIPhiE2->SetYTitle("#Delta Row");
2466  outputContainer->Add(fhBadClusterDeltaIEtaDeltaIPhiE2);
2467 
2468  fhBadClusterDeltaIEtaDeltaIPhiE6 = new TH2F ("hBadClusterDeltaIEtaDeltaIPhiE6"," Cluster size in columns vs rows for E > 6 GeV, #it{n}_{cells} > 3",
2469  50,0,50,50,0,50);
2470  fhBadClusterDeltaIEtaDeltaIPhiE6->SetXTitle("#Delta Column");
2471  fhBadClusterDeltaIEtaDeltaIPhiE6->SetYTitle("#Delta Row");
2472  outputContainer->Add(fhBadClusterDeltaIEtaDeltaIPhiE6);
2473 
2474  fhBadClusterDeltaIA = new TH2F ("hBadClusterDeltaIA"," Cluster *asymmetry* in cell units vs E",
2475  nptbins,ptmin,ptmax,21,-1.05,1.05);
2476  fhBadClusterDeltaIA->SetXTitle("#it{E}_{cluster}");
2477  fhBadClusterDeltaIA->SetYTitle("#it{A}_{cell in cluster}");
2478  outputContainer->Add(fhBadClusterDeltaIA);
2479  }
2480  }
2481 
2482  if(fStudyWeight)
2483  {
2484  fhECellClusterRatio = new TH2F ("hECellClusterRatio"," cell energy / cluster energy vs cluster energy",
2485  nptbins,ptmin,ptmax, 100,0,1.);
2486  fhECellClusterRatio->SetXTitle("#it{E}_{cluster} (GeV) ");
2487  fhECellClusterRatio->SetYTitle("#it{E}_{cell i}/#it{E}_{cluster}");
2488  outputContainer->Add(fhECellClusterRatio);
2489 
2490  fhECellClusterLogRatio = new TH2F ("hECellClusterLogRatio"," Log(cell energy / cluster energy) vs cluster energy",
2491  nptbins,ptmin,ptmax, 100,-10,0);
2492  fhECellClusterLogRatio->SetXTitle("#it{E}_{cluster} (GeV) ");
2493  fhECellClusterLogRatio->SetYTitle("Log(#it{E}_{cell i}/#it{E}_{cluster})");
2494  outputContainer->Add(fhECellClusterLogRatio);
2495 
2496  fhEMaxCellClusterRatio = new TH2F ("hEMaxCellClusterRatio"," max cell energy / cluster energy vs cluster energy",
2497  nptbins,ptmin,ptmax, 100,0,1.);
2498  fhEMaxCellClusterRatio->SetXTitle("#it{E}_{cluster} (GeV) ");
2499  fhEMaxCellClusterRatio->SetYTitle("#it{E}_{max cell}/#it{E}_{cluster}");
2500  outputContainer->Add(fhEMaxCellClusterRatio);
2501 
2502  fhEMaxCellClusterLogRatio = new TH2F ("hEMaxCellClusterLogRatio"," Log(max cell energy / cluster energy) vs cluster energy",
2503  nptbins,ptmin,ptmax, 100,-10,0);
2504  fhEMaxCellClusterLogRatio->SetXTitle("#it{E}_{cluster} (GeV) ");
2505  fhEMaxCellClusterLogRatio->SetYTitle("Log (#it{E}_{max cell}/#it{E}_{cluster})");
2506  outputContainer->Add(fhEMaxCellClusterLogRatio);
2507 
2508  fhECellTotalRatio = new TH2F ("hECellTotalRatio"," cell energy / sum all energy vs all energy",
2509  nptbins*2,ptmin,ptmax*2, 100,0,1.);
2510  fhECellTotalRatio->SetXTitle("#it{E}_{total} (GeV) ");
2511  fhECellTotalRatio->SetYTitle("#it{E}_{cell i}/#it{E}_{total}");
2512  outputContainer->Add(fhECellTotalRatio);
2513 
2514  fhECellTotalLogRatio = new TH2F ("hECellTotalLogRatio"," Log(cell energy / sum all energy) vs all energy",
2515  nptbins*2,ptmin,ptmax*2, 100,-10,0);
2516  fhECellTotalLogRatio->SetXTitle("#it{E}_{total} (GeV) ");
2517  fhECellTotalLogRatio->SetYTitle("Log(#it{E}_{cell i}/#it{E}_{total})");
2518  outputContainer->Add(fhECellTotalLogRatio);
2519 
2520  fhECellTotalRatioMod = new TH2F*[fNModules];
2521  fhECellTotalLogRatioMod = new TH2F*[fNModules];
2522 
2523  for(Int_t imod = 0; imod < fNModules; imod++)
2524  {
2525  fhECellTotalRatioMod[imod] = new TH2F (Form("hECellTotalRatio_Mod%d",imod),
2526  Form("#cell energy / sum all energy vs all energy in Module %d",imod),
2527  nptbins*2,ptmin,ptmax*2, 100,0,1.);
2528  fhECellTotalRatioMod[imod]->SetXTitle("#it{E} (GeV)");
2529  fhECellTotalRatioMod[imod]->SetYTitle("#it{n}_{cells}");
2530  outputContainer->Add(fhECellTotalRatioMod[imod]);
2531 
2532  fhECellTotalLogRatioMod[imod] = new TH2F (Form("hECellTotalLogRatio_Mod%d",imod),
2533  Form("Log(cell energy / sum all energy) vs all energy in Module %d",imod),
2534  nptbins*2,ptmin,ptmax*2, 100,-10,0);
2535  fhECellTotalLogRatioMod[imod]->SetXTitle("#it{E} (GeV)");
2536  fhECellTotalLogRatioMod[imod]->SetYTitle("#it{n}_{cells}");
2537  outputContainer->Add(fhECellTotalLogRatioMod[imod]);
2538  }
2539 
2540  // To be done properly with setters and data members ...
2541  Float_t cellEmin [] = {0.05,0.1,0.15,0.2};
2542  Float_t cellTmin [] = {50.,100.,1000000.};
2543 
2544  for(Int_t iw = 0; iw < 12; iw++)
2545  {
2546  Float_t w0 = 4+0.05*iw; // 3+0.25*iw;
2547  for(Int_t iEmin = 0; iEmin < 4; iEmin++)
2548  {
2549  for(Int_t iTmin = 0; iTmin < 3; iTmin++)
2550  {
2551  fhLambda0ForW0AndCellCuts[iw][iEmin][iTmin] = new TH2F (Form("hLambda0ForW0%d_CellEMin%d_TimeMax%d",iw,iEmin,iTmin),
2552  Form("#lambda^{2}_{0} vs E, w0=%1.2f, cell E>%2.2f MeV, |t|<%2.0f ns",
2553  w0, cellEmin[iEmin], cellTmin[iTmin]),
2554  nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
2555  fhLambda0ForW0AndCellCuts[iw][iEmin][iTmin]->SetXTitle("#it{E}_{cluster}");
2556  fhLambda0ForW0AndCellCuts[iw][iEmin][iTmin]->SetYTitle("#lambda^{2}_{0}");
2557  outputContainer->Add(fhLambda0ForW0AndCellCuts[iw][iEmin][iTmin]);
2558 
2559 
2560 // fhLambda1ForW0AndCellCuts[iw][iEmin][iTmin] = new TH2F (Form("hLambda1ForW0%d_CellEMin%d_TimeMax%d",iw,iEmin,iTmin),
2561 // Form("#lambda^{2}_{1} vs E, w0=%1.2f, cell E>%2.2f MeV, |t|<%2.0f ns"",
2562 // w0, cellEmin[iEmin], cellTmin[iTmin]),
2563 // nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
2564 // fhLambda1ForW0AndCellCuts[iw][iEmin][iTmin]->SetXTitle("#it{E}_{cluster}");
2565 // fhLambda1ForW0AndCellCuts[iw][iEmin][iTmin]->SetYTitle("#lambda^{2}_{1}");
2566 // outputContainer->Add(fhLambda1ForW0AndCellCuts[iw][iEmin][iTmin]);
2567 
2568  fhLambda0ForW0AndCellCutsEta0[iw][iEmin][iTmin] = new TH2F (Form("hLambda0ForW0%d_CellEMin%d_TimeMax%d_Eta0",iw,iEmin,iTmin),
2569  Form("#lambda^{2}_{0} vs E, w0=%1.2f, cell E>%2.2f MeV, |t|<%2.0f ns, |#eta| < 0.15",
2570  w0, cellEmin[iEmin], cellTmin[iTmin]),
2571  nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
2572  fhLambda0ForW0AndCellCutsEta0[iw][iEmin][iTmin]->SetXTitle("#it{E}_{cluster}");
2573  fhLambda0ForW0AndCellCutsEta0[iw][iEmin][iTmin]->SetYTitle("#lambda^{2}_{0}");
2574  outputContainer->Add(fhLambda0ForW0AndCellCutsEta0[iw][iEmin][iTmin]);
2575  }
2576  }
2577 
2578  if(IsDataMC())
2579  {
2580  TString mcnames[] = {"Photon", "Electron","Conversion","Pi0","Hadron"};
2581  for(Int_t imc = 0; imc < 5; imc++)
2582  {
2583  fhLambda0ForW0MC[iw][imc] = new TH2F (Form("hLambda0ForW0%d_MC%s",iw,mcnames[imc].Data()),
2584  Form("shower shape, #lambda^{2}_{0} vs E, w0 = %1.1f, for MC %s",w0,mcnames[imc].Data()),
2585  nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
2586  fhLambda0ForW0MC[iw][imc]->SetXTitle("#it{E}_{cluster}");
2587  fhLambda0ForW0MC[iw][imc]->SetYTitle("#lambda^{2}_{0}");
2588  outputContainer->Add(fhLambda0ForW0MC[iw][imc]);
2589 
2590 // fhLambda1ForW0MC[iw][imc] = new TH2F (Form("hLambda1ForW0%d_MC%s",iw,mcnames[imc].Data()),
2591 // Form("shower shape, #lambda^{2}_{1} vs E, w0 = %1.1f, for MC %s",w0,mcnames[imc].Data()),
2592 // nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
2593 // fhLambda1ForW0MC[iw][imc]->SetXTitle("#it{E}_{cluster}");
2594 // fhLambda1ForW0MC[iw][imc]->SetYTitle("#lambda^{2}_{1}");
2595 // outputContainer->Add(fhLambda1ForW0MC[iw][imc]);
2596  }
2597  }
2598  }
2599  }
2600 
2601  // Track Matching
2602  if(fFillAllTMHisto)
2603  {
2604  Int_t nresetabins = GetHistogramRanges()->GetHistoTrackResidualEtaBins();
2605  Float_t resetamax = GetHistogramRanges()->GetHistoTrackResidualEtaMax();
2606  Float_t resetamin = GetHistogramRanges()->GetHistoTrackResidualEtaMin();
2607  Int_t nresphibins = GetHistogramRanges()->GetHistoTrackResidualPhiBins();
2608  Float_t resphimax = GetHistogramRanges()->GetHistoTrackResidualPhiMax();
2609  Float_t resphimin = GetHistogramRanges()->GetHistoTrackResidualPhiMin();
2610 
2611  fhTrackMatchedDEtaNeg = new TH2F("hTrackMatchedDEtaNeg","d#eta of cluster-negative track vs cluster energy",
2612  nptbins,ptmin,ptmax,nresetabins,resetamin,resetamax);
2613  fhTrackMatchedDEtaNeg->SetYTitle("d#eta");
2614  fhTrackMatchedDEtaNeg->SetXTitle("#it{E}_{cluster} (GeV)");
2615 
2616  fhTrackMatchedDPhiNeg = new TH2F("hTrackMatchedDPhiNeg","d#phi of cluster-negative track vs cluster energy",
2617  nptbins,ptmin,ptmax,nresphibins,resphimin,resphimax);
2618  fhTrackMatchedDPhiNeg->SetYTitle("d#phi (rad)");
2619  fhTrackMatchedDPhiNeg->SetXTitle("#it{E}_{cluster} (GeV)");
2620 
2621  fhTrackMatchedDEtaDPhiNeg = new TH2F("hTrackMatchedDEtaDPhiNeg","d#eta vs d#phi of cluster- negative track vs cluster energy",
2622  nresetabins,resetamin,resetamax,nresphibins,resphimin,resphimax);
2623  fhTrackMatchedDEtaDPhiNeg->SetYTitle("d#phi (rad)");
2624  fhTrackMatchedDEtaDPhiNeg->SetXTitle("d#eta");
2625 
2626  fhTrackMatchedDEtaPos = new TH2F("hTrackMatchedDEtaPos","d#eta of cluster-positive track vs cluster energy",
2627  nptbins,ptmin,ptmax,nresetabins,resetamin,resetamax);
2628  fhTrackMatchedDEtaPos->SetYTitle("d#eta");
2629  fhTrackMatchedDEtaPos->SetXTitle("#it{E}_{cluster} (GeV)");
2630 
2631  fhTrackMatchedDPhiPos = new TH2F("hTrackMatchedDPhiPos","d#phi of cluster-positive track vs cluster energy",
2632  nptbins,ptmin,ptmax,nresphibins,resphimin,resphimax);
2633  fhTrackMatchedDPhiPos->SetYTitle("d#phi (rad)");
2634  fhTrackMatchedDPhiPos->SetXTitle("#it{E}_{cluster} (GeV)");
2635 
2636  fhTrackMatchedDEtaDPhiPos = new TH2F("hTrackMatchedDEtaDPhiPos","d#eta vs d#phi of cluster-positive track vs cluster energy",
2637  nresetabins,resetamin,resetamax,nresphibins,resphimin,resphimax);
2638  fhTrackMatchedDEtaDPhiPos->SetYTitle("d#phi (rad)");
2639  fhTrackMatchedDEtaDPhiPos->SetXTitle("d#eta");
2640 
2641 
2642  fhTrackMatchedDEtaNegMod = new TH2F("hTrackMatchedDEtaNegPerModule","d#eta of cluster-negative track vs module, E > 0.5 GeV",
2643  nresetabins,resetamin,resetamax,fNModules,0,fNModules);
2644  fhTrackMatchedDEtaNegMod->SetXTitle("d#eta");
2645  fhTrackMatchedDEtaNegMod->SetYTitle("Module");
2646 
2647  fhTrackMatchedDPhiNegMod = new TH2F("hTrackMatchedDPhiNegPerModule","d#phi of cluster-negative track vs module, E > 0.5 GeV",
2648  nresetabins,resetamin,resetamax,fNModules,0,fNModules);
2649  fhTrackMatchedDPhiNegMod->SetXTitle("d#phi (rad)");
2650  fhTrackMatchedDPhiNegMod->SetYTitle("Module");
2651 
2652  fhTrackMatchedDEtaPosMod = new TH2F("hTrackMatchedDEtaPosPerModule","d#eta of cluster-positive track vs module, E > 0.5 GeV",
2653  nresetabins,resetamin,resetamax,fNModules,0,fNModules);
2654  fhTrackMatchedDEtaPosMod->SetXTitle("d#eta");
2655  fhTrackMatchedDEtaPosMod->SetYTitle("Module");
2656 
2657  fhTrackMatchedDPhiPosMod = new TH2F("hTrackMatchedDPhiPosPerModule","d#phi of cluster-positive track vs module, E > 0.5 GeV",
2658  nresetabins,resetamin,resetamax,fNModules,0,fNModules);
2659  fhTrackMatchedDPhiPosMod->SetXTitle("d#phi (rad)");
2660  fhTrackMatchedDPhiPosMod->SetYTitle("Module");
2661 
2662  outputContainer->Add(fhTrackMatchedDEtaNeg) ;
2663  outputContainer->Add(fhTrackMatchedDPhiNeg) ;
2664  outputContainer->Add(fhTrackMatchedDEtaPos) ;
2665  outputContainer->Add(fhTrackMatchedDPhiPos) ;
2666  outputContainer->Add(fhTrackMatchedDEtaDPhiNeg) ;
2667  outputContainer->Add(fhTrackMatchedDEtaDPhiPos) ;
2668 
2669  outputContainer->Add(fhTrackMatchedDEtaNegMod) ;
2670  outputContainer->Add(fhTrackMatchedDPhiNegMod) ;
2671  outputContainer->Add(fhTrackMatchedDEtaPosMod) ;
2672  outputContainer->Add(fhTrackMatchedDPhiPosMod) ;
2673 
2674  fhECharged = new TH1F ("hECharged","#it{E} reconstructed clusters, matched with track", nptbins,ptmin,ptmax);
2675  fhECharged->SetXTitle("#it{E} (GeV)");
2676  outputContainer->Add(fhECharged);
2677 
2678  fhPtCharged = new TH1F ("hPtCharged","#it{p}_{T} reconstructed clusters, matched with track", nptbins,ptmin,ptmax);
2679  fhPtCharged->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2680  outputContainer->Add(fhPtCharged);
2681 
2682  fhPhiCharged = new TH1F ("hPhiCharged","#phi reconstructed clusters, matched with track",nphibins,phimin,phimax);
2683  fhPhiCharged->SetXTitle("#phi (rad)");
2684  outputContainer->Add(fhPhiCharged);
2685 
2686  fhEtaCharged = new TH1F ("hEtaCharged","#eta reconstructed clusters, matched with track",netabins,etamin,etamax);
2687  fhEtaCharged->SetXTitle("#eta ");
2688  outputContainer->Add(fhEtaCharged);
2689 
2690  if(fFillAllTH3)
2691  {
2692  fhEtaPhiECharged = new TH3F ("hEtaPhiECharged","#eta vs #phi, reconstructed clusters, matched with track",
2693  netabins,etamin,etamax,nphibins,phimin,phimax,nptbins,ptmin,ptmax);
2694  fhEtaPhiECharged->SetXTitle("#eta ");
2695  fhEtaPhiECharged->SetYTitle("#phi ");
2696  fhEtaPhiECharged->SetZTitle("#it{E} (GeV) ");
2697  outputContainer->Add(fhEtaPhiECharged);
2698  }
2699  else
2700  {
2701  fhEtaPhiCharged = new TH2F ("hEtaPhiCharged","#eta vs #phi for #it{E} > 0.5 GeV, reconstructed clusters, with matched track",
2702  netabins,etamin,etamax,nphibins,phimin,phimax);
2703  fhEtaPhiCharged->SetXTitle("#eta ");
2704  fhEtaPhiCharged->SetYTitle("#phi (rad)");
2705  outputContainer->Add(fhEtaPhiCharged);
2706  }
2707 
2708  fh1EOverP = new TH2F("h1EOverP","TRACK matches #it{E}/#it{p}",nptbins,ptmin,ptmax, nPoverEbins,eOverPmin,eOverPmax);
2709  fh1EOverP->SetYTitle("#it{E}/#it{p}");
2710  fh1EOverP->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2711  outputContainer->Add(fh1EOverP);
2712 
2713  fh2dR = new TH2F("h2dR","TRACK matches #Delta #it{R}",nptbins,ptmin,ptmax,ndRbins,dRmin,dRmax);
2714  fh2dR->SetYTitle("#Delta #it{R} (rad)");
2715  fh2dR->SetXTitle("#it{E} cluster (GeV)");
2716  outputContainer->Add(fh2dR) ;
2717 
2718  fh2MatchdEdx = new TH2F("h2MatchdEdx","#it{dE/dx} vs. #it{p} for all matches",nptbins,ptmin,ptmax,ndedxbins,dedxmin,dedxmax);
2719  fh2MatchdEdx->SetXTitle("p (GeV/#it{c})");
2720  fh2MatchdEdx->SetYTitle("<#it{dE/dx}>");
2721  outputContainer->Add(fh2MatchdEdx);
2722 
2723  fh2EledEdx = new TH2F("h2EledEdx","#it{dE/dx} vs. #it{p} for electrons",nptbins,ptmin,ptmax,ndedxbins,dedxmin,dedxmax);
2724  fh2EledEdx->SetXTitle("p (GeV/#it{c})");
2725  fh2EledEdx->SetYTitle("<#it{dE/dx}>");
2726  outputContainer->Add(fh2EledEdx) ;
2727 
2728  fh1EOverPR02 = new TH2F("h1EOverPR02","TRACK matches #it{E}/#it{p}, all",nptbins,ptmin,ptmax, nPoverEbins,eOverPmin,eOverPmax);
2729  fh1EOverPR02->SetYTitle("#it{E}/#it{p}");
2730  fh1EOverPR02->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2731  outputContainer->Add(fh1EOverPR02);
2732 
2733  fh1EleEOverP = new TH2F("h1EleEOverP","Electron candidates #it{E}/#it{p} (60<#it{dE/dx}<100)",nptbins,ptmin,ptmax, nPoverEbins,eOverPmin,eOverPmax);
2734  fh1EleEOverP->SetYTitle("#it{E}/#it{p}");
2735  fh1EleEOverP->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2736  outputContainer->Add(fh1EleEOverP);
2737 
2738 
2739  // Projections per SM
2740 
2741  fh1EOverPMod = new TH2F("h1EOverP_PerModule","TRACK matches #it{E}/#it{p}, #it{E}_{cl}&#it{p}_{tr}>0.5 Gev/#it{c}",nPoverEbins,eOverPmin,eOverPmax,fNModules,0,fNModules);
2742  fh1EOverPMod->SetXTitle("#it{E}/#it{p}");
2743  fh1EOverPMod->SetYTitle("Module");
2744  outputContainer->Add(fh1EOverPMod);
2745 
2746  fh2dRMod = new TH2F("h2dR_PerModule","TRACK matches #Delta #it{R}, #it{E}_{cl}&#it{p}_{tr}>0.5 Gev/#it{c}",ndRbins,dRmin,dRmax,fNModules,0,fNModules);
2747  fh2dRMod->SetXTitle("#Delta #it{R} (rad)");
2748  fh2dRMod->SetYTitle("Module");
2749  outputContainer->Add(fh2dRMod) ;
2750 
2751  fh2MatchdEdxMod = new TH2F("h2MatchdEdx_PerModule","#it{dE/dx} vs. #it{p} for all matches, #it{E}_{cl}&#it{p}_{tr}>0.5 Gev/#it{c}",ndedxbins,dedxmin,dedxmax,fNModules,0,fNModules);
2752  fh2MatchdEdxMod->SetYTitle("Module");
2753  fh2MatchdEdxMod->SetXTitle("<#it{dE/dx}>");
2754  outputContainer->Add(fh2MatchdEdxMod);
2755 
2756  fh2EledEdxMod = new TH2F("h2EledEdx_PerModule","#it{dE/dx} vs. #it{p} for electrons, #it{E}_{cl}&#it{p}_{tr}>0.5 Gev/#it{c}",ndedxbins,dedxmin,dedxmax,fNModules,0,fNModules);
2757  fh2EledEdxMod->SetYTitle("Module");
2758  fh2EledEdxMod->SetXTitle("<#it{dE/dx}>");
2759  outputContainer->Add(fh2EledEdxMod) ;
2760 
2761  fh1EOverPR02Mod = new TH2F("h1EOverPR02_PerModule","TRACK matches #it{E}/#it{p}, all, #it{E}_{cl}&#it{p}_{tr}>0.5 Gev/#it{c}",nPoverEbins,eOverPmin,eOverPmax,fNModules,0,fNModules);
2762  fh1EOverPR02Mod->SetXTitle("#it{E}/#it{p}");
2763  fh1EOverPR02Mod->SetYTitle("Module");
2764  outputContainer->Add(fh1EOverPR02Mod);
2765 
2766  fh1EleEOverPMod = new TH2F("h1EleEOverP_PerModule","Electron candidates #it{E}/#it{p} (60<#it{dE/dx}<100), #it{E}_{cl}&#it{p}_{tr}>0.5 Gev/#it{c}",nPoverEbins,eOverPmin,eOverPmax,fNModules,0,fNModules);
2767  fh1EleEOverPMod->SetXTitle("#it{E}/#it{p}");
2768  fh1EleEOverPMod->SetYTitle("Module");
2769  outputContainer->Add(fh1EleEOverPMod);
2770  }
2771 
2772  if(fFillAllPi0Histo)
2773  {
2774  fhIM = new TH2F ("hIM","Cluster pairs (EMCAL or PHOS) Invariant mass vs reconstructed pair #it{p}_{T}, ncell > 1",
2775  nptbins,ptmin,ptmax,nmassbins,massmin,massmax);
2776  fhIM->SetXTitle("#it{p}_{T, cluster pairs} (GeV) ");
2777  fhIM->SetYTitle("M_{cluster pairs} (GeV/#it{c}^{2})");
2778  outputContainer->Add(fhIM);
2779 
2780  fhIMDiff = new TH2F ("hIMDiff","Cluster pairs (EMCAL or PHOS) Invariant mass vs reconstructed pair #it{p}_{T}, ncell > 1, different SM",
2781  nptbins,ptmin,ptmax,nmassbins,massmin,massmax);
2782  fhIMDiff->SetXTitle("#it{p}_{T, cluster pairs} (GeV) ");
2783  fhIMDiff->SetYTitle("M_{cluster pairs} (GeV/#it{c}^{2})");
2784  outputContainer->Add(fhIMDiff);
2785 
2786  fhIMSame = new TH2F ("hIMSame","Cluster pairs (EMCAL or PHOS) Invariant mass vs reconstructed pair #it{p}_{T}, ncell > 1, same SM",
2787  nptbins,ptmin,ptmax,nmassbins,massmin,massmax);
2788  fhIMSame->SetXTitle("#it{p}_{T, cluster pairs} (GeV) ");
2789  fhIMSame->SetYTitle("M_{cluster pairs} (GeV/#it{c}^{2})");
2790  outputContainer->Add(fhIMSame);
2791 
2792  if(fNModules > 12 && (GetCalorimeter() == kEMCAL || GetCalorimeter() == kDCAL))
2793  {
2794  fhIMDCAL = new TH2F ("hIMDCAL","Cluster pairs in DCAL Invariant mass vs reconstructed pair #it{p}_{T}, ncell > 1",
2795  nptbins,ptmin,ptmax,nmassbins,massmin,massmax);
2796  fhIMDCAL->SetXTitle("#it{p}_{T, cluster pairs} (GeV) ");
2797  fhIMDCAL->SetYTitle("M_{cluster pairs} (GeV/#it{c}^{2})");
2798  outputContainer->Add(fhIMDCAL);
2799 
2800  fhIMDCALDiff = new TH2F ("hIMDCALDiff","Cluster pairs in DCAL Invariant mass vs reconstructed pair #it{p}_{T}, ncell > 1, different SM",
2801  nptbins,ptmin,ptmax,nmassbins,massmin,massmax);
2802  fhIMDCALDiff->SetXTitle("#it{p}_{T, cluster pairs} (GeV) ");
2803  fhIMDCALDiff->SetYTitle("M_{cluster pairs} (GeV/#it{c}^{2})");
2804  outputContainer->Add(fhIMDCALDiff);
2805 
2806  fhIMDCALSame = new TH2F ("hIMDCALSame","Cluster pairs in DCAL Invariant mass vs reconstructed pair #it{p}_{T}, ncell > 1, same SM",
2807  nptbins,ptmin,ptmax,nmassbins,massmin,massmax);
2808  fhIMDCALSame->SetXTitle("#it{p}_{T, cluster pairs} (GeV) ");
2809  fhIMDCALSame->SetYTitle("M_{cluster pairs} (GeV/#it{c}^{2})");
2810  outputContainer->Add(fhIMDCALSame);
2811 
2812  fhIMDCALPHOS = new TH2F ("hIMDCALPHOS","Cluster pairs in DCAL-PHOS Invariant mass vs reconstructed pair #it{p}_{T}, ncell > 1",
2813  nptbins,ptmin,ptmax,nmassbins,massmin,massmax);
2814  fhIMDCALPHOS->SetXTitle("#it{p}_{T, cluster pairs} (GeV) ");
2815  fhIMDCALPHOS->SetYTitle("M_{cluster pairs} (GeV/#it{c}^{2})");
2816  outputContainer->Add(fhIMDCALPHOS);
2817 
2818  fhIMDCALPHOSSame = new TH2F ("hIMDCALPHOSSame","Cluster pairs in DCAL-PHOS Invariant mass vs reconstructed pair #it{p}_{T}, ncell > 1, same #phi sector",
2819  nptbins,ptmin,ptmax,nmassbins,massmin,massmax);
2820  fhIMDCALPHOSSame->SetXTitle("#it{p}_{T, cluster pairs} (GeV) ");
2821  fhIMDCALPHOSSame->SetYTitle("M_{cluster pairs} (GeV/#it{c}^{2})");
2822  outputContainer->Add(fhIMDCALPHOSSame);
2823  }
2824 
2825  if(fFillInvMassInEMCALWithPHOSDCalAcc)
2826  {
2827  fhIMEMCALPHOS = new TH2F ("hIMEMCALPHOS","Cluster pairs in DCAL-PHOS Invariant mass vs reconstructed pair #it{p}_{T}, ncell > 1",
2828  nptbins,ptmin,ptmax,nmassbins,massmin,massmax);
2829  fhIMEMCALPHOS->SetXTitle("#it{p}_{T, cluster pairs} (GeV) ");
2830  fhIMEMCALPHOS->SetYTitle("M_{cluster pairs} (GeV/#it{c}^{2})");
2831  outputContainer->Add(fhIMEMCALPHOS);
2832 
2833  fhIMEMCALPHOSSame = new TH2F ("hIMEMCALPHOSSame","Cluster pairs in DCAL-PHOS Invariant mass vs reconstructed pair #it{p}_{T}, ncell > 1, same #phi sector",
2834  nptbins,ptmin,ptmax,nmassbins,massmin,massmax);
2835  fhIMEMCALPHOSSame->SetXTitle("#it{p}_{T, cluster pairs} (GeV) ");
2836  fhIMEMCALPHOSSame->SetYTitle("M_{cluster pairs} (GeV/#it{c}^{2})");
2837  outputContainer->Add(fhIMEMCALPHOSSame);
2838  }
2839 
2840  if(fFillPi0PairDiffTime)
2841  {
2842  fhClusterPairDiffTimeEPi0Mass = new TH2F("hClusterPairDiffTimeEPi0Mass","cluster pair time difference vs E, only good clusters",
2843  nptbins,ptmin,ptmax, tdbins,tdmin,tdmax);
2844  fhClusterPairDiffTimeEPi0Mass->SetXTitle("#it{E}_{cluster} (GeV)");
2845  fhClusterPairDiffTimeEPi0Mass->SetYTitle("#Delta #it{t} (ns)");
2846  outputContainer->Add(fhClusterPairDiffTimeEPi0Mass);
2847 
2848  fhClusterPairDiffTimeEPi0MassSame = new TH2F("hClusterPairDiffTimeEPi0MassSame","cluster pair time difference vs E, only good clusters",
2849  nptbins,ptmin,ptmax, tdbins,tdmin,tdmax);
2850  fhClusterPairDiffTimeEPi0MassSame->SetXTitle("#it{E}_{cluster} (GeV)");
2851  fhClusterPairDiffTimeEPi0MassSame->SetYTitle("#Delta #it{t} (ns)");
2852  outputContainer->Add(fhClusterPairDiffTimeEPi0MassSame);
2853 
2854  if(fNModules > 12 && (GetCalorimeter() == kEMCAL || GetCalorimeter() == kDCAL))
2855  {
2856  fhClusterPairDiffTimeEPi0MassDCal = new TH2F("hClusterPairDiffTimeEPi0MassDCal","cluster pair time difference vs E, only good clusters",
2857  nptbins,ptmin,ptmax, tdbins,tdmin,tdmax);
2858  fhClusterPairDiffTimeEPi0MassDCal->SetXTitle("#it{E}_{cluster} (GeV)");
2859  fhClusterPairDiffTimeEPi0MassDCal->SetYTitle("#Delta #it{t} (ns)");
2860  outputContainer->Add(fhClusterPairDiffTimeEPi0MassDCal);
2861 
2862  fhClusterPairDiffTimeEPi0MassDCalSame = new TH2F("hClusterPairDiffTimeEPi0MassDCalSame","cluster pair time difference vs E, only good clusters",
2863  nptbins,ptmin,ptmax, tdbins,tdmin,tdmax);
2864  fhClusterPairDiffTimeEPi0MassDCalSame->SetXTitle("#it{E}_{cluster} (GeV)");
2865  fhClusterPairDiffTimeEPi0MassDCalSame->SetYTitle("#Delta #it{t} (ns)");
2866  outputContainer->Add(fhClusterPairDiffTimeEPi0MassDCalSame);
2867  }
2868  }
2869 
2870  fhAsym = new TH2F ("hAssym","Cluster pairs Asymmetry vs reconstructed pair energy",nptbins,ptmin,ptmax,nasymbins,asymmin,asymmax);
2871  fhAsym->SetXTitle("#it{p}_{T, cluster pairs} (GeV) ");
2872  fhAsym->SetYTitle("#it{Asymmetry}");
2873  outputContainer->Add(fhAsym);
2874 
2875  if(fFillInvMassOpenAngle)
2876  {
2877  fhOpAngle = new TH2F ("hOpeningAngle","Cluster pairs opening angle vs reconstructed pair #it{p}_{T}, ncell > 1",
2878  nptbins,ptmin,ptmax, 180,0,TMath::Pi());
2879  fhOpAngle->SetXTitle("#it{p}_{T, cluster pairs} (GeV) ");
2880  fhOpAngle->SetYTitle("Opening angle (degrees)");
2881  outputContainer->Add(fhOpAngle);
2882 
2883  Int_t nBinOpAngle = TMath::Nint(fInvMassMaxOpenAngle*TMath::RadToDeg()); // bin of 1 degree size
2884  fhIMvsOpAngle = new TH2F ("hIMvsOpAngle","Cluster pairs Invariant mass vs reconstructed pair opening angle, ncell > 1",
2885  nBinOpAngle,0.,fInvMassMaxOpenAngle,nmassbins,massmin,massmax);
2886  fhIMvsOpAngle->SetXTitle("Opening angle (degrees)");
2887  fhIMvsOpAngle->SetYTitle("M_{cluster pairs} (GeV/#it{c}^{2})");
2888  outputContainer->Add(fhIMvsOpAngle);
2889  }
2890  }
2891 
2892  if(fFillAllPosHisto2)
2893  {
2894  fhXYZ = new TH3F ("hXYZ","Cluster: #it{x} vs #it{y} vs #it{z}",xbins,xmin,xmax,ybins,ymin,ymax,zbins,zmin,zmax);
2895  fhXYZ->SetXTitle("#it{x} (cm)");
2896  fhXYZ->SetYTitle("#it{y} (cm)");
2897  fhXYZ->SetZTitle("#it{z} (cm) ");
2898  outputContainer->Add(fhXYZ);
2899 
2900  fhXE = new TH2F ("hXE","Cluster X position vs cluster energy",xbins,xmin,xmax,nptbins,ptmin,ptmax);
2901  fhXE->SetXTitle("#it{x} (cm)");
2902  fhXE->SetYTitle("#it{E} (GeV)");
2903  outputContainer->Add(fhXE);
2904 
2905  fhYE = new TH2F ("hYE","Cluster Y position vs cluster energy",ybins,ymin,ymax,nptbins,ptmin,ptmax);
2906  fhYE->SetXTitle("#it{y} (cm)");
2907  fhYE->SetYTitle("#it{E} (GeV)");
2908  outputContainer->Add(fhYE);
2909 
2910  fhZE = new TH2F ("hZE","Cluster Z position vs cluster energy",zbins,zmin,zmax,nptbins,ptmin,ptmax);
2911  fhZE->SetXTitle("#it{z} (cm)");
2912  fhZE->SetYTitle("#it{E} (GeV)");
2913  outputContainer->Add(fhZE);
2914 
2915  fhRE = new TH2F ("hRE","Cluster R position vs cluster energy",rbins,rmin,rmax,nptbins,ptmin,ptmax);
2916  fhRE->SetXTitle("r = #sqrt{x^{2}+y^{2}} (cm)");
2917  fhRE->SetYTitle("#it{E} (GeV)");
2918  outputContainer->Add(fhRE);
2919 
2920  fhXNCells = new TH2F ("hXNCells","Cluster X position vs N Cells per Cluster",xbins,xmin,xmax,nceclbins,nceclmin,nceclmax);
2921  fhXNCells->SetXTitle("#it{x} (cm)");
2922  fhXNCells->SetYTitle("N cells per cluster");
2923  outputContainer->Add(fhXNCells);
2924 
2925  fhYNCells = new TH2F ("hYNCells","Cluster Y position vs N Cells per Cluster",ybins,ymin,ymax,nceclbins,nceclmin,nceclmax);
2926  fhYNCells->SetXTitle("#it{y} (cm)");
2927  fhYNCells->SetYTitle("N cells per cluster");
2928  outputContainer->Add(fhYNCells);
2929 
2930  fhZNCells = new TH2F ("hZNCells","Cluster Z position vs N Cells per Cluster",zbins,zmin,zmax,nceclbins,nceclmin,nceclmax);
2931  fhZNCells->SetXTitle("#it{z} (cm)");
2932  fhZNCells->SetYTitle("N cells per cluster");
2933  outputContainer->Add(fhZNCells);
2934 
2935  fhRNCells = new TH2F ("hRNCells","Cluster R position vs N Cells per Cluster",rbins,rmin,rmax,nceclbins,nceclmin,nceclmax);
2936  fhRNCells->SetXTitle("r = #sqrt{x^{2}+y^{2}} (cm)");
2937  fhRNCells->SetYTitle("N cells per cluster");
2938  outputContainer->Add(fhRNCells);
2939  }
2940 
2941  if(fFillAllPosHisto)
2942  {
2943  fhRCellE = new TH2F ("hRCellE","Cell R position vs cell energy",rbins,rmin,rmax,nptbins,ptmin,ptmax);
2944  fhRCellE->SetXTitle("r = #sqrt{x^{2}+y^{2}} (cm)");
2945  fhRCellE->SetYTitle("#it{E} (GeV)");
2946  outputContainer->Add(fhRCellE);
2947 
2948  fhXCellE = new TH2F ("hXCellE","Cell X position vs cell energy",xbins,xmin,xmax,nptbins,ptmin,ptmax);
2949  fhXCellE->SetXTitle("#it{x} (cm)");
2950  fhXCellE->SetYTitle("#it{E} (GeV)");
2951  outputContainer->Add(fhXCellE);
2952 
2953  fhYCellE = new TH2F ("hYCellE","Cell Y position vs cell energy",ybins,ymin,ymax,nptbins,ptmin,ptmax);
2954  fhYCellE->SetXTitle("#it{y} (cm)");
2955  fhYCellE->SetYTitle("#it{E} (GeV)");
2956  outputContainer->Add(fhYCellE);
2957 
2958  fhZCellE = new TH2F ("hZCellE","Cell Z position vs cell energy",zbins,zmin,zmax,nptbins,ptmin,ptmax);
2959  fhZCellE->SetXTitle("#it{z} (cm)");
2960  fhZCellE->SetYTitle("#it{E} (GeV)");
2961  outputContainer->Add(fhZCellE);
2962 
2963  fhXYZCell = new TH3F ("hXYZCell","Cell : #it{x} vs #it{y} vs #it{z}",xbins,xmin,xmax,ybins,ymin,ymax,zbins,zmin,zmax);
2964  fhXYZCell->SetXTitle("#it{x} (cm)");
2965  fhXYZCell->SetYTitle("#it{y} (cm)");
2966  fhXYZCell->SetZTitle("#it{z} (cm)");
2967  outputContainer->Add(fhXYZCell);
2968 
2969  Float_t dx = TMath::Abs(xmin)+TMath::Abs(xmax);
2970  Float_t dy = TMath::Abs(ymin)+TMath::Abs(ymax);
2971  Float_t dz = TMath::Abs(zmin)+TMath::Abs(zmax);
2972  Float_t dr = TMath::Abs(rmin)+TMath::Abs(rmax);
2973 
2974  fhDeltaCellClusterRNCells = new TH2F ("hDeltaCellClusterRNCells","Cluster-Cell R position vs N Cells per Cluster",rbins*2,-dr,dr,nceclbins,nceclmin,nceclmax);
2975  fhDeltaCellClusterRNCells->SetXTitle("#it{r} = #sqrt{x^{2}+y^{2}}, #it{r}_{clus}-#it{r}_{cell} (cm)");
2976  fhDeltaCellClusterRNCells->SetYTitle("#it{n}_{cells per cluster}");
2977  outputContainer->Add(fhDeltaCellClusterRNCells);
2978 
2979  fhDeltaCellClusterXNCells = new TH2F ("hDeltaCellClusterXNCells","Cluster-Cell X position vs N Cells per Cluster",xbins*2,-dx,dx,nceclbins,nceclmin,nceclmax);
2980  fhDeltaCellClusterXNCells->SetXTitle("#it{x}_{clus}-#it{x}_{cell} (cm)");
2981  fhDeltaCellClusterXNCells->SetYTitle("#it{n}_{cells per cluster}");
2982  outputContainer->Add(fhDeltaCellClusterXNCells);
2983 
2984  fhDeltaCellClusterYNCells = new TH2F ("hDeltaCellClusterYNCells","Cluster-Cell Y position vs N Cells per Cluster",ybins*2,-dy,dy,nceclbins,nceclmin,nceclmax);
2985  fhDeltaCellClusterYNCells->SetXTitle("#it{y}_{clus}-#it{y}_{cell} (cm)");
2986  fhDeltaCellClusterYNCells->SetYTitle("N cells per cluster");
2987  outputContainer->Add(fhDeltaCellClusterYNCells);
2988 
2989  fhDeltaCellClusterZNCells = new TH2F ("hDeltaCellClusterZNCells","Cluster-Cell Z position vs N Cells per Cluster",zbins*2,-dz,dz,nceclbins,nceclmin,nceclmax);
2990  fhDeltaCellClusterZNCells->SetXTitle("#it{z}_{clus}-#it{z}_{cell} (cm)");
2991  fhDeltaCellClusterZNCells->SetYTitle("#it{n}_{cells per cluster}");
2992  outputContainer->Add(fhDeltaCellClusterZNCells);
2993 
2994  fhDeltaCellClusterRE = new TH2F ("hDeltaCellClusterRE","Cluster-Cell R position vs cluster energy",rbins*2,-dr,dr,nptbins,ptmin,ptmax);
2995  fhDeltaCellClusterRE->SetXTitle("#it{r} = #sqrt{x^{2}+y^{2}}, #it{r}_{clus}-#it{r}_{cell} (cm)");
2996  fhDeltaCellClusterRE->SetYTitle("#it{E} (GeV)");
2997  outputContainer->Add(fhDeltaCellClusterRE);
2998 
2999  fhDeltaCellClusterXE = new TH2F ("hDeltaCellClusterXE","Cluster-Cell X position vs cluster energy",xbins*2,-dx,dx,nptbins,ptmin,ptmax);
3000  fhDeltaCellClusterXE->SetXTitle("#it{x}_{clus}-#it{x}_{cell} (cm)");
3001  fhDeltaCellClusterXE->SetYTitle("#it{E} (GeV)");
3002  outputContainer->Add(fhDeltaCellClusterXE);
3003 
3004  fhDeltaCellClusterYE = new TH2F ("hDeltaCellClusterYE","Cluster-Cell Y position vs cluster energy",ybins*2,-dy,dy,nptbins,ptmin,ptmax);
3005  fhDeltaCellClusterYE->SetXTitle("#it{y}_{clus}-#it{y}_{cell} (cm)");
3006  fhDeltaCellClusterYE->SetYTitle("#it{E} (GeV)");
3007  outputContainer->Add(fhDeltaCellClusterYE);
3008 
3009  fhDeltaCellClusterZE = new TH2F ("hDeltaCellClusterZE","Cluster-Cell Z position vs cluster energy",zbins*2,-dz,dz,nptbins,ptmin,ptmax);
3010  fhDeltaCellClusterZE->SetXTitle("#it{z}_{clus}-#it{z}_{cell} (cm)");
3011  fhDeltaCellClusterZE->SetYTitle("#it{E} (GeV)");
3012  outputContainer->Add(fhDeltaCellClusterZE);
3013 
3014  if(fFillAllTH3)
3015  {
3016  fhEtaPhiAmpCell = new TH3F ("hEtaPhiAmpCell","Cell #eta vs cell #phi vs cell energy",
3017  netabins,etamin,etamax,nphibins,phimin,phimax,nptbins,ptmin,ptmax);
3018  fhEtaPhiAmpCell->SetXTitle("#eta ");
3019  fhEtaPhiAmpCell->SetYTitle("#phi (rad)");
3020  fhEtaPhiAmpCell->SetZTitle("#it{E} (GeV) ");
3021  outputContainer->Add(fhEtaPhiAmpCell);
3022  }
3023  else
3024  {
3025  fhEtaPhiCell = new TH2F ("hEtaPhiCell","Cell #eta vs cell #phi vs cell energy",
3026  netabins,etamin,etamax,nphibins,phimin,phimax);
3027  fhEtaPhiCell->SetXTitle("#eta ");
3028  fhEtaPhiCell->SetYTitle("#phi (rad)");
3029  outputContainer->Add(fhEtaPhiCell);
3030  }
3031  }
3032 
3033  // Calorimeter cells
3034 
3035  fhNCells = new TH1F ("hNCells","# cells", ncebins,ncemin+0.5,ncemax);
3036  fhNCells->SetXTitle("#it{n}_{cells}");
3037  outputContainer->Add(fhNCells);
3038 
3039  fhNCellsCutAmpMin = new TH1F ("hNCellsCutAmpMin",Form("# cells amp > %1.2f-%1.2f",fEMCALCellAmpMin,fPHOSCellAmpMin), ncebins,ncemin+0.5,ncemax);
3040  fhNCellsCutAmpMin->SetXTitle("#it{n}_{cells}");
3041  outputContainer->Add(fhNCellsCutAmpMin);
3042 
3043  fhAmplitude = new TH1F ("hAmplitude","#it{E}_{cell}", nptbins*2,ptmin,ptmax);
3044  fhAmplitude->SetXTitle("#it{E}_{cell} (GeV)");
3045  outputContainer->Add(fhAmplitude);
3046 
3047  fhAmpId = new TH2F ("hAmpId","#it{E}_{cell}", nfineptbins,ptfinemin,ptfinemax,fNMaxRows*fNMaxCols*fNModules,0,fNMaxRows*fNMaxCols*fNModules);
3048  fhAmpId->SetXTitle("#it{E}_{cell} (GeV)");
3049  outputContainer->Add(fhAmpId);
3050 
3051  fhAmpIdLowGain = new TH2F ("hAmpIdLG","Low gain: #it{E}_{cell}", nfineptbins,15,ptfinemax+15,fNMaxRows*fNMaxCols*fNModules,0,fNMaxRows*fNMaxCols*fNModules);
3052  fhAmpIdLowGain->SetXTitle("#it{E}_{cell} (GeV)");
3053  outputContainer->Add(fhAmpIdLowGain);
3054 
3055  if(fFillAllCellTimeHisto)
3056  {
3057  fhCellTimeSpreadRespectToCellMax = new TH2F ("hCellTimeSpreadRespectToCellMax","t_{cell max}-t_{cell i} per cluster", nptbins,ptmin,ptmax,tdbins,tdmin,tdmax);
3058  fhCellTimeSpreadRespectToCellMax->SetXTitle("#it{E} (GeV)");
3059  fhCellTimeSpreadRespectToCellMax->SetYTitle("#Delta #it{t}_{cell max-i} (ns)");
3060  outputContainer->Add(fhCellTimeSpreadRespectToCellMax);
3061 
3062 // fhClusterMaxCellDiffAverageTime = new TH2F ("hClusterMaxCellDiffAverageTime","t_{cell max}-t_{average} per cluster", nptbins,ptmin,ptmax, tdbins,tdmin,tdmax);
3063 // fhClusterMaxCellDiffAverageTime->SetXTitle("#it{E} (GeV)");
3064 // fhClusterMaxCellDiffAverageTime->SetYTitle("#Delta #it{t}_{cell max - average} (ns)");
3065 // outputContainer->Add(fhClusterMaxCellDiffAverageTime);
3066 //
3067 // fhClusterMaxCellDiffWeightedTime = new TH2F ("hClusterMaxCellDiffWeightedTime","t_{cell max}-t_{weighted} per cluster", nptbins,ptmin,ptmax, tdbins,tdmin,tdmax);
3068 // fhClusterMaxCellDiffWeightedTime->SetXTitle("#it{E} (GeV)");
3069 // fhClusterMaxCellDiffWeightedTime->SetYTitle("#Delta #it{t}_{cell max - weighted} (ns)");
3070 // outputContainer->Add(fhClusterMaxCellDiffWeightedTime);
3071 
3072  fhCellIdCellLargeTimeSpread= new TH1F ("hCellIdCellLargeTimeSpread","Cells with time 100 ns larger than cell max in cluster ",
3073  fNMaxCols*fNMaxRows*fNModules,0,fNMaxCols*fNMaxRows*fNModules);
3074  fhCellIdCellLargeTimeSpread->SetXTitle("Absolute Cell Id");
3075  outputContainer->Add(fhCellIdCellLargeTimeSpread);
3076 
3077  fhTime = new TH1F ("hTime","#it{t}_{cell}",ntimebins,timemin,timemax);
3078  fhTime->SetXTitle("#it{t}_{cell} (ns)");
3079  outputContainer->Add(fhTime);
3080 
3081 // fhTimeVz = new TH2F ("hTimeVz","#it{t}_{cell} vs vertex, amplitude > 0.5 GeV",100, 0, 50,ntimebins,timemin,timemax);
3082 // fhTimeVz->SetXTitle("|v_{z}| (cm)");
3083 // fhTimeVz->SetYTitle("#it{t}_{cell} (ns)");
3084 // outputContainer->Add(fhTimeVz);
3085 
3086  fhTimeId = new TH2F ("hTimeId","#it{t}_{cell} vs Absolute Id",
3087  ntimebins,timemin,timemax,fNMaxRows*fNMaxCols*fNModules,0,fNMaxRows*fNMaxCols*fNModules);
3088  fhTimeId->SetXTitle("#it{t}_{cell} (ns)");
3089  fhTimeId->SetYTitle("Cell Absolute Id");
3090  outputContainer->Add(fhTimeId);
3091 
3092  fhTimeAmp = new TH2F ("hTimeAmp","#it{t}_{cell} vs #it{E}_{cell}",nptbins*2,ptmin,ptmax,ntimebins,timemin,timemax);
3093  fhTimeAmp->SetYTitle("#it{t}_{cell} (ns)");
3094  fhTimeAmp->SetXTitle("#it{E}_{cell} (GeV)");
3095  outputContainer->Add(fhTimeAmp);
3096 
3097  for(Int_t bc = 0; bc < 4; bc++)
3098  {
3099  fhTimePerSMPerBC[bc] = new TH2F (Form("hTimePerSM_BC%d",bc),
3100  Form("#it{t}_{cell} vs super-module, for BC/4=%d",bc),
3101  ntimebins,timemin,timemax,fNModules,0,fNModules);
3102  fhTimePerSMPerBC[bc]->SetXTitle("#it{t}_{cell} (ns)");
3103  fhTimePerSMPerBC[bc]->SetYTitle("Module");
3104  outputContainer->Add(fhTimePerSMPerBC[bc]);
3105  }
3106 
3107  fhTimeIdLowGain = new TH2F ("hTimeIdLG","Low gain: #it{t}_{cell} vs Absolute Id",
3108  ntimebins,timemin,timemax,fNMaxRows*fNMaxCols*fNModules,0,fNMaxRows*fNMaxCols*fNModules);
3109  fhTimeIdLowGain->SetXTitle("#it{t}_{cell} (ns)");
3110  fhTimeIdLowGain->SetYTitle("Cell Absolute Id");
3111  outputContainer->Add(fhTimeIdLowGain);
3112 
3113  fhTimeAmpLowGain = new TH2F ("hTimeAmpLG","Low gain: #it{t}_{cell} vs #it{E}_{cell}",nptbins*2,ptmin,ptmax,ntimebins,timemin,timemax);
3114  fhTimeAmpLowGain->SetYTitle("#it{t}_{cell} (ns)");
3115  fhTimeAmpLowGain->SetXTitle("#it{E}_{cell} (GeV)");
3116  outputContainer->Add(fhTimeAmpLowGain);
3117  }
3118 
3119  fhCellECross = new TH2F ("hCellECross","1 - Energy in cross around cell / cell energy",
3120  nptbins,ptmin,ptmax, 400,-1,1.);
3121  fhCellECross->SetXTitle("#it{E}_{cell} (GeV) ");
3122  fhCellECross->SetYTitle("1- #it{E}_{cross}/#it{E}_{cell}");
3123  outputContainer->Add(fhCellECross);
3124 
3125 
3126  if(fCorrelate)
3127  {
3128  // PHOS vs EMCAL
3129  fhEMCALPHOSCorrNClusters = new TH2F ("hEMCALPHOSCorrNClusters","# clusters in EMCAL vs PHOS", nclbins,nclmin,nclmax,nclbins,nclmin,nclmax);
3130  fhEMCALPHOSCorrNClusters->SetXTitle("number of clusters in EMCAL");
3131  fhEMCALPHOSCorrNClusters->SetYTitle("number of clusters in PHOS");
3132  outputContainer->Add(fhEMCALPHOSCorrNClusters);
3133 
3134  fhEMCALPHOSCorrEClusters = new TH2F ("hEMCALPHOSCorrEClusters","summed energy of clusters in EMCAL vs PHOS", nptbins,ptmin,ptmax*2,nptbins,ptmin,ptmax*2);
3135  fhEMCALPHOSCorrEClusters->SetXTitle("#Sigma #it{E} of clusters in EMCAL (GeV)");
3136  fhEMCALPHOSCorrEClusters->SetYTitle("#Sigma #it{E} of clusters in PHOS (GeV)");
3137  outputContainer->Add(fhEMCALPHOSCorrEClusters);
3138 
3139  fhEMCALPHOSCorrNCells = new TH2F ("hEMCALPHOSCorrNCells","# Cells in EMCAL vs PHOS", ncebins,ncemin,ncemax, ncebins,ncemin,ncemax);
3140  fhEMCALPHOSCorrNCells->SetXTitle("number of Cells in EMCAL");
3141  fhEMCALPHOSCorrNCells->SetYTitle("number of Cells in PHOS");
3142  outputContainer->Add(fhEMCALPHOSCorrNCells);
3143 
3144  fhEMCALPHOSCorrECells = new TH2F ("hEMCALPHOSCorrECells","summed energy of Cells in EMCAL vs PHOS", nptbins*2,ptmin,ptmax*4,nptbins*2,ptmin,ptmax*4);
3145  fhEMCALPHOSCorrECells->SetXTitle("#Sigma #it{E} of Cells in EMCAL (GeV)");
3146  fhEMCALPHOSCorrECells->SetYTitle("#Sigma #it{E} of Cells in PHOS (GeV)");
3147  outputContainer->Add(fhEMCALPHOSCorrECells);
3148 
3149  // DCal vs EMCAL
3150  fhEMCALDCALCorrNClusters = new TH2F ("hEMCALDCALCorrNClusters","# clusters in EMCAL vs DCAL", nclbins,nclmin,nclmax,nclbins,nclmin,nclmax);
3151  fhEMCALDCALCorrNClusters->SetXTitle("number of clusters in EMCAL");
3152  fhEMCALDCALCorrNClusters->SetYTitle("number of clusters in DCAL");
3153  outputContainer->Add(fhEMCALDCALCorrNClusters);
3154 
3155  fhEMCALDCALCorrEClusters = new TH2F ("hEMCALDCALCorrEClusters","summed energy of clusters in EMCAL vs DCAL", nptbins,ptmin,ptmax*2,nptbins,ptmin,ptmax*2);
3156  fhEMCALDCALCorrEClusters->SetXTitle("#Sigma #it{E} of clusters in EMCAL (GeV)");
3157  fhEMCALDCALCorrEClusters->SetYTitle("#Sigma #it{E} of clusters in DCAL (GeV)");
3158  outputContainer->Add(fhEMCALDCALCorrEClusters);
3159 
3160  fhEMCALDCALCorrNCells = new TH2F ("hEMCALDCALCorrNCells","# Cells in EMCAL vs DCAL", ncebins,ncemin,ncemax, ncebins,ncemin,ncemax);
3161  fhEMCALDCALCorrNCells->SetXTitle("number of Cells in EMCAL");
3162  fhEMCALDCALCorrNCells->SetYTitle("number of Cells in DCAL");
3163  outputContainer->Add(fhEMCALDCALCorrNCells);
3164 
3165  fhEMCALDCALCorrECells = new TH2F ("hEMCALDCALCorrECells","summed energy of Cells in EMCAL vs DCAL", nptbins*2,ptmin,ptmax*4,nptbins*2,ptmin,ptmax*4);
3166  fhEMCALDCALCorrECells->SetXTitle("#Sigma #it{E} of Cells in EMCAL (GeV)");
3167  fhEMCALDCALCorrECells->SetYTitle("#Sigma #it{E} of Cells in DCAL (GeV)");
3168  outputContainer->Add(fhEMCALDCALCorrECells);
3169 
3170 
3171  // DCAL vs PHOS
3172  fhDCALPHOSCorrNClusters = new TH2F ("hDCALPHOSCorrNClusters","# clusters in DCAL vs PHOS", nclbins,nclmin,nclmax,nclbins,nclmin,nclmax);
3173  fhDCALPHOSCorrNClusters->SetXTitle("number of clusters in DCAL");
3174  fhDCALPHOSCorrNClusters->SetYTitle("number of clusters in PHOS");
3175  outputContainer->Add(fhDCALPHOSCorrNClusters);
3176 
3177  fhDCALPHOSCorrEClusters = new TH2F ("hDCALPHOSCorrEClusters","summed energy of clusters in DCAL vs PHOS", nptbins,ptmin,ptmax*2,nptbins,ptmin,ptmax*2);
3178  fhDCALPHOSCorrEClusters->SetXTitle("#Sigma #it{E} of clusters in DCAL (GeV)");
3179  fhDCALPHOSCorrEClusters->SetYTitle("#Sigma #it{E} of clusters in PHOS (GeV)");
3180  outputContainer->Add(fhDCALPHOSCorrEClusters);
3181 
3182  fhDCALPHOSCorrNCells = new TH2F ("hDCALPHOSCorrNCells","# Cells in DCAL vs PHOS", ncebins,ncemin,ncemax, ncebins,ncemin,ncemax);
3183  fhDCALPHOSCorrNCells->SetXTitle("number of Cells in DCAL");
3184  fhDCALPHOSCorrNCells->SetYTitle("number of Cells in PHOS");
3185  outputContainer->Add(fhDCALPHOSCorrNCells);
3186 
3187  fhDCALPHOSCorrECells = new TH2F ("hDCALPHOSCorrECells","summed energy of Cells in DCAL vs PHOS", nptbins*2,ptmin,ptmax*4,nptbins*2,ptmin,ptmax*4);
3188  fhDCALPHOSCorrECells->SetXTitle("#Sigma #it{E} of Cells in DCAL (GeV)");
3189  fhDCALPHOSCorrECells->SetYTitle("#Sigma #it{E} of Cells in PHOS (GeV)");
3190  outputContainer->Add(fhDCALPHOSCorrECells);
3191 
3192  // Calorimeter vs. V0 signal
3193 
3194  fhCaloV0SCorrNClusters = new TH2F ("hCaloV0SNClusters",Form("# clusters in %s vs V0 signal",GetCalorimeterString().Data()), nv0sbins,nv0smin,nv0smax,nclbins,nclmin,nclmax);
3195  fhCaloV0SCorrNClusters->SetXTitle("V0 signal");
3196  fhCaloV0SCorrNClusters->SetYTitle(Form("number of clusters in %s",GetCalorimeterString().Data()));
3197  outputContainer->Add(fhCaloV0SCorrNClusters);
3198 
3199  fhCaloV0SCorrEClusters = new TH2F ("hCaloV0SEClusters",Form("summed energy of clusters in %s vs V0 signal",GetCalorimeterString().Data()), nv0sbins,nv0smin,nv0smax,nptbins,ptmin,ptmax*2);
3200  fhCaloV0SCorrEClusters->SetXTitle("V0 signal");
3201  fhCaloV0SCorrEClusters->SetYTitle(Form("#Sigma #it{E} of clusters in %s (GeV)",GetCalorimeterString().Data()));
3202  outputContainer->Add(fhCaloV0SCorrEClusters);
3203 
3204  fhCaloV0SCorrNCells = new TH2F ("hCaloV0SNCells",Form("# Cells in %s vs V0 signal",GetCalorimeterString().Data()), nv0sbins,nv0smin,nv0smax, ncebins,ncemin,ncemax);
3205  fhCaloV0SCorrNCells->SetXTitle("V0 signal");
3206  fhCaloV0SCorrNCells->SetYTitle(Form("number of Cells in %s",GetCalorimeterString().Data()));
3207  outputContainer->Add(fhCaloV0SCorrNCells);
3208 
3209  fhCaloV0SCorrECells = new TH2F ("hCaloV0SECells",Form("summed energy of Cells in %s vs V0 signal",GetCalorimeterString().Data()), nv0sbins,nv0smin,nv0smax,nptbins,ptmin,ptmax*2);
3210  fhCaloV0SCorrECells->SetXTitle("V0 signal");
3211  fhCaloV0SCorrECells->SetYTitle(Form("#Sigma #it{E} of Cells in %s (GeV)",GetCalorimeterString().Data()));
3212  outputContainer->Add(fhCaloV0SCorrECells);
3213 
3214  // Calorimeter vs V0 multiplicity
3215 
3216  fhCaloV0MCorrNClusters = new TH2F ("hCaloV0MNClusters",Form("# clusters in %s vs V0 signal",GetCalorimeterString().Data()), nv0mbins,nv0mmin,nv0mmax,nclbins,nclmin,nclmax);
3217  fhCaloV0MCorrNClusters->SetXTitle("V0 signal");
3218  fhCaloV0MCorrNClusters->SetYTitle(Form("number of clusters in %s",GetCalorimeterString().Data()));
3219  outputContainer->Add(fhCaloV0MCorrNClusters);
3220 
3221  fhCaloV0MCorrEClusters = new TH2F ("hCaloV0MEClusters",Form("summed energy of clusters in %s vs V0 signal",GetCalorimeterString().Data()), nv0mbins,nv0mmin,nv0mmax,nptbins,ptmin,ptmax*2);
3222  fhCaloV0MCorrEClusters->SetXTitle("V0 signal");
3223  fhCaloV0MCorrEClusters->SetYTitle(Form("#Sigma #it{E} of clusters in %s (GeV)",GetCalorimeterString().Data()));
3224  outputContainer->Add(fhCaloV0MCorrEClusters);
3225 
3226  fhCaloV0MCorrNCells = new TH2F ("hCaloV0MNCells",Form("# Cells in %s vs V0 signal",GetCalorimeterString().Data()), nv0mbins,nv0mmin,nv0mmax, ncebins,ncemin,ncemax);
3227  fhCaloV0MCorrNCells->SetXTitle("V0 signal");
3228  fhCaloV0MCorrNCells->SetYTitle(Form("number of Cells in %s",GetCalorimeterString().Data()));
3229  outputContainer->Add(fhCaloV0MCorrNCells);
3230 
3231  fhCaloV0MCorrECells = new TH2F ("hCaloV0MECells",Form("summed energy of Cells in %s vs V0 signal",GetCalorimeterString().Data()), nv0mbins,nv0mmin,nv0mmax,nptbins,ptmin,ptmax*2);
3232  fhCaloV0MCorrECells->SetXTitle("V0 signal");
3233  fhCaloV0MCorrECells->SetYTitle(Form("#Sigma #it{E} of Cells in %s (GeV)",GetCalorimeterString().Data()));
3234  outputContainer->Add(fhCaloV0MCorrECells);
3235 
3236  //Calorimeter VS Track multiplicity
3237  fhCaloTrackMCorrNClusters = new TH2F ("hCaloTrackMNClusters",Form("# clusters in %s vs # tracks",GetCalorimeterString().Data()), ntrmbins,ntrmmin,ntrmmax,nclbins,nclmin,nclmax);
3238  fhCaloTrackMCorrNClusters->SetXTitle("# tracks");
3239  fhCaloTrackMCorrNClusters->SetYTitle(Form("number of clusters in %s",GetCalorimeterString().Data()));
3240  outputContainer->Add(fhCaloTrackMCorrNClusters);
3241 
3242  fhCaloTrackMCorrEClusters = new TH2F ("hCaloTrackMEClusters",Form("summed energy of clusters in %s vs # tracks",GetCalorimeterString().Data()), ntrmbins,ntrmmin,ntrmmax,nptbins,ptmin,ptmax*2);
3243  fhCaloTrackMCorrEClusters->SetXTitle("# tracks");
3244  fhCaloTrackMCorrEClusters->SetYTitle(Form("#Sigma #it{E} of clusters in %s (GeV)",GetCalorimeterString().Data()));
3245  outputContainer->Add(fhCaloTrackMCorrEClusters);
3246 
3247  fhCaloTrackMCorrNCells = new TH2F ("hCaloTrackMNCells",Form("# Cells in %s vs # tracks",GetCalorimeterString().Data()), ntrmbins,ntrmmin,ntrmmax, ncebins,ncemin,ncemax);
3248  fhCaloTrackMCorrNCells->SetXTitle("# tracks");
3249  fhCaloTrackMCorrNCells->SetYTitle(Form("number of Cells in %s",GetCalorimeterString().Data()));
3250  outputContainer->Add(fhCaloTrackMCorrNCells);
3251 
3252  fhCaloTrackMCorrECells = new TH2F ("hCaloTrackMECells",Form("summed energy of Cells in %s vs # tracks",GetCalorimeterString().Data()), ntrmbins,ntrmmin,ntrmmax,nptbins,ptmin,ptmax*2);
3253  fhCaloTrackMCorrECells->SetXTitle("# tracks");
3254  fhCaloTrackMCorrECells->SetYTitle(Form("#Sigma #it{E} of Cells in %s (GeV)",GetCalorimeterString().Data()));
3255  outputContainer->Add(fhCaloTrackMCorrECells);
3256 
3257  fhCaloCenNClusters = new TH2F ("hCaloCenNClusters","# clusters in calorimeter vs centrality",100,0,100,nclbins,nclmin,nclmax);
3258  fhCaloCenNClusters->SetYTitle("number of clusters in calorimeter");
3259  fhCaloCenNClusters->SetXTitle("Centrality");
3260  outputContainer->Add(fhCaloCenNClusters);
3261 
3262  fhCaloCenEClusters = new TH2F ("hCaloCenEClusters","summed energy of clusters in calorimeter vs centrality",100,0,100,nptbins,ptmin,ptmax*2);
3263  fhCaloCenEClusters->SetYTitle("#Sigma #it{E} of clusters in calorimeter (GeV)");
3264  fhCaloCenEClusters->SetXTitle("Centrality");
3265  outputContainer->Add(fhCaloCenEClusters);
3266 
3267  fhCaloCenNCells = new TH2F ("hCaloCenNCells","# Cells in calorimeter vs centrality",100,0,100,ncebins,ncemin,ncemax);
3268  fhCaloCenNCells->SetYTitle("number of Cells in calorimeter");
3269  fhCaloCenNCells->SetXTitle("Centrality");
3270  outputContainer->Add(fhCaloCenNCells);
3271 
3272  fhCaloCenECells = new TH2F ("hCaloCenECells","summed energy of Cells in calorimeter vs centrality",100,0,100,nptbins*2,ptmin,ptmax*4);
3273  fhCaloCenECells->SetYTitle("#Sigma #it{E} of Cells in calorimeter (GeV)");
3274  fhCaloCenECells->SetXTitle("Centrality");
3275  outputContainer->Add(fhCaloCenECells);
3276 
3277  fhCaloEvPNClusters = new TH2F ("hCaloEvPNClusters","# clusters in calorimeter vs event plane angle",100,0,TMath::Pi(),nclbins,nclmin,nclmax);
3278  fhCaloEvPNClusters->SetYTitle("number of clusters in calorimeter");
3279  fhCaloEvPNClusters->SetXTitle("Event plane angle (rad)");
3280  outputContainer->Add(fhCaloEvPNClusters);
3281 
3282  fhCaloEvPEClusters = new TH2F ("hCaloEvPEClusters","summed energy of clusters in calorimeter vs event plane angle",100,0,TMath::Pi(),nptbins,ptmin,ptmax*2);
3283  fhCaloEvPEClusters->SetYTitle("#Sigma #it{E} of clusters in calorimeter (GeV)");
3284  fhCaloEvPEClusters->SetXTitle("Event plane angle (rad)");
3285  outputContainer->Add(fhCaloEvPEClusters);
3286 
3287  fhCaloEvPNCells = new TH2F ("hCaloEvPNCells","# Cells in calorimeter vs event plane angle",100,0,TMath::Pi(),ncebins,ncemin,ncemax);
3288  fhCaloEvPNCells->SetYTitle("number of Cells in calorimeter");
3289  fhCaloEvPNCells->SetXTitle("Event plane angle (rad)");
3290  outputContainer->Add(fhCaloEvPNCells);
3291 
3292  fhCaloEvPECells = new TH2F ("hCaloEvPECells","summed energy of Cells in calorimeter vs event plane angle",100,0,TMath::Pi(),nptbins*2,ptmin,ptmax*4);
3293  fhCaloEvPECells->SetYTitle("#Sigma #it{E} of Cells in calorimeter (GeV)");
3294  fhCaloEvPECells->SetXTitle("Event plane angle (rad)");
3295  outputContainer->Add(fhCaloEvPECells);
3296  } // correlate calorimeters
3297 
3298  // Module histograms
3299 
3300  fhEMod = new TH2F ("hE_Mod","Cluster reconstructed Energy in each present Module",nptbins,ptmin,ptmax,fNModules,0,fNModules);
3301  fhEMod->SetXTitle("#it{E} (GeV)");
3302  fhEMod->SetYTitle("Module");
3303  outputContainer->Add(fhEMod);
3304 
3305  fhAmpMod = new TH2F ("hAmp_Mod","Cell energy in each present Module",nptbins,ptmin,ptmax,fNModules,0,fNModules);
3306  fhAmpMod->SetXTitle("#it{E} (GeV)");
3307  fhAmpMod->SetYTitle("Module");
3308  outputContainer->Add(fhAmpMod);
3309 
3310  fhEWeirdMod = new TH2F ("hEWeird_Mod","Cluster reconstructed Energy in each present Module, ridiculously large E",200,0,10000,fNModules,0,fNModules);
3311  fhEWeirdMod->SetXTitle("#it{E} (GeV)");
3312  fhEWeirdMod->SetYTitle("Module");
3313  outputContainer->Add(fhEWeirdMod);
3314 
3315  fhAmpWeirdMod = new TH2F ("hAmpWeird_Mod","Cell energy in each present Module, ridiculously large E",200,0,10000,fNModules,0,fNModules);
3316  fhAmpWeirdMod->SetXTitle("#it{E} (GeV)");
3317  fhAmpWeirdMod->SetYTitle("Module");
3318  outputContainer->Add(fhAmpWeirdMod);
3319 
3320  if(fFillAllCellTimeHisto)
3321  {
3322  fhTimeMod = new TH2F ("hTime_Mod","Cell time in each present Module",ntimebins,timemin,timemax,fNModules,0,fNModules);
3323  fhTimeMod->SetXTitle("t (ns)");
3324  fhTimeMod->SetYTitle("Module");
3325  outputContainer->Add(fhTimeMod);
3326  }
3327 
3328  fhNClustersMod = new TH2F ("hNClusters_Mod","# clusters vs Module", nclbins,nclmin+0.5,nclmax,fNModules,0,fNModules);
3329  fhNClustersMod->SetXTitle("number of clusters");
3330  fhNClustersMod->SetYTitle("Module");
3331  outputContainer->Add(fhNClustersMod);
3332 
3333  fhNCellsMod = new TH2F ("hNCells_Mod","# cells vs Module", ncebins,ncemin+0.5,ncemax,fNModules,0,fNModules);
3334  fhNCellsMod->SetXTitle("#it{n}_{cells}");
3335  fhNCellsMod->SetYTitle("Module");
3336  outputContainer->Add(fhNCellsMod);
3337 
3338  fhSumClustersEnergyMod = new TH2F ("hSumClustersEnergy_Mod","# clusters vs Module", 1000, 0, 2000,fNModules,0,fNModules);
3339  fhSumClustersEnergyMod->SetXTitle("#Sigma_{clusters} #it{E} (GeV)");
3340  fhSumClustersEnergyMod->SetYTitle("Module");
3341  outputContainer->Add(fhSumClustersEnergyMod);
3342 
3343  fhSumCellsAmpMod = new TH2F ("hSumCellsAmp_Mod","# cells vs Module", 1000, 0, 2000,fNModules,0,fNModules);
3344  fhSumCellsAmpMod->SetXTitle("#Sigma_{cells} #it{Amp} (GeV)");
3345  fhSumCellsAmpMod->SetYTitle("Module");
3346  outputContainer->Add(fhSumCellsAmpMod);
3347 
3348  Int_t colmaxs = fNMaxCols;
3349  Int_t rowmaxs = fNMaxRows;
3350  if(GetCalorimeter()==kEMCAL)
3351  {
3352  colmaxs=2*fNMaxCols;
3353  rowmaxs=Int_t(fNModules/2)*fNMaxRows;
3354  }
3355  else
3356  {
3357  rowmaxs=fNModules*fNMaxRows;
3358  }
3359 
3360  fhGridCells = new TH2F ("hGridCells",Form("Entries in grid of cells"),
3361  colmaxs+2,-1.5,colmaxs+0.5, rowmaxs+2,-1.5,rowmaxs+0.5);
3362  fhGridCells->SetYTitle("row (phi direction)");
3363  fhGridCells->SetXTitle("column (eta direction)");
3364  outputContainer->Add(fhGridCells);
3365 
3366  fhGridCellsE = new TH2F ("hGridCellsE","Accumulated energy in grid of cells",
3367  colmaxs+2,-1.5,colmaxs+0.5, rowmaxs+2,-1.5,rowmaxs+0.5);
3368  fhGridCellsE->SetYTitle("row (phi direction)");
3369  fhGridCellsE->SetXTitle("column (eta direction)");
3370  outputContainer->Add(fhGridCellsE);
3371 
3372  fhGridCellsLowGain = new TH2F ("hGridCellsLG",Form("Low gain: Entries in grid of cells"),
3373  colmaxs+2,-1.5,colmaxs+0.5, rowmaxs+2,-1.5,rowmaxs+0.5);
3374  fhGridCellsLowGain->SetYTitle("row (phi direction)");
3375  fhGridCellsLowGain->SetXTitle("column (eta direction)");
3376  outputContainer->Add(fhGridCellsLowGain);
3377 
3378  fhGridCellsELowGain = new TH2F ("hGridCellsELG","Low gain: Accumulated energy in grid of cells",
3379  colmaxs+2,-1.5,colmaxs+0.5, rowmaxs+2,-1.5,rowmaxs+0.5);
3380  fhGridCellsELowGain->SetYTitle("row (phi direction)");
3381  fhGridCellsELowGain->SetXTitle("column (eta direction)");
3382  outputContainer->Add(fhGridCellsELowGain);
3383 
3384  if(fFillAllCellTimeHisto)
3385  {
3386  fhGridCellsTime = new TH2F ("hGridCellsTime","Accumulated time in grid of cells",
3387  colmaxs+2,-1.5,colmaxs+0.5, rowmaxs+2,-1.5,rowmaxs+0.5);
3388  fhGridCellsTime->SetYTitle("row (phi direction)");
3389  fhGridCellsTime->SetXTitle("column (eta direction)");
3390  outputContainer->Add(fhGridCellsTime);
3391 
3392  fhGridCellsTimeLowGain = new TH2F ("hGridCellsTimeLG","Low gain: Accumulated time in grid of cells",
3393  colmaxs+2,-1.5,colmaxs+0.5, rowmaxs+2,-1.5,rowmaxs+0.5);
3394  fhGridCellsTimeLowGain->SetYTitle("row (phi direction)");
3395  fhGridCellsTimeLowGain->SetXTitle("column (eta direction)");
3396  outputContainer->Add(fhGridCellsTimeLowGain);
3397  }
3398 
3399  fhNCellsSumAmpPerMod = new TH2F*[fNModules];
3400 
3401  fhNClustersSumEnergyPerMod = new TH2F*[fNModules];
3402 
3403  fhIMMod = new TH2F*[fNModules];
3404 
3405  fhNCellsPerClusterMod = new TH2F*[fNModules];
3406 
3407  if(fStudyBadClusters)
3408  fhNCellsPerClusterModNoCut = new TH2F*[fNModules];
3409 
3410  if(fFillAllCellTimeHisto)
3411  fhTimeAmpPerRCU = new TH2F*[fNModules*fNRCU];
3412 
3413  for(Int_t imod = 0; imod < fNModules; imod++)
3414  {
3415  fhNCellsSumAmpPerMod[imod] = new TH2F (Form("hNCellsSumAmp_Mod%d",imod),
3416  Form("# cells in SM vs sum of cells energy in Module %d",imod),
3417  nptbins,ptmin,ptmax*4, ncebins,ncemin,ncemax);
3418  fhNCellsSumAmpPerMod[imod]->SetXTitle("#Sigma #it{Amplitude} (GeV)");
3419  fhNCellsSumAmpPerMod[imod]->SetYTitle("#Sigma #it{n}_{cells}");
3420  outputContainer->Add(fhNCellsSumAmpPerMod[imod]);
3421 
3422  fhNClustersSumEnergyPerMod[imod] = new TH2F (Form("hNClustersSumEnergy_Mod%d",imod),
3423  Form("# clusters in SM vs sum of clusters energy in Module %d",imod),
3424  nptbins,ptmin,ptmax*4, nclbins,nclmin,nclmax);
3425  fhNClustersSumEnergyPerMod[imod]->SetXTitle("#Sigma #it{E} (GeV)");
3426  fhNClustersSumEnergyPerMod[imod]->SetYTitle("#Sigma #it{n}_{clusters}");
3427  outputContainer->Add(fhNClustersSumEnergyPerMod[imod]);
3428 
3429  fhNCellsPerClusterMod[imod] = new TH2F (Form("hNCellsPerCluster_Mod%d",imod),
3430  Form("# cells per cluster vs cluster energy in Module %d",imod),
3431  nptbins,ptmin,ptmax, nceclbins,nceclmin,nceclmax);
3432  fhNCellsPerClusterMod[imod]->SetXTitle("#it{E} (GeV)");
3433  fhNCellsPerClusterMod[imod]->SetYTitle("#it{n}_{cells per cluster}");
3434  outputContainer->Add(fhNCellsPerClusterMod[imod]);
3435 
3436  if(fStudyBadClusters)
3437  {
3438  fhNCellsPerClusterModNoCut[imod] = new TH2F (Form("hNCellsPerClusterNoCut_Mod%d",imod),
3439  Form("# cells per cluster vs cluster energy in Module %d, no cut",imod),
3440  nptbins,ptmin,ptmax, nceclbins,nceclmin,nceclmax);
3441  fhNCellsPerClusterModNoCut[imod]->SetXTitle("#it{E} (GeV)");
3442  fhNCellsPerClusterModNoCut[imod]->SetYTitle("#it{n}_{cells per cluster}");
3443  outputContainer->Add(fhNCellsPerClusterModNoCut[imod]);
3444  }
3445 
3446  if(fFillAllCellTimeHisto)
3447  {
3448  for(Int_t ircu = 0; ircu < fNRCU; ircu++)
3449  {
3450  fhTimeAmpPerRCU[imod*fNRCU+ircu] = new TH2F (Form("hTimeAmp_Mod%d_RCU%d",imod,ircu),
3451  Form("#it{E}_{cell} vs #it{t}_{cell} in Module %d, RCU %d ",imod,ircu),
3452  nptbins,ptmin,ptmax,ntimebins,timemin,timemax);
3453  fhTimeAmpPerRCU[imod*fNRCU+ircu]->SetXTitle("#it{E} (GeV)");
3454  fhTimeAmpPerRCU[imod*fNRCU+ircu]->SetYTitle("#it{t} (ns)");
3455  outputContainer->Add(fhTimeAmpPerRCU[imod*fNRCU+ircu]);
3456  }
3457  }
3458 
3459  if(fFillAllPi0Histo)
3460  {
3461  fhIMMod[imod] = new TH2F (Form("hIM_Mod%d",imod),
3462  Form("Cluster pairs Invariant mass vs reconstructed pair energy in Module %d, n cell > 1",imod),
3463  nptbins,ptmin,ptmax,nmassbins,massmin,massmax);
3464  fhIMMod[imod]->SetXTitle("#it{p}_{T, cluster pairs} (GeV) ");
3465  fhIMMod[imod]->SetYTitle("#it{M}_{cluster pairs} (GeV/#it{c}^{2})");
3466  outputContainer->Add(fhIMMod[imod]);
3467  }
3468  }
3469 
3470  // Monte Carlo Histograms
3471 
3472  TString particleName[] = {
3473  "Photon", "Pi0", "Eta",
3474  "Electron", "PhotonConv",
3475  "NeutralHadron", "ChargedHadron" };
3476 
3477  if(IsDataMC())
3478  {
3479  for(Int_t iPart = 0; iPart < 7; iPart++)
3480  {
3481  for(Int_t iCh = 0; iCh < 2; iCh++)
3482  {
3483  fhRecoMCRatioE[iPart][iCh] = new TH2F (Form("hRecoMCRatioE_%s_Match%d",particleName[iPart].Data(),iCh),
3484  Form("Reconstructed/Generated E, %s, Matched %d",particleName[iPart].Data(),iCh),
3485  nptbins, ptmin, ptmax, 200,0,2);
3486  fhRecoMCRatioE[iPart][iCh]->SetYTitle("#it{E}_{reconstructed}/#it{E}_{generated}");
3487  fhRecoMCRatioE[iPart][iCh]->SetXTitle("#it{E}_{reconstructed} (GeV)");
3488  outputContainer->Add(fhRecoMCRatioE[iPart][iCh]);
3489 
3490 
3491  fhRecoMCDeltaE[iPart][iCh] = new TH2F (Form("hRecoMCDeltaE_%s_Match%d",particleName[iPart].Data(),iCh),
3492  Form("Generated - Reconstructed E, %s, Matched %d",particleName[iPart].Data(),iCh),
3493  nptbins, ptmin, ptmax, nptbins*2,-ptmax,ptmax);
3494  fhRecoMCDeltaE[iPart][iCh]->SetYTitle("#Delta #it{E} (GeV)");
3495  fhRecoMCDeltaE[iPart][iCh]->SetXTitle("#it{E}_{reconstructed} (GeV)");
3496  outputContainer->Add(fhRecoMCDeltaE[iPart][iCh]);
3497 
3498  fhRecoMCDeltaPhi[iPart][iCh] = new TH2F (Form("hRecoMCDeltaPhi_%s_Match%d",particleName[iPart].Data(),iCh),
3499  Form("Generated - Reconstructed #phi, %s, Matched %d",particleName[iPart].Data(),iCh),
3500  nptbins, ptmin, ptmax, nphibins*2,-phimax,phimax);
3501  fhRecoMCDeltaPhi[iPart][iCh]->SetYTitle("#Delta #phi (rad)");
3502  fhRecoMCDeltaPhi[iPart][iCh]->SetXTitle("#it{E}_{reconstructed} (GeV)");
3503  outputContainer->Add(fhRecoMCDeltaPhi[iPart][iCh]);
3504 
3505  fhRecoMCDeltaEta[iPart][iCh] = new TH2F (Form("hRecoMCDeltaEta_%s_Match%d",particleName[iPart].Data(),iCh),
3506  Form("Generated - Reconstructed #eta, %s, Matched %d",particleName[iPart].Data(),iCh),
3507  nptbins, ptmin, ptmax,netabins*2,-etamax,etamax);
3508  fhRecoMCDeltaEta[iPart][iCh]->SetYTitle("#Delta #eta ");
3509  fhRecoMCDeltaEta[iPart][iCh]->SetXTitle("#it{E}_{reconstructed} (GeV)");
3510  outputContainer->Add(fhRecoMCDeltaEta[iPart][iCh]);
3511 
3512  fhRecoMCE[iPart][iCh] = new TH2F (Form("hRecoMCE_%s_Match%d",particleName[iPart].Data(),iCh),
3513  Form("#it{E} distribution, reconstructed vs generated, %s, Matched %d",particleName[iPart].Data(),iCh),
3514  nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
3515  fhRecoMCE[iPart][iCh]->SetXTitle("#it{E}_{rec} (GeV)");
3516  fhRecoMCE[iPart][iCh]->SetYTitle("#it{E}_{gen} (GeV)");
3517  outputContainer->Add(fhRecoMCE[iPart][iCh]);
3518 
3519  fhRecoMCPhi[iPart][iCh] = new TH2F (Form("hRecoMCPhi_%s_Match%d",particleName[iPart].Data(),iCh),
3520  Form("#phi distribution, reconstructed vs generated, %s, Matched %d",particleName[iPart].Data(),iCh),
3521  nphibins,phimin,phimax, nphibins,phimin,phimax);
3522  fhRecoMCPhi[iPart][iCh]->SetXTitle("#phi_{reconstructed} (rad)");
3523  fhRecoMCPhi[iPart][iCh]->SetYTitle("#phi_{generated} (rad)");
3524  outputContainer->Add(fhRecoMCPhi[iPart][iCh]);
3525 
3526  fhRecoMCEta[iPart][iCh] = new TH2F (Form("hRecoMCEta_%s_Match%d",particleName[iPart].Data(),iCh),
3527  Form("#eta distribution, reconstructed vs generated, %s, Matched %d",particleName[iPart].Data(),iCh),
3528  netabins,etamin,etamax,netabins,etamin,etamax);
3529  fhRecoMCEta[iPart][iCh]->SetXTitle("#eta_{reconstructed} ");
3530  fhRecoMCEta[iPart][iCh]->SetYTitle("#eta_{generated} ");
3531  outputContainer->Add(fhRecoMCEta[iPart][iCh]);
3532  }
3533  }
3534 
3535  // Pure MC
3536 
3537  for(Int_t iPart = 0; iPart < 4; iPart++)
3538  {
3539  fhGenMCE [iPart] = new TH1F(Form("hGenMCE_%s",particleName[iPart].Data()) ,
3540  Form("#it{E} of generated %s",particleName[iPart].Data()),
3541  nptbins,ptmin,ptmax);
3542 
3543  fhGenMCPt[iPart] = new TH1F(Form("hGenMCPt_%s",particleName[iPart].Data()) ,
3544  Form("#it{p}_{T} of generated %s",particleName[iPart].Data()),
3545  nptbins,ptmin,ptmax);
3546 
3547  fhGenMCEtaPhi[iPart] = new TH2F(Form("hGenMCEtaPhi_%s",particleName[iPart].Data()),
3548  Form("Y vs #phi of generated %s",particleName[iPart].Data()),
3549  200,-1,1,360,0,TMath::TwoPi());
3550 
3551  fhGenMCE [iPart] ->SetXTitle("#it{E} (GeV)");
3552  fhGenMCPt[iPart] ->SetXTitle("#it{p}_{T} (GeV/#it{c})");
3553  fhGenMCEtaPhi[iPart]->SetXTitle("#eta");
3554  fhGenMCEtaPhi[iPart]->SetYTitle("#phi (rad)");
3555 
3556  outputContainer->Add(fhGenMCE [iPart]);
3557  outputContainer->Add(fhGenMCPt [iPart]);
3558  outputContainer->Add(fhGenMCEtaPhi[iPart]);
3559 
3560 
3561  fhGenMCAccE [iPart] = new TH1F(Form("hGenMCAccE_%s",particleName[iPart].Data()) ,
3562  Form("#it{E} of generated %s",particleName[iPart].Data()),
3563  nptbins,ptmin,ptmax);
3564  fhGenMCAccPt[iPart] = new TH1F(Form("hGenMCAccPt_%s",particleName[iPart].Data()) ,
3565  Form("#it{p}_{T} of generated %s",particleName[iPart].Data()),
3566  nptbins,ptmin,ptmax);
3567  fhGenMCAccEtaPhi[iPart] = new TH2F(Form("hGenMCAccEtaPhi_%s",particleName[iPart].Data()),
3568  Form("Y vs #phi of generated %s",particleName[iPart].Data()),
3569  netabins,etamin,etamax,nphibins,phimin,phimax);
3570 
3571  fhGenMCAccE [iPart] ->SetXTitle("#it{E} (GeV)");
3572  fhGenMCAccPt[iPart] ->SetXTitle("#it{p}_{T} (GeV/#it{c})");
3573  fhGenMCAccEtaPhi[iPart]->SetXTitle("#eta");
3574  fhGenMCAccEtaPhi[iPart]->SetYTitle("#phi (rad)");
3575 
3576  outputContainer->Add(fhGenMCAccE [iPart]);
3577  outputContainer->Add(fhGenMCAccPt [iPart]);
3578  outputContainer->Add(fhGenMCAccEtaPhi[iPart]);
3579 
3580  }
3581 
3582  // Vertex of generated particles
3583 
3584  fhEMVxyz = new TH2F ("hEMVxyz","Production vertex of reconstructed ElectroMagnetic particles",nvdistbins,vdistmin,vdistmax,nvdistbins,vdistmin,vdistmax);//,100,0,500);
3585  fhEMVxyz->SetXTitle("#it{v}_{x}");
3586  fhEMVxyz->SetYTitle("#it{v}_{y}");
3587  //fhEMVxyz->SetZTitle("v_{z}");
3588  outputContainer->Add(fhEMVxyz);
3589 
3590  fhHaVxyz = new TH2F ("hHaVxyz","Production vertex of reconstructed hadrons",nvdistbins,vdistmin,vdistmax,nvdistbins,vdistmin,vdistmax);//,100,0,500);
3591  fhHaVxyz->SetXTitle("#it{v}_{x}");
3592  fhHaVxyz->SetYTitle("#it{v}_{y}");
3593  //fhHaVxyz->SetZTitle("v_{z}");
3594  outputContainer->Add(fhHaVxyz);
3595 
3596  fhEMR = new TH2F ("hEMR","Distance to production vertex of reconstructed ElectroMagnetic particles vs E rec",nptbins,ptmin,ptmax,nvdistbins,vdistmin,vdistmax);
3597  fhEMR->SetXTitle("#it{E} (GeV)");
3598  fhEMR->SetYTitle("TMath::Sqrt(v_{x}^{2}+v_{y}^{2})");
3599  outputContainer->Add(fhEMR);
3600 
3601  fhHaR = new TH2F ("hHaR","Distance to production vertex of reconstructed Hadrons vs E rec",nptbins,ptmin,ptmax,nvdistbins,vdistmin,vdistmax);
3602  fhHaR->SetXTitle("#it{E} (GeV)");
3603  fhHaR->SetYTitle("TMath::Sqrt(v_{x}^{2}+v_{y}^{2})");
3604  outputContainer->Add(fhHaR);
3605 
3606  // Track Matching
3607 
3608  fhMCEle1EOverP = new TH2F("hMCEle1EOverP","TRACK matches #it{E}/#it{p}, MC electrons",nptbins,ptmin,ptmax, nPoverEbins,eOverPmin,eOverPmax);
3609  fhMCEle1EOverP->SetYTitle("#it{E}/#it{p}");
3610  fhMCEle1EOverP->SetXTitle("#it{p}_{T} (GeV/#it{c})");
3611  outputContainer->Add(fhMCEle1EOverP);
3612 
3613  fhMCEle1dR = new TH1F("hMCEle1dR","TRACK matches dR, MC electrons",ndRbins,dRmin,dRmax);
3614  fhMCEle1dR->SetXTitle("#Delta #it{R} (rad)");
3615  outputContainer->Add(fhMCEle1dR) ;
3616 
3617  fhMCEle2MatchdEdx = new TH2F("hMCEle2MatchdEdx","#it{dE/dx} vs. #it{p} for all matches, MC electrons",nptbins,ptmin,ptmax,ndedxbins,dedxmin,dedxmax);
3618  fhMCEle2MatchdEdx->SetXTitle("#it{p} (GeV/#it{c})");
3619  fhMCEle2MatchdEdx->SetYTitle("<#it{dE/dx}>");
3620  outputContainer->Add(fhMCEle2MatchdEdx);
3621 
3622  fhMCChHad1EOverP = new TH2F("hMCChHad1EOverP","TRACK matches #it{E}/#it{p}, MC charged hadrons",nptbins,ptmin,ptmax, nPoverEbins,eOverPmin,eOverPmax);
3623  fhMCChHad1EOverP->SetYTitle("#it{E}/#it{p}");
3624  fhMCChHad1EOverP->SetXTitle("#it{p}_{T} (GeV/#it{c})");
3625  outputContainer->Add(fhMCChHad1EOverP);
3626 
3627  fhMCChHad1dR = new TH1F("hMCChHad1dR","TRACK matches dR, MC charged hadrons",ndRbins,dRmin,dRmax);
3628  fhMCChHad1dR->SetXTitle("#Delta R (rad)");
3629  outputContainer->Add(fhMCChHad1dR) ;
3630 
3631  fhMCChHad2MatchdEdx = new TH2F("hMCChHad2MatchdEdx","#it{dE/dx} vs. #it{p} for all matches, MC charged hadrons",nptbins,ptmin,ptmax,ndedxbins,dedxmin,dedxmax);
3632  fhMCChHad2MatchdEdx->SetXTitle("#it{p} (GeV/#it{c})");
3633  fhMCChHad2MatchdEdx->SetYTitle("#it{dE/dx}>");
3634  outputContainer->Add(fhMCChHad2MatchdEdx);
3635 
3636  fhMCNeutral1EOverP = new TH2F("hMCNeutral1EOverP","TRACK matches #it{E}/#it{p}, MC neutrals",nptbins,ptmin,ptmax, nPoverEbins,eOverPmin,eOverPmax);
3637  fhMCNeutral1EOverP->SetYTitle("#it{E}/#it{p}");
3638  fhMCNeutral1EOverP->SetXTitle("#it{p}_{T} (GeV/#it{c})");
3639  outputContainer->Add(fhMCNeutral1EOverP);
3640 
3641  fhMCNeutral1dR = new TH1F("hMCNeutral1dR","TRACK matches dR, MC neutrals",ndRbins,dRmin,dRmax);
3642  fhMCNeutral1dR->SetXTitle("#Delta #it{R} (rad)");
3643  outputContainer->Add(fhMCNeutral1dR) ;
3644 
3645  fhMCNeutral2MatchdEdx = new TH2F("hMCNeutral2MatchdEdx","#it{dE/dx} vs. #it{p} for all matches, MC neutrals",nptbins,ptmin,ptmax,ndedxbins,dedxmin,dedxmax);
3646  fhMCNeutral2MatchdEdx->SetXTitle("#it{p} (GeV/#it{c})");
3647  fhMCNeutral2MatchdEdx->SetYTitle("#it{dE/dx}>");
3648  outputContainer->Add(fhMCNeutral2MatchdEdx);
3649 
3650  fhMCEle1EOverPR02 = new TH2F("hMCEle1EOverPR02","TRACK matches #it{E}/#it{p}, MC electrons",nptbins,ptmin,ptmax, nPoverEbins,eOverPmin,eOverPmax);
3651  fhMCEle1EOverPR02->SetYTitle("#it{E}/#it{p}");
3652  fhMCEle1EOverPR02->SetXTitle("#it{p}_{T} (GeV/#it{c})");
3653  outputContainer->Add(fhMCEle1EOverPR02);
3654 
3655  fhMCChHad1EOverPR02 = new TH2F("hMCChHad1EOverPR02","TRACK matches #it{E}/#it{p}, MC charged hadrons",nptbins,ptmin,ptmax, nPoverEbins,eOverPmin,eOverPmax);
3656  fhMCChHad1EOverPR02->SetYTitle("#it{E}/#it{p}");
3657  fhMCChHad1EOverPR02->SetXTitle("#it{p}_{T} (GeV/#it{c})");
3658  outputContainer->Add(fhMCChHad1EOverPR02);
3659 
3660  fhMCNeutral1EOverPR02 = new TH2F("hMCNeutral1EOverPR02","TRACK matches #it{E}/#it{p}, MC neutrals",nptbins,ptmin,ptmax, nPoverEbins,eOverPmin,eOverPmax);
3661  fhMCNeutral1EOverPR02->SetYTitle("#it{E}/#it{p}");
3662  fhMCNeutral1EOverPR02->SetXTitle("#it{p}_{T} (GeV/#it{c})");
3663  outputContainer->Add(fhMCNeutral1EOverPR02);
3664 
3665  fhMCEle1EleEOverP = new TH2F("hMCEle1EleEOverP","Electron candidates #it{E}/#it{p} (60<dEdx<100), MC electrons",nptbins,ptmin,ptmax, nPoverEbins,eOverPmin,eOverPmax);
3666  fhMCEle1EleEOverP->SetYTitle("#it{E}/#it{p}");
3667  fhMCEle1EleEOverP->SetXTitle("#it{p}_{T} (GeV/#it{c})");
3668  outputContainer->Add(fhMCEle1EleEOverP);
3669 
3670  fhMCChHad1EleEOverP = new TH2F("hMCEle1EleEOverP","Electron candidates #it{E}/#it{p} (60<dEdx<100), MC charged hadrons",nptbins,ptmin,ptmax, nPoverEbins,eOverPmin,eOverPmax);
3671  fhMCChHad1EleEOverP->SetYTitle("#it{E}/#it{p}");
3672  fhMCChHad1EleEOverP->SetXTitle("#it{p}_{T} (GeV/#it{c})");
3673  outputContainer->Add(fhMCChHad1EleEOverP);
3674 
3675  fhMCNeutral1EleEOverP = new TH2F("hMCNeutral1EleEOverP","Electron candidates #it{E}/#it{p} (60<dEdx<100), MC neutrals",nptbins,ptmin,ptmax, nPoverEbins,eOverPmin,eOverPmax);
3676  fhMCNeutral1EleEOverP->SetYTitle("#it{E}/#it{p}");
3677  fhMCNeutral1EleEOverP->SetXTitle("#it{p}_{T} (GeV/#it{c})");
3678  outputContainer->Add(fhMCNeutral1EleEOverP);
3679  }
3680 
3681 
3682  if(fFillEBinAcceptanceHisto)
3683  {
3684  for(Int_t ie=0; ie<fNEBinCuts; ie++)
3685  {
3686  fhEBinClusterEtaPhi[ie] = new TH2F
3687  (Form("hEBin%d_Cluster_EtaPhi",ie),
3688  Form("#eta vs #phi, cluster, %2.2f<#it{p}_{T}<%2.2f GeV/#it{c}",fEBinCuts[ie],fEBinCuts[ie+1]),
3689  netabins,etamin,etamax,nphibins,phimin,phimax);
3690  fhEBinClusterEtaPhi[ie]->SetYTitle("#phi (rad)");
3691  fhEBinClusterEtaPhi[ie]->SetXTitle("#eta");
3692  outputContainer->Add(fhEBinClusterEtaPhi[ie]) ;
3693 
3694  fhEBinClusterColRow[ie] = new TH2F
3695  (Form("hEBin%d_Cluster_ColRow",ie),
3696  Form("column vs row, cluster max E cell, %2.2f<#it{p}_{T}<%2.2f GeV/#it{c}",fEBinCuts[ie],fEBinCuts[ie+1]),
3697  96+2,-1.5,96+0.5,(8*24+2*8)+2,-1.5,(8*24+2*8)+0.5);
3698  fhEBinClusterColRow[ie]->SetYTitle("row");
3699  fhEBinClusterColRow[ie]->SetXTitle("column");
3700  outputContainer->Add(fhEBinClusterColRow[ie]) ;
3701 
3702  fhEBinCellColRow[ie] = new TH2F
3703  (Form("hEBin%d_Cell_ColRow",ie),
3704  Form("column vs row, cell, %2.2f<#it{p}_{T}<%2.2f GeV/#it{c}",fEBinCuts[ie],fEBinCuts[ie+1]),
3705  96+2,-1.5,96+0.5,(8*24+2*8)+2,-1.5,(8*24+2*8)+0.5);
3706  fhEBinCellColRow[ie]->SetYTitle("row");
3707  fhEBinCellColRow[ie]->SetXTitle("column");
3708  outputContainer->Add(fhEBinCellColRow[ie]) ;
3709  }
3710  }
3711 
3712 
3713  // for(Int_t i = 0; i < outputContainer->GetEntries() ; i++)
3714  // printf("i=%d, name= %s\n",i,outputContainer->At(i)->GetName());
3715 
3716  return outputContainer;
3717 }
3718 
3719 //______________________________________________________________________________________
3722 //______________________________________________________________________________________
3723 Float_t AliAnaCalorimeterQA::GetECross(Int_t absID, AliVCaloCells* cells, Float_t dtcut)
3724 {
3725  Int_t icol =-1, irow=-1,iRCU = -1;
3726  Int_t imod = GetModuleNumberCellIndexes(absID, GetCalorimeter(), icol, irow, iRCU);
3727 
3728  if ( GetCalorimeter() == kEMCAL)
3729  {
3730  // Get close cells index, energy and time, not in corners
3731 
3732  Int_t absID1 = -1;
3733  Int_t absID2 = -1;
3734 
3735  if( irow < AliEMCALGeoParams::fgkEMCALRows-1) absID1 = GetCaloUtils()->GetEMCALGeometry()->GetAbsCellIdFromCellIndexes(imod, irow+1, icol);
3736  if( irow > 0 ) absID2 = GetCaloUtils()->GetEMCALGeometry()->GetAbsCellIdFromCellIndexes(imod, irow-1, icol);
3737 
3738  // In case of cell in eta = 0 border, depending on SM shift the cross cell index
3739  Int_t absID3 = -1;
3740  Int_t absID4 = -1;
3741 
3742  if ( icol == AliEMCALGeoParams::fgkEMCALCols - 1 && !(imod%2) )
3743  {
3744  absID3 = GetCaloUtils()->GetEMCALGeometry()-> GetAbsCellIdFromCellIndexes(imod+1, irow, 0);
3745  absID4 = GetCaloUtils()->GetEMCALGeometry()-> GetAbsCellIdFromCellIndexes(imod , irow, icol-1);
3746  }
3747  else if( icol == 0 && imod%2 )
3748  {
3749  absID3 = GetCaloUtils()->GetEMCALGeometry()-> GetAbsCellIdFromCellIndexes(imod , irow, icol+1);
3750  absID4 = GetCaloUtils()->GetEMCALGeometry()-> GetAbsCellIdFromCellIndexes(imod-1, irow, AliEMCALGeoParams::fgkEMCALCols-1);
3751  }
3752  else
3753  {
3754  if( icol < AliEMCALGeoParams::fgkEMCALCols-1 )
3755  absID3 = GetCaloUtils()->GetEMCALGeometry()-> GetAbsCellIdFromCellIndexes(imod, irow, icol+1);
3756  if( icol > 0 )
3757  absID4 = GetCaloUtils()->GetEMCALGeometry()-> GetAbsCellIdFromCellIndexes(imod, irow, icol-1);
3758  }
3759 
3760  // Recalibrate cell energy if needed
3761  //Float_t ecell = cells->GetCellAmplitude(absID);
3762  //GetCaloUtils()->RecalibrateCellAmplitude(ecell,GetCalorimeter(), absID);
3763  Double_t tcell = cells->GetCellTime(absID);
3764  GetCaloUtils()->RecalibrateCellTime(tcell, GetCalorimeter(), absID,GetReader()->GetInputEvent()->GetBunchCrossNumber());
3765 
3766  Float_t ecell1 = 0, ecell2 = 0, ecell3 = 0, ecell4 = 0;
3767  Double_t tcell1 = 0, tcell2 = 0, tcell3 = 0, tcell4 = 0;
3768 
3769  if(absID1 > 0 )
3770  {
3771  ecell1 = cells->GetCellAmplitude(absID1);
3772  GetCaloUtils()->RecalibrateCellAmplitude(ecell1, GetCalorimeter(), absID1);
3773  tcell1 = cells->GetCellTime(absID1);
3774  GetCaloUtils()->RecalibrateCellTime (tcell1, GetCalorimeter(), absID1,GetReader()->GetInputEvent()->GetBunchCrossNumber());
3775  }
3776 
3777  if(absID2 > 0 )
3778  {
3779  ecell2 = cells->GetCellAmplitude(absID2);
3780  GetCaloUtils()->RecalibrateCellAmplitude(ecell2, GetCalorimeter(), absID2);
3781  tcell2 = cells->GetCellTime(absID2);
3782  GetCaloUtils()->RecalibrateCellTime (tcell2, GetCalorimeter(), absID2, GetReader()->GetInputEvent()->GetBunchCrossNumber());
3783  }
3784 
3785  if(absID3 > 0 )
3786  {
3787  ecell3 = cells->GetCellAmplitude(absID3);
3788  GetCaloUtils()->RecalibrateCellAmplitude(ecell3, GetCalorimeter(), absID3);
3789  tcell3 = cells->GetCellTime(absID3);
3790  GetCaloUtils()->RecalibrateCellTime (tcell3, GetCalorimeter(), absID3, GetReader()->GetInputEvent()->GetBunchCrossNumber());
3791  }
3792 
3793  if(absID4 > 0 )
3794  {
3795  ecell4 = cells->GetCellAmplitude(absID4);
3796  GetCaloUtils()->RecalibrateCellAmplitude(ecell4, GetCalorimeter(), absID4);
3797  tcell4 = cells->GetCellTime(absID4);
3798  GetCaloUtils()->RecalibrateCellTime (tcell4, GetCalorimeter(), absID4, GetReader()->GetInputEvent()->GetBunchCrossNumber());
3799  }
3800 
3801  if(TMath::Abs(tcell-tcell1)*1.e9 > dtcut) ecell1 = 0 ;
3802  if(TMath::Abs(tcell-tcell2)*1.e9 > dtcut) ecell2 = 0 ;
3803  if(TMath::Abs(tcell-tcell3)*1.e9 > dtcut) ecell3 = 0 ;
3804  if(TMath::Abs(tcell-tcell4)*1.e9 > dtcut) ecell4 = 0 ;
3805 
3806  return ecell1+ecell2+ecell3+ecell4;
3807  }
3808  else // PHOS
3809  {
3810  Int_t absId1 = -1, absId2 = -1, absId3 = -1, absId4 = -1;
3811 
3812  Int_t relId1[] = { imod+1, 0, irow+1, icol };
3813  Int_t relId2[] = { imod+1, 0, irow-1, icol };
3814  Int_t relId3[] = { imod+1, 0, irow , icol+1 };
3815  Int_t relId4[] = { imod+1, 0, irow , icol-1 };
3816 
3817  GetCaloUtils()->GetPHOSGeometry()->RelToAbsNumbering(relId1, absId1);
3818  GetCaloUtils()->GetPHOSGeometry()->RelToAbsNumbering(relId2, absId2);
3819  GetCaloUtils()->GetPHOSGeometry()->RelToAbsNumbering(relId3, absId3);
3820  GetCaloUtils()->GetPHOSGeometry()->RelToAbsNumbering(relId4, absId4);
3821 
3822  Float_t ecell1 = 0, ecell2 = 0, ecell3 = 0, ecell4 = 0;
3823 
3824  if(absId1 > 0 ) ecell1 = cells->GetCellAmplitude(absId1);
3825  if(absId2 > 0 ) ecell2 = cells->GetCellAmplitude(absId2);
3826  if(absId3 > 0 ) ecell3 = cells->GetCellAmplitude(absId3);
3827  if(absId4 > 0 ) ecell4 = cells->GetCellAmplitude(absId4);
3828 
3829  return ecell1+ecell2+ecell3+ecell4;
3830  }
3831 }
3832 
3833 //___________________________________________________________________________________________________________
3839 //___________________________________________________________________________________________________________
3840 void AliAnaCalorimeterQA::InvariantMassHistograms(Int_t iclus, Int_t nModule, const TObjArray* caloClusters,
3841  AliVCaloCells * cells)
3842 {
3843  AliDebug(1,"Start");
3844 
3845  //Get vertex for photon momentum calculation and event selection
3846  Double_t v[3] = {0,0,0}; //vertex ;
3847  //GetReader()->GetVertex(v);
3848 
3849  Int_t nModule2 = -1;
3850  Int_t nCaloClusters = caloClusters->GetEntriesFast();
3851 
3852  Float_t phi1 = fClusterMomentum.Phi();
3853  if(phi1 < 0) phi1 += TMath::TwoPi();
3854 
3855  Double_t tof1 = ((AliVCluster*) caloClusters->At(iclus))->GetTOF()*1.e9;
3856  if(tof1>400) tof1-=fConstantTimeShift;
3857 
3858  for(Int_t jclus = iclus + 1 ; jclus < nCaloClusters ; jclus++)
3859  {
3860  AliVCluster* clus2 = (AliVCluster*) caloClusters->At(jclus);
3861 
3862  Float_t maxCellFraction = 0.;
3863  Int_t absIdMax = GetCaloUtils()->GetMaxEnergyCell(cells, clus2, maxCellFraction);
3864 
3865  Double_t tof2 = clus2->GetTOF()*1.e9;
3866  if(tof2>400) tof2-=fConstantTimeShift;
3867 
3868  Double_t diffTof = tof1-tof2;
3869 
3870  // Try to reduce background with a mild shower shape cut and no more
3871  // than 1 local maximum in cluster and remove low energy clusters
3872 
3873  if( !IsGoodCluster(absIdMax, clus2->GetM02(), clus2->GetNCells(), cells)
3874  || GetCaloUtils()->GetNumberOfLocalMaxima(clus2,cells) > 1
3875  || clus2->GetM02() > fInvMassMaxM02Cut
3876  || clus2->GetM02() < fInvMassMinM02Cut
3877  || clus2->E() < fInvMassMinECut
3878  || clus2->E() > fInvMassMaxECut
3879  || TMath::Abs(diffTof) > fInvMassMaxTimeDifference
3880  ) continue;
3881 
3882  // Get cluster kinematics
3883  clus2->GetMomentum(fClusterMomentum2,v);
3884 
3885  // Check only certain regions
3886  Bool_t in2 = kTRUE;
3887  if(IsFiducialCutOn()) in2 = GetFiducialCut()->IsInFiducialCut(fClusterMomentum2.Eta(),fClusterMomentum2.Phi(),GetCalorimeter()) ;
3888  if(!in2) continue;
3889 
3890  Float_t pairPt = (fClusterMomentum+fClusterMomentum2).Pt();
3891 
3892  // Opening angle cut, avoid combination of DCal and EMCal clusters
3893  Double_t angle = fClusterMomentum.Angle(fClusterMomentum2.Vect());
3894 
3895  if ( fFillInvMassOpenAngle ) fhOpAngle->Fill(pairPt, angle, GetEventWeight()) ;
3896 
3897  if( angle > fInvMassMaxOpenAngle ) continue;
3898 
3899  // Get module of cluster
3900  nModule2 = GetModuleNumber(clus2);
3901 
3902  // Fill histograms
3903  Float_t mass = (fClusterMomentum+fClusterMomentum2).M ();
3904  Float_t asym = TMath::Abs( fClusterMomentum.E() - fClusterMomentum2.E() ) /
3905  ( fClusterMomentum.E() + fClusterMomentum2.E() );
3906 
3907  // All modules
3908  // Combine EMCal or PHOS clusters
3909 
3910  Float_t phi2 = fClusterMomentum2.Phi();
3911  if(phi2 < 0) phi2 += TMath::TwoPi();
3912 
3913  AliDebug(1,Form("Selected pair: pT %f, mass %f, asym %f, angle %f, diffTof %f, SM1 %d, SM2 %d\n",pairPt,mass,asym,angle,diffTof,nModule,nModule2));
3914 
3915  Bool_t inPi0Window = kFALSE;
3916  if(mass < 0.18 && mass > 0.1) inPi0Window = kTRUE ;
3917 
3918  if ( nModule < 12 && nModule2 < 12 )
3919  {
3920  fhIM ->Fill(pairPt, mass, GetEventWeight());
3921 
3922  if( fFillPi0PairDiffTime && inPi0Window )
3923  fhClusterPairDiffTimeEPi0Mass->Fill(pairPt, diffTof, GetEventWeight());
3924 
3925  if ( nModule == nModule2 )
3926  {
3927  fhIMSame->Fill(pairPt, mass, GetEventWeight());
3928 
3929  if( fFillPi0PairDiffTime && inPi0Window )
3930  fhClusterPairDiffTimeEPi0MassSame->Fill(pairPt, diffTof, GetEventWeight());
3931  }
3932  else
3933  {
3934  fhIMDiff->Fill(pairPt, mass, GetEventWeight());
3935  }
3936  }
3937  // Combine DCal clusters
3938  else if ( ( GetCalorimeter() == kEMCAL || GetCalorimeter() == kDCAL ) &&
3939  nModule > 11 && nModule2 > 11 && fNModules > 12 )
3940  {
3941  fhIMDCAL->Fill(pairPt, mass, GetEventWeight());
3942 
3943  if( fFillPi0PairDiffTime && inPi0Window )
3944  fhClusterPairDiffTimeEPi0MassDCal->Fill(pairPt, diffTof, GetEventWeight());
3945 
3946  if ( nModule == nModule2 )
3947  {
3948  fhIMDCALSame->Fill(pairPt, mass, GetEventWeight());
3949 
3950  if( fFillPi0PairDiffTime && inPi0Window )
3951  fhClusterPairDiffTimeEPi0MassDCalSame->Fill(pairPt, diffTof, GetEventWeight());
3952  }
3953  else
3954  {
3955  fhIMDCALDiff->Fill(pairPt, mass, GetEventWeight());
3956  }
3957  }
3958 
3959  if ( fFillInvMassOpenAngle ) fhIMvsOpAngle->Fill(mass, angle, GetEventWeight());
3960 
3961  // Single module
3962  if(nModule == nModule2 && nModule >= 0 && nModule < fNModules)
3963  fhIMMod[nModule]->Fill(pairPt, mass, GetEventWeight());
3964 
3965  // Asymetry histograms
3966  fhAsym->Fill(pairPt, asym, GetEventWeight());
3967  } // 2nd cluster loop
3968 
3969  //
3970  // Combine PHOS and DCal
3971  //
3972 
3973  if( ( GetCalorimeter() == kEMCAL || GetCalorimeter() == kDCAL ) &&
3974  fNModules > 12 && nModule > 11)
3975  {
3976  AliDebug(1,"Check DCal-PHOS pairs\n");
3977 
3978  Int_t sector1 = -1;
3979  Int_t sector2 = -1;
3980 
3981  if(phi1 >= 260*TMath::DegToRad() && phi1 < 280) sector1 = 0;
3982  if(phi1 >= 280*TMath::DegToRad() && phi1 < 300) sector1 = 1;
3983  if(phi1 >= 300*TMath::DegToRad() && phi1 < 320) sector1 = 2;
3984 
3985  for(Int_t jclus = 0 ; jclus < GetPHOSClusters()->GetEntriesFast() ; jclus++)
3986  {
3987  AliVCluster* clus2 = (AliVCluster*) GetPHOSClusters()->At(jclus);
3988 
3989  Float_t maxCellFraction = 0.;
3990  Int_t absIdMax = GetCaloUtils()->GetMaxEnergyCell(cells, clus2, maxCellFraction);
3991 
3992  // Try to reduce background, remove low energy clusters
3993  if( !IsGoodCluster(absIdMax, clus2->GetM02(), clus2->GetNCells(), cells)
3994  || clus2->E() < fInvMassMinECut
3995  || clus2->E() > fInvMassMaxECut
3996  ) continue;
3997 
3998  // Get cluster kinematics
3999  clus2->GetMomentum(fClusterMomentum2,v);
4000 
4001  // Fill histograms
4002 
4003  Float_t mass = (fClusterMomentum+fClusterMomentum2).M ();
4004  Float_t pairPt = (fClusterMomentum+fClusterMomentum2).Pt();
4005  //Float_t asym = TMath::Abs( fClusterMomentum.E() - fClusterMomentum2.E() ) /
4006  //( fClusterMomentum.E() + fClusterMomentum2.E() );
4007 
4008  fhIMDCALPHOS->Fill(pairPt, mass, GetEventWeight());
4009 
4010  Float_t phiPHOS = fClusterMomentum2.Phi();
4011  if(phiPHOS < 0) phiPHOS += TMath::TwoPi();
4012 
4013  if(phiPHOS >= 260*TMath::DegToRad() && phiPHOS < 280) sector2 = 0;
4014  if(phiPHOS >= 280*TMath::DegToRad() && phiPHOS < 300) sector2 = 1;
4015  if(phiPHOS >= 300*TMath::DegToRad() && phiPHOS < 320) sector2 = 2;
4016 
4017  if(sector1 == sector2) fhIMDCALPHOSSame->Fill(pairPt, mass, GetEventWeight());
4018 
4019  } // PHOS cluster loop
4020  } // DCal-PHOS combination
4021 
4022 
4023  // Select EMCal clusters in DCal eta acceptance
4024  // Cross check combination of DCal-PHOS pairs
4025  if( fFillInvMassInEMCALWithPHOSDCalAcc && TMath::Abs(fClusterMomentum.Eta() > 0.22))
4026  {
4027  AliDebug(1,"Check EMCAL(DCal)-EMCAL(PHOS) pairs\n");
4028 
4029  Int_t sector1 = -1;
4030  Int_t sector2 = -1;
4031 
4032  if(phi1 >= 80*TMath::DegToRad() && phi1 < 100) sector1 = 0;
4033  if(phi1 >= 100*TMath::DegToRad() && phi1 < 120) sector1 = 1;
4034  if(phi1 >= 120*TMath::DegToRad() && phi1 < 140) sector1 = 2;
4035  if(phi1 >= 140*TMath::DegToRad() && phi1 < 160) sector1 = 3;
4036  if(phi1 >= 160*TMath::DegToRad() && phi1 < 180) sector1 = 4;
4037  if(phi1 >= 180*TMath::DegToRad() && phi1 < 190) sector1 = 5;
4038 
4039  for(Int_t jclus = 0 ; jclus < caloClusters->GetEntriesFast() ; jclus++)
4040  {
4041  AliVCluster* clus2 = (AliVCluster*) caloClusters->At(jclus);
4042 
4043  Float_t maxCellFraction = 0.;
4044  Int_t absIdMax = GetCaloUtils()->GetMaxEnergyCell(cells, clus2, maxCellFraction);
4045 
4046  Double_t tof2 = clus2->GetTOF()*1.e9;
4047  if(tof2>400) tof2-=fConstantTimeShift;
4048 
4049  Double_t diffTof = TMath::Abs(tof1-tof2);
4050 
4051  // Try to reduce background with a mild shower shape cut and no more
4052  // than 1 local maximum in cluster and remove low energy clusters
4053  if( !IsGoodCluster(absIdMax, clus2->GetM02(), clus2->GetNCells(), cells)
4054  || GetCaloUtils()->GetNumberOfLocalMaxima(clus2,cells) > 1
4055  || clus2->GetM02() > fInvMassMaxM02Cut
4056  || clus2->GetM02() < fInvMassMinM02Cut
4057  || clus2->E() < fInvMassMinECut
4058  || clus2->E() > fInvMassMaxECut
4059  || TMath::Abs(diffTof) > fInvMassMaxTimeDifference
4060  ) continue;
4061 
4062  // Get cluster kinematics
4063  clus2->GetMomentum(fClusterMomentum2,v);
4064 
4065  // Select EMCal clusters in PHOS acceptance
4066  if(TMath::Abs(fClusterMomentum2.Eta() > 0.13)) continue ;
4067 
4068  // Fill histograms
4069 
4070  Float_t mass = (fClusterMomentum+fClusterMomentum2).M ();
4071  Float_t pairPt = (fClusterMomentum+fClusterMomentum2).Pt();
4072  //Float_t asym = TMath::Abs( fClusterMomentum.E() - fClusterMomentum2.E() ) /
4073  //( fClusterMomentum.E() + fClusterMomentum2.E() );
4074 
4075  fhIMEMCALPHOS->Fill(pairPt, mass, GetEventWeight());
4076 
4077  Float_t phiPHOS = fClusterMomentum2.Phi();
4078  if(phiPHOS < 0) phiPHOS += TMath::TwoPi();
4079 
4080  if(phiPHOS >= 80*TMath::DegToRad() && phiPHOS < 100) sector2 = 0;
4081  if(phiPHOS >= 100*TMath::DegToRad() && phiPHOS < 120) sector2 = 1;
4082  if(phiPHOS >= 120*TMath::DegToRad() && phiPHOS < 140) sector2 = 2;
4083  if(phiPHOS >= 140*TMath::DegToRad() && phiPHOS < 160) sector2 = 3;
4084  if(phiPHOS >= 160*TMath::DegToRad() && phiPHOS < 180) sector2 = 4;
4085  if(phiPHOS >= 180*TMath::DegToRad() && phiPHOS < 190) sector2 = 5;
4086 
4087  if(sector1 == sector2) fhIMEMCALPHOSSame->Fill(pairPt, mass, GetEventWeight());
4088 
4089  } // PHOS cluster loop
4090  } // EMCal(DCal)-EMCAL(PHOS) combination
4091 
4092  AliDebug(1,"End");
4093 }
4094 
4095 //______________________________
4096 // Check if the calorimeter setting is ok, if not abort.
4097 //______________________________
4098 void AliAnaCalorimeterQA::Init()
4099 {
4100  if(GetCalorimeter() != kPHOS && GetCalorimeter() !=kEMCAL)
4101  AliFatal(Form("Wrong calorimeter name <%s>", GetCalorimeterString().Data()));
4102 
4103  //if(GetReader()->GetDataType()== AliCaloTrackReader::kMC)
4104  // AliFatal("Analysis of reconstructed data, MC reader not aplicable");
4105 }
4106 
4107 //________________________________________
4109 //________________________________________
4110 void AliAnaCalorimeterQA::InitParameters()
4111 {
4112  AddToHistogramsName("AnaCaloQA_");
4113 
4114  fNModules = 22; // set maximum to maximum number of EMCAL modules
4115  fNRCU = 2; // set maximum number of RCU in EMCAL per SM
4116 
4117  fTimeCutMin = -9999999;
4118  fTimeCutMax = 9999999;
4119 
4120  fEMCALCellAmpMin = 0.2; // 200 MeV
4121  fPHOSCellAmpMin = 0.2; // 200 MeV
4122  fCellAmpMin = 0.2; // 200 MeV
4123 
4124  fEMCALClusterNCellMin = 2; // at least 2
4125  fPHOSClusterNCellMin = 3; // at least 3
4126 
4127  fInvMassMinECut = 0.5; // 500 MeV
4128  fInvMassMaxECut = 10;
4129 
4130  fInvMassMinM02Cut = 0.1;
4131  fInvMassMaxM02Cut = 0.4;
4132 
4133  fInvMassMaxTimeDifference = 200; // ns, quite open
4134 
4135  fInvMassMaxOpenAngle = 100*TMath::DegToRad(); // 100 degrees
4136 
4137  // Exotic studies
4138  fExoNECrossCuts = 10 ;
4139  fExoNDTimeCuts = 4 ;
4140 
4141  fExoDTimeCuts [0] = 1.e4 ; fExoDTimeCuts [1] = 50.0 ; fExoDTimeCuts [2] = 25.0 ; fExoDTimeCuts [3] = 10.0 ;
4142  fExoECrossCuts[0] = 0.80 ; fExoECrossCuts[1] = 0.85 ; fExoECrossCuts[2] = 0.90 ; fExoECrossCuts[3] = 0.92 ; fExoECrossCuts[4] = 0.94 ;
4143  fExoECrossCuts[5] = 0.95 ; fExoECrossCuts[6] = 0.96 ; fExoECrossCuts[7] = 0.97 ; fExoECrossCuts[8] = 0.98 ; fExoECrossCuts[9] = 0.99 ;
4144 
4145  fNEBinCuts = 14;
4146  fEBinCuts[0] = 0.; fEBinCuts[1] = 0.3; fEBinCuts[2] = 0.5;
4147  fEBinCuts[3] = 1.; fEBinCuts[4] = 2. ; fEBinCuts[5] = 3. ;
4148  fEBinCuts[6] = 4.; fEBinCuts[7] = 5. ; fEBinCuts[8] = 7. ;
4149  fEBinCuts[9] = 9.; fEBinCuts[10]= 12.; fEBinCuts[11]= 15.;
4150  fEBinCuts[12]= 20.; fEBinCuts[13]= 50.; fEBinCuts[14]= 100.;
4151  for(Int_t i = fNEBinCuts; i < 15; i++) fEBinCuts[i] = 1000.;
4152 
4153 
4154 }
4155 
4156 //_____________________________________________________________________________
4158 //_____________________________________________________________________________
4159 Bool_t AliAnaCalorimeterQA::IsGoodCluster(Int_t absIdMax, Float_t m02,
4160  Int_t nCellsPerCluster, AliVCaloCells* cells)
4161 {
4162  //if(!fStudyBadClusters) return kTRUE;
4163 
4164  if(GetCalorimeter() == kEMCAL)
4165  {
4166 
4167  if( m02 < 0.05 || nCellsPerCluster < fEMCALClusterNCellMin ) return kFALSE ; // mild shower shape cut for exotics
4168 
4169  if(!GetCaloUtils()->GetEMCALRecoUtils()->IsRejectExoticCluster())
4170  {
4171  return !(GetCaloUtils()->GetEMCALRecoUtils()->IsExoticCell(absIdMax,cells,(GetReader()->GetInputEvent())->GetBunchCrossNumber()));
4172  }
4173  else
4174  {
4175  return kTRUE;
4176  }
4177  }
4178  else // PHOS
4179  {
4180  if( m02 < 0.05 || nCellsPerCluster < fPHOSClusterNCellMin ) return kFALSE ; // mild shower shape cut for exotics
4181 
4182  Float_t ampMax = cells->GetCellAmplitude(absIdMax);
4183  GetCaloUtils()->RecalibrateCellAmplitude(ampMax, GetCalorimeter(), absIdMax);
4184 
4185  if(ampMax < 0.01) return kFALSE;
4186 
4187  if(1-GetECross(absIdMax,cells)/ampMax > 0.95) return kFALSE;
4188  else return kTRUE;
4189  }
4190 }
4191 
4192 //_________________________________________________________
4194 //_________________________________________________________
4195 void AliAnaCalorimeterQA::Print(const Option_t * opt) const
4196 {
4197  if(! opt)
4198  return;
4199 
4200  printf("**** Print %s %s ****\n", GetName(), GetTitle() ) ;
4201  AliAnaCaloTrackCorrBaseClass::Print(" ");
4202 
4203  printf("Select Calorimeter %s \n",GetCalorimeterString().Data());
4204  printf("Time Cut: %3.1f < TOF < %3.1f\n" , fTimeCutMin, fTimeCutMax);
4205  printf("EMCAL Min Amplitude : %2.1f GeV/c\n", fEMCALCellAmpMin) ;
4206  printf("PHOS Min Amplitude : %2.1f GeV/c\n", fPHOSCellAmpMin) ;
4207  printf("EMCAL Min n cells : %d\n" , fEMCALClusterNCellMin) ;
4208  printf("PHOS Min n cells : %d\n" , fPHOSClusterNCellMin) ;
4209 
4210  printf("Inv. Mass %2.1f < E_clus < %2.1f GeV/c\n" , fInvMassMinECut , fInvMassMaxECut ) ;
4211  printf("Inv. Mass %2.1f < M02_clus < %2.1f GeV/c\n", fInvMassMinM02Cut, fInvMassMaxM02Cut) ;
4212  printf("Inv. Mass open angle : %2.1f deg\n" , fInvMassMaxOpenAngle*TMath::RadToDeg()) ;
4213  printf("Inv. Mass time difference: %2.1f ns\n" , fInvMassMaxTimeDifference) ;
4214 }
4215 
4216 //_____________________________________________________
4218 //_____________________________________________________
4219 void AliAnaCalorimeterQA::MakeAnalysisFillHistograms()
4220 {
4221  AliDebug(1,"Start");
4222 
4223  //Print("");
4224 
4225  // Play with the MC stack if available
4226  if(IsDataMC()) MCHistograms();
4227 
4228  // Correlate Calorimeters and V0 and track Multiplicity
4229  if(fCorrelate) Correlate();
4230 
4231  // Get List with CaloClusters , calo Cells, init min amplitude
4232  TObjArray * caloClusters = NULL;
4233  AliVCaloCells * cells = 0x0;
4234  if (GetCalorimeter() == kPHOS)
4235  {
4236  fCellAmpMin = fPHOSCellAmpMin;
4237  caloClusters = GetPHOSClusters();
4238  cells = GetPHOSCells();
4239  }
4240  else if (GetCalorimeter() == kEMCAL)
4241  {
4242  fCellAmpMin = fEMCALCellAmpMin;
4243  caloClusters = GetEMCALClusters();
4244  cells = GetEMCALCells();
4245  }
4246  else
4247  AliFatal(Form("AliAnaCalorimeterQA::MakeAnalysisFillHistograms() - Wrong calorimeter name <%s>, END", GetCalorimeterString().Data()));
4248 
4249  if( !caloClusters || !cells )
4250  {
4251  AliFatal(Form("AliAnaCalorimeterQA::MakeAnalysisFillHistograms() - No CaloClusters or CaloCells available"));
4252  return; // trick coverity
4253  }
4254 
4255  if(caloClusters->GetEntriesFast() == 0) return ;
4256 
4257  //printf("QA: N cells %d, N clusters %d \n",cells->GetNumberOfCells(),caloClusters->GetEntriesFast());
4258 
4259  // Clusters
4260  ClusterLoopHistograms(caloClusters,cells);
4261 
4262  // Cells
4263  CellHistograms(cells);
4264 
4265  AliDebug(1,"End");
4266 }
4267 
4268 //______________________________________
4272 //______________________________________
4273 void AliAnaCalorimeterQA::MCHistograms()
4274 {
4275  Int_t pdg = 0 ;
4276  Int_t status = 0 ;
4277  Int_t nprim = 0 ;
4278 
4279  TParticle * primStack = 0;
4280  AliAODMCParticle * primAOD = 0;
4281 
4282  // Get the MC arrays and do some checks before filling MC histograms
4283 
4284  // Get the ESD MC particles container
4285  AliStack * stack = 0;
4286  if( GetReader()->ReadStack() )
4287  {
4288  stack = GetMCStack();
4289  if(!stack ) return;
4290  nprim = stack->GetNtrack();
4291  }
4292 
4293  // Get the AOD MC particles container
4294  TClonesArray * mcparticles = 0;
4295  if( GetReader()->ReadAODMCParticles() )
4296  {
4297  mcparticles = GetReader()->GetAODMCParticles();
4298  if( !mcparticles ) return;
4299  nprim = mcparticles->GetEntriesFast();
4300  }
4301 
4302  //printf("N primaries %d\n",nprim);
4303  for(Int_t i=0 ; i < nprim; i++)
4304  {
4305  if ( !GetReader()->AcceptParticleMCLabel( i ) ) continue ;
4306 
4307  // Get the generated particles, check that it is primary (not parton, resonance)
4308  // and get its momentum. Different way to recover from ESD or AOD
4309  if(GetReader()->ReadStack())
4310  {
4311  primStack = stack->Particle(i) ;
4312  if(!primStack)
4313  {
4314  AliWarning("ESD primaries pointer not available!!");
4315  continue;
4316  }
4317 
4318  pdg = primStack->GetPdgCode();
4319  status = primStack->GetStatusCode();
4320 
4321  //printf("Input: i %d, %s, pdg %d, status %d \n",i, primStack->GetName(), pdg, status);
4322 
4323  if ( status > 11 ) continue; //Working for PYTHIA and simple generators, check for HERWIG, HIJING?
4324 
4325  // Protection against floating point exception
4326  if ( primStack->Energy() == TMath::Abs(primStack->Pz()) ||
4327  (primStack->Energy() - primStack->Pz()) < 1e-3 ||
4328  (primStack->Energy() + primStack->Pz()) < 0 ) continue ;
4329 
4330  //printf("Take : i %d, %s, pdg %d, status %d \n",i, primStack->GetName(), pdg, status);
4331 
4332  // Photon kinematics
4333  primStack->Momentum(fPrimaryMomentum);
4334  }
4335  else
4336  {
4337  primAOD = (AliAODMCParticle *) mcparticles->At(i);
4338  if(!primAOD)
4339  {
4340  AliWarning("AOD primaries pointer not available!!");
4341  continue;
4342  }
4343 
4344  pdg = primAOD->GetPdgCode();
4345  status = primAOD->GetStatus();
4346 
4347  //printf("Input: i %d, %s, pdg %d, status %d \n",i, primAOD->GetName(), pdg, status);
4348 
4349  if (!primAOD->IsPrimary()) continue; //accept all which is not MC transport generated. Don't know how to avoid partons
4350 
4351  if ( status > 11 ) continue; //Working for PYTHIA and simple generators, check for HERWIG
4352 
4353  //Protection against floating point exception
4354  if ( primAOD->E() == TMath::Abs(primAOD->Pz()) ||
4355  (primAOD->E() - primAOD->Pz()) < 1e-3 ||
4356  (primAOD->E() + primAOD->Pz()) < 0 ) continue ;
4357 
4358  //printf("Take : i %d, %s, pdg %d, status %d \n",i, primAOD->GetName(), pdg, status);
4359 
4360  // Kinematics
4361  fPrimaryMomentum.SetPxPyPzE(primAOD->Px(),primAOD->Py(),primAOD->Pz(),primAOD->E());
4362  }
4363 
4364  Float_t eMC = fPrimaryMomentum.E();
4365  if(eMC < 0.2) continue;
4366 
4367  Float_t ptMC = fPrimaryMomentum.E();
4368 
4369  Float_t etaMC = fPrimaryMomentum.Eta();
4370 
4371  // Only particles in |eta| < 1
4372  if (TMath::Abs(etaMC) > 1) continue;
4373 
4374  Float_t phiMC = fPrimaryMomentum.Phi();
4375  if(phiMC < 0)
4376  phiMC += TMath::TwoPi();
4377 
4378  Int_t mcIndex = -1;
4379  if (pdg==22) mcIndex = kmcPhoton;
4380  else if (pdg==111) mcIndex = kmcPi0;
4381  else if (pdg==221) mcIndex = kmcEta;
4382  else if (TMath::Abs(pdg)==11) mcIndex = kmcElectron;
4383 
4384  if( mcIndex >=0 )
4385  {
4386  fhGenMCE [mcIndex]->Fill( eMC, GetEventWeight());
4387  fhGenMCPt[mcIndex]->Fill(ptMC, GetEventWeight());
4388  if(eMC > 0.5) fhGenMCEtaPhi[mcIndex]->Fill(etaMC, phiMC, GetEventWeight());
4389 
4390  Bool_t inacceptance = kTRUE;
4391  // Check same fidutial borders as in data analysis on top of real acceptance if real was requested.
4392  if( IsFiducialCutOn() && !GetFiducialCut()->IsInFiducialCut(fPrimaryMomentum.Eta(),fPrimaryMomentum.Phi(),GetCalorimeter()) )
4393  inacceptance = kFALSE ;
4394 
4395  if(IsRealCaloAcceptanceOn()) // defined on base class
4396  {
4397  if(GetReader()->ReadStack() &&
4398  !GetCaloUtils()->IsMCParticleInCalorimeterAcceptance(GetCalorimeter(), primStack)) inacceptance = kFALSE ;
4399  if(GetReader()->ReadAODMCParticles() &&
4400  !GetCaloUtils()->IsMCParticleInCalorimeterAcceptance(GetCalorimeter(), primAOD )) inacceptance = kFALSE ;
4401  }
4402 
4403  if(!inacceptance) continue;
4404 
4405  fhGenMCAccE [mcIndex]->Fill( eMC, GetEventWeight());
4406  fhGenMCAccPt[mcIndex]->Fill(ptMC, GetEventWeight());
4407  if(eMC > 0.5) fhGenMCAccEtaPhi[mcIndex]->Fill(etaMC, phiMC, GetEventWeight());
4408  }
4409  }
4410 }
4411 
4412 //_________________________________________________________________________________
4415 //_________________________________________________________________________________
4416 void AliAnaCalorimeterQA::WeightHistograms(AliVCluster *clus, AliVCaloCells* cells)
4417 {
4418  // First recalculate energy in case non linearity was applied
4419  Float_t energy = 0;
4420  Float_t ampMax = 0;
4421  Int_t bc = GetReader()->GetInputEvent()->GetBunchCrossNumber();
4422  Float_t enOrg = clus->E();
4423 
4424  // Do study when there are enough cells in cluster
4425  if(clus->GetNCells() < 3) return ;
4426 
4427  for (Int_t ipos = 0; ipos < clus->GetNCells(); ipos++)
4428  {
4429  Int_t id = clus->GetCellsAbsId()[ipos];
4430 
4431  //Recalibrate cell energy if needed
4432  Float_t amp = cells->GetCellAmplitude(id);
4433  GetCaloUtils()->RecalibrateCellAmplitude(amp, GetCalorimeter(), id);
4434 
4435  energy += amp;
4436 
4437  if ( amp > ampMax )
4438  ampMax = amp;
4439 
4440  } // energy loop
4441 
4442  if ( energy <=0 )
4443  {
4444  AliWarning(Form("Wrong calculated energy %f",energy));
4445  return;
4446  }
4447 
4448  // Remove non linearity correction
4449  clus->SetE(energy);
4450 
4451  fhEMaxCellClusterRatio ->Fill(energy, ampMax/energy , GetEventWeight());
4452  fhEMaxCellClusterLogRatio->Fill(energy, TMath::Log(ampMax/energy), GetEventWeight());
4453 
4454  // Get the ratio and log ratio to all cells in cluster
4455  for (Int_t ipos = 0; ipos < clus->GetNCells(); ipos++)
4456  {
4457  Int_t id = clus->GetCellsAbsId()[ipos];
4458 
4459  // Recalibrate cell energy if needed
4460  Float_t amp = cells->GetCellAmplitude(id);
4461  GetCaloUtils()->RecalibrateCellAmplitude(amp, GetCalorimeter(), id);
4462 
4463  fhECellClusterRatio ->Fill(energy, amp/energy , GetEventWeight());
4464  fhECellClusterLogRatio->Fill(energy, TMath::Log(amp/energy), GetEventWeight());
4465  }
4466 
4467  //
4468  // Recalculate shower shape for different W0 and cell cuts
4469  //
4470  if ( GetCalorimeter() == kEMCAL )
4471  {
4472  Float_t l0org = clus->GetM02();
4473  Float_t l1org = clus->GetM20();
4474  Float_t dorg = clus->GetDispersion();
4475  Float_t w0org = GetCaloUtils()->GetEMCALRecoUtils()->GetW0();
4476 
4477  Int_t tagMC = -1;
4478  if(IsDataMC() && clus->GetNLabels() > 0)
4479  {
4480  Int_t tag = GetMCAnalysisUtils()->CheckOrigin(clus->GetLabels(),clus->GetNLabels(), GetReader(),GetCalorimeter());
4481 
4482  if( GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCPhoton) &&
4483  !GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCPi0) &&
4484  !GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCEta) &&
4485  !GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCConversion) ){
4486  tagMC = 0;
4487  } // Pure Photon
4488  else if( GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCElectron) &&
4489  !GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCConversion) ){
4490  tagMC = 1;
4491  } // Electron
4492  else if( GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCConversion) ){
4493  tagMC = 2;
4494  } // Conversion
4495  else if( GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCPi0) ){
4496  tagMC = 3;
4497  }// Pi0
4498  else if(!GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCEta) &&
4499  !GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCPhoton) ){
4500  tagMC = 4;
4501  } // Hadron
4502  } // Is MC
4503 
4504  // To be done properly with setters and data members ...
4505  Float_t cellEmin [] = {0.05,0.1,0.15,0.2};
4506  Float_t cellTmin [] = {50.,100.,1000000.};
4507 
4508  for(Int_t iw = 0; iw < 12; iw++)
4509  {
4510  for(Int_t iEmin = 0; iEmin < 4; iEmin++)
4511  {
4512  for(Int_t iTmin = 0; iTmin < 3; iTmin++)
4513  {
4514  GetCaloUtils()->GetEMCALRecoUtils()->SetW0(4+iw*0.05);
4515 
4516  Float_t newEnergy = 0;
4517  Float_t l0 = 0, l1 = 0;
4518  Float_t disp = 0, dEta = 0, dPhi = 0;
4519  Float_t sEta = 0, sPhi = 0, sEtaPhi = 0;
4520 
4521  GetCaloUtils()->GetEMCALRecoUtils()->RecalculateClusterShowerShapeParametersWithCellCuts(GetEMCALGeometry(),
4522  GetReader()->GetEMCALCells(), clus,
4523  cellEmin[iEmin], cellTmin[iTmin], bc,
4524  newEnergy, l0, l1, disp, dEta, dPhi,
4525  sEta, sPhi, sEtaPhi);
4526 
4527 
4528 
4529  fhLambda0ForW0AndCellCuts[iw][iEmin][iTmin]->Fill(enOrg, l0, GetEventWeight());
4530 // fhLambda1ForW0AndCellCuts[iw][iEmin][iTmin]->Fill(enOrg, l0, GetEventWeight());
4531 
4532  if(TMath::Abs(fClusterMomentum.Eta()) < 0.15)
4533  fhLambda0ForW0AndCellCutsEta0[iw][iEmin][iTmin]->Fill(enOrg, l0, GetEventWeight());
4534 
4535  } // E cell loop
4536  } // Time loop
4537 
4538  if(IsDataMC() && tagMC >= 0)
4539  {
4540  fhLambda0ForW0MC[iw][tagMC]->Fill(energy, clus->GetM02(), GetEventWeight());
4541 // fhLambda1ForW0MC[iw][tagMC]->Fill(energy, clus->GetM20(), GetEventWeight());
4542  }
4543  } // w0 loop
4544 
4545  // Set the original values back
4546  clus->SetM02(l0org);
4547  clus->SetM20(l1org);
4548  clus->SetDispersion(dorg);
4549 
4550  GetCaloUtils()->GetEMCALRecoUtils()->SetW0(w0org);
4551 
4552  } // EMCAL
4553 
4554  clus->SetE(enOrg);
4555 }
4556 
4557 
4564 //void AliAnaCalorimeterQA::RecalculateClusterShowerShapeParametersWithCellCut(const AliEMCALGeometry * geom,
4565 // AliVCaloCells* cells, AliVCluster * cluster,
4566 // Float_t eCellMin, Float_t & newEnergy,
4567 // Float_t & l0, Float_t & l1,
4568 // Float_t & disp, Float_t & dEta, Float_t & dPhi,
4569 // Float_t & sEta, Float_t & sPhi, Float_t & sEtaPhi)
4570 //{
4571 // if(!cluster)
4572 // {
4573 // AliWarning("Cluster pointer null!");
4574 // return;
4575 // }
4576 //
4577 // Double_t eCell = 0.;
4578 // Float_t fraction = 1.;
4579 // Float_t recalFactor = 1.;
4580 //
4581 // Int_t iSupMod = -1;
4582 // Int_t iTower = -1;
4583 // Int_t iIphi = -1;
4584 // Int_t iIeta = -1;
4585 // Int_t iphi = -1;
4586 // Int_t ieta = -1;
4587 // Double_t etai = -1.;
4588 // Double_t phii = -1.;
4589 //
4590 // Int_t nstat = 0 ;
4591 // Float_t wtot = 0.;
4592 // Double_t w = 0.;
4593 // Double_t etaMean = 0.;
4594 // Double_t phiMean = 0.;
4595 //
4596 // newEnergy = 0;
4597 //
4598 // Bool_t shared = GetCaloUtils()-> IsClusterSharedByTwoSuperModules(geom,cluster);
4599 //
4600 // Float_t energy = GetCaloUtils()->RecalibrateClusterEnergy(cluster, cells);
4601 //
4608 //
4609 // // Loop on cells, get weighted parameters
4610 // for(Int_t iDigit=0; iDigit < cluster->GetNCells(); iDigit++)
4611 // {
4612 // // Get from the absid the supermodule, tower and eta/phi numbers
4613 // geom->GetCellIndex(cluster->GetCellAbsId(iDigit),iSupMod,iTower,iIphi,iIeta);
4614 // geom->GetCellPhiEtaIndexInSModule(iSupMod,iTower,iIphi,iIeta, iphi,ieta);
4615 //
4616 // // Get the cell energy, if recalibration is on, apply factors
4617 // fraction = cluster->GetCellAmplitudeFraction(iDigit);
4618 // if(fraction < 1e-4) fraction = 1.; // in case unfolding is off
4619 //
4620 // if(GetCaloUtils()->GetEMCALRecoUtils()->IsRecalibrationOn())
4621 // {
4622 // recalFactor = GetCaloUtils()->GetEMCALRecoUtils()->GetEMCALChannelRecalibrationFactor(iSupMod,ieta,iphi);
4623 // }
4624 //
4625 // eCell = cells->GetCellAmplitude(cluster->GetCellAbsId(iDigit))*fraction*recalFactor;
4626 //
4627 // // In case of a shared cluster, index of SM in C side, columns start at 48 and ends at 48*2
4628 // // C Side impair SM, nSupMod%2=1; A side pair SM, nSupMod%2=0
4629 // if(shared && iSupMod%2) ieta+=AliEMCALGeoParams::fgkEMCALCols;
4630 //
4631 // if(energy > 0 && eCell > eCellMin)
4632 // {
4633 // if(GetCaloUtils()->IsMCECellClusFracCorrectionOn())
4634 // eCell*=GetCaloUtils()->GetMCECellClusFracCorrection(eCell,energy);// /simuTotWeight;
4635 //
4636 // newEnergy+=eCell;
4637 //
4638 // w = GetCaloUtils()->GetEMCALRecoUtils()->GetCellWeight(eCell,energy);
4639 //
4640 // // correct weight, ONLY in simulation
4641 // //w *= (fWSimu[0] - fWSimu[1] * w );
4642 //
4643 // etai=(Double_t)ieta;
4644 // phii=(Double_t)iphi;
4645 //
4646 // if(w > 0.0)
4647 // {
4648 // wtot += w ;
4649 // nstat++;
4650 // //Shower shape
4651 // sEta += w * etai * etai ;
4652 // etaMean += w * etai ;
4653 // sPhi += w * phii * phii ;
4654 // phiMean += w * phii ;
4655 // sEtaPhi += w * etai * phii ;
4656 // }
4657 // }
4658 // else if(energy == 0 || (eCellMin <0.01 && eCell == 0)) AliError(Form("Wrong energy %f and/or amplitude %f", eCell, energy));
4659 // }//cell loop
4660 //
4661 // // Normalize to the weight
4662 // if (wtot > 0)
4663 // {
4664 // etaMean /= wtot ;
4665 // phiMean /= wtot ;
4666 // }
4667 // //else
4668 // // AliError(Form("Wrong weight %f", wtot));
4669 //
4670 // // Calculate dispersion
4671 // for(Int_t iDigit=0; iDigit < cluster->GetNCells(); iDigit++)
4672 // {
4673 // // Get from the absid the supermodule, tower and eta/phi numbers
4674 // geom->GetCellIndex(cluster->GetCellAbsId(iDigit),iSupMod,iTower,iIphi,iIeta);
4675 // geom->GetCellPhiEtaIndexInSModule(iSupMod,iTower,iIphi,iIeta, iphi,ieta);
4676 //
4677 // // Get the cell energy, if recalibration is on, apply factors
4678 // fraction = cluster->GetCellAmplitudeFraction(iDigit);
4679 // if(fraction < 1e-4) fraction = 1.; // in case unfolding is off
4680 // if (GetCaloUtils()->GetEMCALRecoUtils()->IsRecalibrationOn())
4681 // {
4682 // recalFactor = GetCaloUtils()->GetEMCALRecoUtils()->GetEMCALChannelRecalibrationFactor(iSupMod,ieta,iphi);
4683 // }
4684 //
4685 // eCell = cells->GetCellAmplitude(cluster->GetCellAbsId(iDigit))*fraction*recalFactor;
4686 //
4687 // // In case of a shared cluster, index of SM in C side, columns start at 48 and ends at 48*2
4688 // // C Side impair SM, nSupMod%2=1; A side pair SM, nSupMod%2=0
4689 // if(shared && iSupMod%2) ieta+=AliEMCALGeoParams::fgkEMCALCols;
4690 //
4691 // if(energy > 0 && eCell > eCellMin)
4692 // {
4693 // if(GetCaloUtils()->IsMCECellClusFracCorrectionOn())
4694 // eCell*=GetCaloUtils()->GetMCECellClusFracCorrection(eCell,energy); // /simuTotWeight;
4695 //
4696 // w = GetCaloUtils()->GetEMCALRecoUtils()->GetCellWeight(eCell,energy);
4697 //
4698 // //correct weight, ONLY in simulation
4699 // //w *= (fWSimu[0] - fWSimu[1] * w );
4700 //
4701 // etai=(Double_t)ieta;
4702 // phii=(Double_t)iphi;
4703 // if(w > 0.0)
4704 // {
4705 // disp += w *((etai-etaMean)*(etai-etaMean)+(phii-phiMean)*(phii-phiMean));
4706 // dEta += w * (etai-etaMean)*(etai-etaMean) ;
4707 // dPhi += w * (phii-phiMean)*(phii-phiMean) ;
4708 // }
4709 // }
4710 // else if(energy == 0 || (eCellMin <0.01 && eCell == 0)) AliError(Form("Wrong energy %f and/or amplitude %f", eCell, energy));
4711 // } // cell loop
4712 //
4713 // // Normalize to the weigth and set shower shape parameters
4714 // if (wtot > 0 && nstat > 1)
4715 // {
4716 // disp /= wtot ;
4717 // dEta /= wtot ;
4718 // dPhi /= wtot ;
4719 // sEta /= wtot ;
4720 // sPhi /= wtot ;
4721 // sEtaPhi /= wtot ;
4722 //
4723 // sEta -= etaMean * etaMean ;
4724 // sPhi -= phiMean * phiMean ;
4725 // sEtaPhi -= etaMean * phiMean ;
4726 //
4727 // l0 = (0.5 * (sEta + sPhi) + TMath::Sqrt( 0.25 * (sEta - sPhi) * (sEta - sPhi) + sEtaPhi * sEtaPhi ));
4728 // l1 = (0.5 * (sEta + sPhi) - TMath::Sqrt( 0.25 * (sEta - sPhi) * (sEta - sPhi) + sEtaPhi * sEtaPhi ));
4729 // }
4730 // else
4731 // {
4732 // l0 = 0. ;
4733 // l1 = 0. ;
4734 // dEta = 0. ; dPhi = 0. ; disp = 0. ;
4735 // sEta = 0. ; sPhi = 0. ; sEtaPhi = 0. ;
4736 // }
4737 //}
4738 
4739 
4740 
charge
Int_t charge
Definition: MuonTrackingEfficiency.C:54
AliHistogramRanges::GetHistoZMax
Float_t GetHistoZMax() const
Definition: AliHistogramRanges.h:261
AliAnaCalorimeterQA::fhCaloV0MCorrEClusters
TH2F * fhCaloV0MCorrEClusters
! Calo vs V0 multiplicity, total measured cluster energy
Definition: AliAnaCalorimeterQA.h:441
AliHistogramRanges::GetHistoPtMax
Float_t GetHistoPtMax() const
Definition: AliHistogramRanges.h:39
AliAnaCalorimeterQA::fhTimeIdLowGain
TH2F * fhTimeIdLowGain
! Time vs Absolute cell Id, low gain
Definition: AliAnaCalorimeterQA.h:412
AliAnaCalorimeterQA::fhDeltaCellClusterRE
TH2F * fhDeltaCellClusterRE
! R cluster - R cell distribution (cm) vs cluster energy
Definition: AliAnaCalorimeterQA.h:391
AliAnaCalorimeterQA::IsGoodCluster
Bool_t IsGoodCluster(Int_t absIdMax, Float_t m02, Int_t nCellsPerCluster, AliVCaloCells *cells)
Identify cluster as exotic or not.
Definition: AliAnaCalorimeterQA.cxx:4159
AliAnaCalorimeterQA::fhDCALPHOSCorrNCells
TH2F * fhDCALPHOSCorrNCells
! DCAL vs PHOS, number of cells
Definition: AliAnaCalorimeterQA.h:431
pdg
Int_t pdg
Definition: charmCutsOptimization.C:34
AliAnaCalorimeterQA::fhNCellsPerClusterNoCut
TH2F * fhNCellsPerClusterNoCut
! N cells per cluster vs cluster energy vs eta of cluster
Definition: AliAnaCalorimeterQA.h:310
AliAnaCalorimeterQA::fhLambda0ForW0AndCellCutsEta0
TH2F * fhLambda0ForW0AndCellCutsEta0[12][4][3]
! L0 for different w0 and cell cuts, |eta| < 0.15
Definition: AliAnaCalorimeterQA.h:500
AliAnaCalorimeterQA::fhMCNeutral1EOverPR02
TH2F * fhMCNeutral1EOverPR02
! p/E for track-cluster matches, dR < 0.2, MC neutral
Definition: AliAnaCalorimeterQA.h:581
AliAnaCalorimeterQA::fhNCellsMod
TH2F * fhNCellsMod
! Number of towers/crystals with signal for different module, Reco
Definition: AliAnaCalorimeterQA.h:474
AliAnaCalorimeterQA::fhXYZ
TH3F * fhXYZ
! cluster X vs Y vs Z (cm)
Definition: AliAnaCalorimeterQA.h:378
AliAnaCalorimeterQA::fh1EOverPMod
TH2F * fh1EOverPMod
! p/E for track-cluster matches, per SM
Definition: AliAnaCalorimeterQA.h:560
AliAnaCalorimeterQA::fh1EOverPR02Mod
TH2F * fh1EOverPR02Mod
! p/E for track-cluster matches, dR < 0.2, per SM
Definition: AliAnaCalorimeterQA.h:564
AliHistogramRanges::GetHistoNClusterCellMin
Int_t GetHistoNClusterCellMin() const
Definition: AliHistogramRanges.h:159
AliAnaCalorimeterQA::fhDeltaCellClusterYNCells
TH2F * fhDeltaCellClusterYNCells
! Y cluster - Y cell distribution (cm) vs N cells in cluster
Definition: AliAnaCalorimeterQA.h:388
AliAnaCalorimeterQA::fhPt
TH1F * fhPt
! pT distribution, Reco
Definition: AliAnaCalorimeterQA.h:275
AliHistogramRanges::GetHistoPtMin
Float_t GetHistoPtMin() const
Definition: AliHistogramRanges.h:38
AliHistogramRanges::GetHistoFinePtBins
Int_t GetHistoFinePtBins() const
Definition: AliHistogramRanges.h:51
AliAnaCalorimeterQA::fhClusterMaxCellECross
TH2F * fhClusterMaxCellECross
! 1 - Energy in cross around max energy cell / max energy cell vs cluster energy, good clusters ...
Definition: AliAnaCalorimeterQA.h:327
ymax
const Double_t ymax
Definition: AddTaskCFDStar.C:7
AliAnaCalorimeterQA::fFillPi0PairDiffTime
Bool_t fFillPi0PairDiffTime
Fill time difference histograms of cluster pairs in pi0 mass window, only if fFillAllPi0Histo=kTRUE.
Definition: AliAnaCalorimeterQA.h:217
Double_t
double Double_t
Definition: External.C:58
AliCaloTrackReader.h
AliHistogramRanges::GetHistoShowerShapeBins
Int_t GetHistoShowerShapeBins() const
Definition: AliHistogramRanges.h:268
TH3F
Definition: External.C:260
AliAnaCalorimeterQA::fhGridCellsE
TH2F * fhGridCellsE
! Cells ordered in column/row for different module, weighted with energy, Reco
Definition: AliAnaCalorimeterQA.h:484