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