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