]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGGA/CaloTrackCorrelations/AliAnaCalorimeterQA.cxx
Merge branch 'master' of https://git.cern.ch/reps/AliRoot
[u/mrichter/AliRoot.git] / PWGGA / CaloTrackCorrelations / AliAnaCalorimeterQA.cxx
CommitLineData
9725fd2a 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 **************************************************************************/
9725fd2a 15
a6f26052 16//_________________________________________________________________________
17// Class to check results from simulations or reconstructed real data.
18// Fill few histograms and do some checking plots
19//
20//-- Author: Gustavo Conesa (INFN-LNF)
21//_________________________________________________________________________
9725fd2a 22
23
a6f26052 24// --- ROOT system ---
d55bb5e1 25#include <TObjArray.h>
26#include <TParticle.h>
27#include <TDatabasePDG.h>
28#include <TH3F.h>
0c1383b5 29#include <TObjString.h>
9725fd2a 30
a6f26052 31//---- AliRoot system ----
9725fd2a 32#include "AliAnaCalorimeterQA.h"
33#include "AliCaloTrackReader.h"
34#include "AliStack.h"
c8fe2783 35#include "AliVCaloCells.h"
ff45398a 36#include "AliFiducialCut.h"
c8fe2783 37#include "AliVCluster.h"
d55bb5e1 38#include "AliVTrack.h"
c8fe2783 39#include "AliVEvent.h"
902aa95c 40#include "AliVEventHandler.h"
902aa95c 41#include "AliAODMCParticle.h"
42#include "AliMCAnalysisUtils.h"
9725fd2a 43
c5693f62 44// --- Detectors ---
45#include "AliPHOSGeoUtils.h"
46#include "AliEMCALGeometry.h"
47
9725fd2a 48ClassImp(AliAnaCalorimeterQA)
c8fe2783 49
649b825d 50//________________________________________
c8fe2783 51AliAnaCalorimeterQA::AliAnaCalorimeterQA() :
0cea6003 52AliAnaCaloTrackCorrBaseClass(),
649b825d 53
54//Switches
e6fec6f5 55fFillAllCellTimeHisto(kTRUE),
45769d5b 56fFillAllPosHisto(kFALSE), fFillAllPosHisto2(kTRUE),
57fFillAllTH3(kFALSE),
35c71d5c 58fFillAllTMHisto(kTRUE), fFillAllPi0Histo(kTRUE),
649b825d 59fCorrelate(kTRUE), fStudyBadClusters(kFALSE),
f1538a5f 60fStudyClustersAsymmetry(kFALSE), fStudyExotic(kFALSE),
61fStudyWeight(kFALSE),
649b825d 62
63//Parameters and cuts
35c71d5c 64fNModules(12), fNRCU(2),
65fNMaxCols(48), fNMaxRows(24),
f15c25da 66fTimeCutMin(-10000), fTimeCutMax(10000),
07e4c878 67fCellAmpMin(0), fEMCALCellAmpMin(0),
68fPHOSCellAmpMin(0), fMinInvMassECut(0),
649b825d 69
f1538a5f 70// Exotic
71fExoNECrossCuts(0), fExoECrossCuts(),
72fExoNDTimeCuts(0), fExoDTimeCuts(),
73
f4282b09 74fClusterMomentum(), fClusterMomentum2(),
75fPrimaryMomentum(),
649b825d 76//Histograms
9e9f04cb 77fhE(0), fhPt(0),
78fhPhi(0), fhEta(0), fhEtaPhiE(0),
79fhECharged(0), fhPtCharged(0),
80fhPhiCharged(0), fhEtaCharged(0), fhEtaPhiECharged(0),
521636d2 81
82//Invariant mass
9e9f04cb 83fhIM(0 ), fhAsym(0),
a82b4462 84
85fhNCellsPerCluster(0), fhNCellsPerClusterNoCut(0), fhNClusters(0),
3129a79e 86
e1e62b89 87//Timing
9e9f04cb 88fhClusterTimeEnergy(0), fhCellTimeSpreadRespectToCellMax(0),
9e9f04cb 89fhCellIdCellLargeTimeSpread(0), fhClusterPairDiffTimeE(0),
9e9f04cb 90fhClusterMaxCellCloseCellRatio(0), fhClusterMaxCellCloseCellDiff(0),
91fhClusterMaxCellDiff(0), fhClusterMaxCellDiffNoCut(0),
a82b4462 92fhClusterMaxCellDiffAverageTime(0), fhClusterMaxCellDiffWeightedTime(0),
1a72f6c5 93fhClusterMaxCellECross(0),
649b825d 94fhLambda0(0), fhLambda1(0), fhDispersion(0),
715fd81f 95
649b825d 96//bad clusters
97fhBadClusterEnergy(0), fhBadClusterTimeEnergy(0),
98fhBadClusterPairDiffTimeE(0), fhBadCellTimeSpreadRespectToCellMax(0),
9e9f04cb 99fhBadClusterMaxCellCloseCellRatio(0), fhBadClusterMaxCellCloseCellDiff(0), fhBadClusterMaxCellDiff(0),
a82b4462 100fhBadClusterMaxCellDiffAverageTime(0), fhBadClusterMaxCellDiffWeightedTime(0),
1a72f6c5 101fhBadClusterMaxCellECross(0),
1a83b960 102fhBadClusterDeltaIEtaDeltaIPhiE0(0), fhBadClusterDeltaIEtaDeltaIPhiE2(0),
103fhBadClusterDeltaIEtaDeltaIPhiE6(0), fhBadClusterDeltaIA(0),
9e9f04cb 104
521636d2 105//Position
9e9f04cb 106fhRNCells(0), fhXNCells(0),
107fhYNCells(0), fhZNCells(0),
108fhRE(0), fhXE(0),
109fhYE(0), fhZE(0),
521636d2 110fhXYZ(0),
9e9f04cb 111fhRCellE(0), fhXCellE(0),
112fhYCellE(0), fhZCellE(0),
521636d2 113fhXYZCell(0),
9e9f04cb 114fhDeltaCellClusterRNCells(0), fhDeltaCellClusterXNCells(0),
115fhDeltaCellClusterYNCells(0), fhDeltaCellClusterZNCells(0),
116fhDeltaCellClusterRE(0), fhDeltaCellClusterXE(0),
117fhDeltaCellClusterYE(0), fhDeltaCellClusterZE(0),
649b825d 118
521636d2 119// Cells
701cbf54 120fhNCells(0), fhNCellsCutAmpMin(0),
121fhAmplitude(0), fhAmpId(0), fhEtaPhiAmp(0),
1a72f6c5 122fhTime(0), fhTimeVz(0),
123fhTimeId(0), fhTimeAmp(0),
638916c4 124fhAmpIdLowGain(0), fhTimeIdLowGain(0), fhTimeAmpLowGain(0),
125
1a72f6c5 126fhCellECross(0),
9e9f04cb 127fhCaloCorrNClusters(0), fhCaloCorrEClusters(0),
128fhCaloCorrNCells(0), fhCaloCorrECells(0),
129fhCaloV0SCorrNClusters(0), fhCaloV0SCorrEClusters(0),
130fhCaloV0SCorrNCells(0), fhCaloV0SCorrECells(0),
131fhCaloV0MCorrNClusters(0), fhCaloV0MCorrEClusters(0),
132fhCaloV0MCorrNCells(0), fhCaloV0MCorrECells(0),
133fhCaloTrackMCorrNClusters(0), fhCaloTrackMCorrEClusters(0),
134fhCaloTrackMCorrNCells(0), fhCaloTrackMCorrECells(0),
653aed3c 135fhCaloCenNClusters(0), fhCaloCenEClusters(0),
136fhCaloCenNCells(0), fhCaloCenECells(0),
137fhCaloEvPNClusters(0), fhCaloEvPEClusters(0),
138fhCaloEvPNCells(0), fhCaloEvPECells(0),
521636d2 139//Super-Module dependent histgrams
649b825d 140fhEMod(0), fhAmpMod(0), fhTimeMod(0),
141fhNClustersMod(0), fhNCellsMod(0),
142fhNCellsPerClusterMod(0), fhNCellsPerClusterModNoCut(0),
143
638916c4 144fhGridCells(0), fhGridCellsE(0), fhGridCellsTime(0),
145fhGridCellsLowGain(0), fhGridCellsELowGain(0), fhGridCellsTimeLowGain(0),
146fhTimeAmpPerRCU(0), fhIMMod(0),
649b825d 147
148// Weight studies
149fhECellClusterRatio(0), fhECellClusterLogRatio(0),
150fhEMaxCellClusterRatio(0), fhEMaxCellClusterLogRatio(0),
701cbf54 151fhECellTotalRatio(0), fhECellTotalLogRatio(0),
152fhECellTotalRatioMod(0), fhECellTotalLogRatioMod(0),
715fd81f 153
765206a5 154fhExoL0ECross(0), fhExoL1ECross(0),
155
715fd81f 156// MC and reco
649b825d 157fhRecoMCE(), fhRecoMCPhi(), fhRecoMCEta(),
158fhRecoMCDeltaE(), fhRecoMCRatioE(),
159fhRecoMCDeltaPhi(), fhRecoMCDeltaEta(),
715fd81f 160
521636d2 161// MC only
95aee5e1 162fhGenMCE(), fhGenMCPt(), fhGenMCEtaPhi(),
163fhGenMCAccE(), fhGenMCAccPt(), fhGenMCAccEtaPhi(),
35c71d5c 164
165//matched MC
649b825d 166fhEMVxyz(0), fhEMR(0),
167fhHaVxyz(0), fhHaR(0),
d55bb5e1 168fh1EOverP(0), fh2dR(0),
649b825d 169fh2EledEdx(0), fh2MatchdEdx(0),
d55bb5e1 170fhMCEle1EOverP(0), fhMCEle1dR(0), fhMCEle2MatchdEdx(0),
171fhMCChHad1EOverP(0), fhMCChHad1dR(0), fhMCChHad2MatchdEdx(0),
a054a582 172fhMCNeutral1EOverP(0), fhMCNeutral1dR(0), fhMCNeutral2MatchdEdx(0), fh1EOverPR02(0),
173fhMCEle1EOverPR02(0), fhMCChHad1EOverPR02(0), fhMCNeutral1EOverPR02(0),
653aed3c 174fh1EleEOverP(0), fhMCEle1EleEOverP(0),
175fhMCChHad1EleEOverP(0), fhMCNeutral1EleEOverP(0),
176fhTrackMatchedDEta(0), fhTrackMatchedDPhi(0), fhTrackMatchedDEtaDPhi(0),
177fhTrackMatchedDEtaPos(0), fhTrackMatchedDPhiPos(0), fhTrackMatchedDEtaDPhiPos(0)
9725fd2a 178{
a6f26052 179 //Default Ctor
afabc52f 180
649b825d 181 //Weight studies
701cbf54 182 for(Int_t i =0; i < 12; i++){
649b825d 183 fhLambda0ForW0[i] = 0;
1a72f6c5 184 //fhLambda1ForW0[i] = 0;
649b825d 185
186 for(Int_t j = 0; j < 5; j++){
187 fhLambda0ForW0MC[i][j] = 0;
1a72f6c5 188 //fhLambda1ForW0MC[i][j] = 0;
649b825d 189 }
190
191 }
c8fe2783 192
649b825d 193 //Cluster size
194 fhDeltaIEtaDeltaIPhiE0[0] = 0 ; fhDeltaIEtaDeltaIPhiE2[0] = 0; fhDeltaIEtaDeltaIPhiE6[0] = 0;
195 fhDeltaIEtaDeltaIPhiE0[1] = 0 ; fhDeltaIEtaDeltaIPhiE2[1] = 0; fhDeltaIEtaDeltaIPhiE6[1] = 0;
196 fhDeltaIA[0] = 0 ; fhDeltaIAL0[0] = 0; fhDeltaIAL1[0] = 0;
197 fhDeltaIA[1] = 0 ; fhDeltaIAL0[1] = 0; fhDeltaIAL1[1] = 0;
198 fhDeltaIANCells[0] = 0 ; fhDeltaIANCells[1] = 0;
199 fhDeltaIAMC[0] = 0 ; fhDeltaIAMC[1] = 0;
200 fhDeltaIAMC[2] = 0 ; fhDeltaIAMC[3] = 0;
2302a644 201
f1538a5f 202 // Exotic
203 for (Int_t ie = 0; ie < 10 ; ie++)
204 {
205 fhExoDTime[ie] = 0;
206 for (Int_t idt = 0; idt < 5 ; idt++)
207 {
208 fhExoNCell [ie][idt] = 0;
209 fhExoL0 [ie][idt] = 0;
765206a5 210 fhExoL1 [ie][idt] = 0;
f1538a5f 211 fhExoECross [ie][idt] = 0;
212 fhExoTime [ie][idt] = 0;
765206a5 213 fhExoL0NCell [ie][idt] = 0;
214 fhExoL1NCell [ie][idt] = 0;
f1538a5f 215 }
216 }
217
649b825d 218 // MC
c8fe2783 219
6aba6683 220 for(Int_t i = 0; i < 7; i++)
45769d5b 221 {
649b825d 222 fhRecoMCE[i][0] = 0; fhRecoMCE[i][1] = 0;
223 fhRecoMCPhi[i][0] = 0; fhRecoMCPhi[i][1] = 0;
224 fhRecoMCEta[i][0] = 0; fhRecoMCEta[i][1] = 0;
225 fhRecoMCDeltaE[i][0] = 0; fhRecoMCDeltaE[i][1] = 0;
226 fhRecoMCRatioE[i][0] = 0; fhRecoMCRatioE[i][1] = 0;
227 fhRecoMCDeltaPhi[i][0] = 0; fhRecoMCDeltaPhi[i][1] = 0;
45769d5b 228 fhRecoMCDeltaEta[i][0] = 0; fhRecoMCDeltaEta[i][1] = 0;
649b825d 229 }
2302a644 230
95aee5e1 231 for(Int_t i = 0; i < 4; i++)
232 {
233 fhGenMCE[i] = 0;
234 fhGenMCPt[i] = 0;
235 fhGenMCEtaPhi[i] = 0;
236 fhGenMCAccE[i] = 0;
237 fhGenMCAccPt[i] = 0;
238 fhGenMCAccEtaPhi[i] = 0;
239 }
240
649b825d 241 //Initialize parameters
242 InitParameters();
243}
3f5990d6 244
b94e038e 245//______________________________________________________________________________________________________________________
c5693f62 246void AliAnaCalorimeterQA::BadClusterHistograms(AliVCluster* clus, const TObjArray *caloClusters, AliVCaloCells * cells,
b94e038e 247 Int_t absIdMax, Double_t maxCellFraction, Float_t eCrossFrac,
45769d5b 248 Double_t tmax)
649b825d 249{
250 //Bad cluster histograms
b0114dba 251
252 // printf("AliAnaCalorimeterQA::BadClusterHistograms() - Event %d - Calorimeter %s \n \t E %f, n cells %d, max cell absId %d, maxCellFrac %f\n",
f3da05f3 253 // GetReader()->GetEventNumber(), GetCalorimeterString().Data(),
b0114dba 254 // clus->E(),clus->GetNCells(),absIdMax,maxCellFraction);
1a72f6c5 255
649b825d 256 fhBadClusterEnergy ->Fill(clus->E());
257 Double_t tof = clus->GetTOF()*1.e9;
1a83b960 258 fhBadClusterTimeEnergy ->Fill(clus->E(),tof);
259 fhBadClusterMaxCellDiff ->Fill(clus->E(),maxCellFraction);
260 fhBadClusterMaxCellECross->Fill(clus->E(),eCrossFrac);
261
262 if(fStudyClustersAsymmetry) ClusterAsymmetryHistograms(clus,absIdMax,kFALSE);
1a72f6c5 263
a82b4462 264 //Clusters in event time differencem bad minus good
2302a644 265
45769d5b 266 for(Int_t iclus2 = 0; iclus2 < caloClusters->GetEntriesFast(); iclus2++ )
267 {
649b825d 268 AliVCluster* clus2 = (AliVCluster*)caloClusters->At(iclus2);
a6f26052 269
649b825d 270 if(clus->GetID()==clus2->GetID()) continue;
715fd81f 271
a82b4462 272 Float_t maxCellFraction2 = 0.;
273 Int_t absIdMax2 = GetCaloUtils()->GetMaxEnergyCell(cells, clus2,maxCellFraction2);
45769d5b 274 if(IsGoodCluster(absIdMax2,cells))
275 {
a82b4462 276 Double_t tof2 = clus2->GetTOF()*1.e9;
649b825d 277 fhBadClusterPairDiffTimeE ->Fill(clus->E(), (tof-tof2));
649b825d 278 }
a82b4462 279
649b825d 280 } // loop
924e319f 281
649b825d 282 // Max cell compared to other cells in cluster
e6fec6f5 283 if(fFillAllCellTimeHisto)
2747966a 284 {
45769d5b 285 // Get some time averages
286 Double_t timeAverages[2] = {0.,0.};
287 CalculateAverageTime(clus, cells, timeAverages);
288
649b825d 289 fhBadClusterMaxCellDiffAverageTime ->Fill(clus->E(),tmax-timeAverages[0]);
649b825d 290 fhBadClusterMaxCellDiffWeightedTime ->Fill(clus->E(),tmax-timeAverages[1]);
649b825d 291 }
715fd81f 292
2747966a 293 for (Int_t ipos = 0; ipos < clus->GetNCells(); ipos++)
294 {
649b825d 295 Int_t absId = clus->GetCellsAbsId()[ipos];
2747966a 296 if(absId!=absIdMax && cells->GetCellAmplitude(absIdMax) > 0.01)
297 {
649b825d 298 Float_t frac = cells->GetCellAmplitude(absId)/cells->GetCellAmplitude(absIdMax);
299
300 fhBadClusterMaxCellCloseCellRatio->Fill(clus->E(),frac);
301 fhBadClusterMaxCellCloseCellDiff ->Fill(clus->E(),cells->GetCellAmplitude(absIdMax)-cells->GetCellAmplitude(absId));
302
e6fec6f5 303 if(fFillAllCellTimeHisto)
dbba06ca 304 {
649b825d 305 Double_t time = cells->GetCellTime(absId);
0cea6003 306 GetCaloUtils()->RecalibrateCellTime(time, GetCalorimeter(), absId,GetReader()->GetInputEvent()->GetBunchCrossNumber());
dbba06ca 307
649b825d 308 Float_t diff = (tmax-time*1e9);
309 fhBadCellTimeSpreadRespectToCellMax->Fill(clus->E(), diff);
310
e6fec6f5 311 }
649b825d 312 }// Not max
313 }//loop
715fd81f 314
649b825d 315}
715fd81f 316
dbf54f1e 317//______________________________________________________________________
318void AliAnaCalorimeterQA::CalculateAverageTime(AliVCluster *clus,
319 AliVCaloCells* cells,
320 Double_t timeAverages[2])
649b825d 321{
322 // Calculate time averages and weights
a95eac90 323
649b825d 324 // First recalculate energy in case non linearity was applied
325 Float_t energy = 0;
326 Float_t ampMax = 0, amp = 0;
ecdde216 327// Int_t absIdMax =-1;
2747966a 328 for (Int_t ipos = 0; ipos < clus->GetNCells(); ipos++)
329 {
649b825d 330 Int_t id = clus->GetCellsAbsId()[ipos];
e1e62b89 331
649b825d 332 //Recalibrate cell energy if needed
333 amp = cells->GetCellAmplitude(id);
0cea6003 334 GetCaloUtils()->RecalibrateCellAmplitude(amp,GetCalorimeter(), id);
4c8f7c2e 335
649b825d 336 energy += amp;
e1e62b89 337
2747966a 338 if(amp> ampMax)
339 {
649b825d 340 ampMax = amp;
ecdde216 341// absIdMax = id;
649b825d 342 }
649b825d 343 } // energy loop
344
345 // Calculate average time of cells in cluster and weighted average
dbf54f1e 346 Double_t aTime = 0;
347 Double_t wTime = 0;
348 Float_t wTot = 0;
349 Double_t time = 0;
350 Int_t id =-1;
351 Double_t w = 0;
352 Int_t ncells = clus->GetNCells();
45769d5b 353
2747966a 354 for (Int_t ipos = 0; ipos < ncells; ipos++)
355 {
dbf54f1e 356 id = clus ->GetCellsAbsId()[ipos];
357 amp = cells->GetCellAmplitude(id);
358 time = cells->GetCellTime(id);
649b825d 359
360 //Recalibrate energy and time
0cea6003 361 GetCaloUtils()->RecalibrateCellAmplitude(amp , GetCalorimeter(), id);
362 GetCaloUtils()->RecalibrateCellTime (time, GetCalorimeter(), id, GetReader()->GetInputEvent()->GetBunchCrossNumber());
dbba06ca 363
dbf54f1e 364 w = GetCaloUtils()->GetEMCALRecoUtils()->GetCellWeight(cells->GetCellAmplitude(id),energy);
649b825d 365 aTime += time*1e9;
366 wTime += time*1e9 * w;
367 wTot += w;
35c71d5c 368
649b825d 369 }
f3138ecf 370
dbf54f1e 371 if(ncells > 0) aTime /= ncells;
372 else aTime = 0;
649b825d 373
dbf54f1e 374 if(wTot > 0) wTime /= wTot;
f3138ecf 375 else wTime = 0;
376
dbf54f1e 377 timeAverages[0] = aTime;
378 timeAverages[1] = wTime;
39de6caa 379
649b825d 380}
381
382//____________________________________________________________
383void AliAnaCalorimeterQA::CellHistograms(AliVCaloCells *cells)
384{
385 // Plot histograms related to cells only
9e9f04cb 386
649b825d 387 Int_t ncells = cells->GetNumberOfCells();
701cbf54 388 if( ncells > 0 ) fhNCells->Fill(ncells) ;
389
390 Int_t ncellsCut = 0;
391 Float_t ecellsCut = 0;
39de6caa 392
2db10729 393 AliDebug(1,Form("%s cell entries %d", GetCalorimeterString().Data(), ncells));
649b825d 394
395 //Init arrays and used variables
701cbf54 396 Int_t *nCellsInModule = new Int_t [fNModules];
397 Float_t *eCellsInModule = new Float_t[fNModules];
398
399 for(Int_t imod = 0; imod < fNModules; imod++ )
400 {
401 nCellsInModule[imod] = 0 ;
402 eCellsInModule[imod] = 0.;
403 }
a82b4462 404
649b825d 405 Int_t icol = -1;
406 Int_t irow = -1;
407 Int_t iRCU = -1;
408 Float_t amp = 0.;
409 Double_t time = 0.;
410 Int_t id = -1;
638916c4 411 Bool_t highG = kFALSE;
649b825d 412 Float_t recalF = 1.;
a82b4462 413 Int_t bc = GetReader()->GetInputEvent()->GetBunchCrossNumber();
649b825d 414
45769d5b 415 for (Int_t iCell = 0; iCell < cells->GetNumberOfCells(); iCell++)
416 {
2db10729 417 AliDebug(2,Form("Cell : amp %f, absId %d", cells->GetAmplitude(iCell), cells->GetCellNumber(iCell)));
45769d5b 418
0cea6003 419 Int_t nModule = GetModuleNumberCellIndexes(cells->GetCellNumber(iCell),GetCalorimeter(), icol, irow, iRCU);
2db10729 420
421 AliDebug(2,Form("\t module %d, column %d, row %d", nModule,icol,irow));
649b825d 422
f1538a5f 423 if(nModule < fNModules)
424 {
649b825d 425 //Check if the cell is a bad channel
45769d5b 426 if(GetCaloUtils()->IsBadChannelsRemovalSwitchedOn())
427 {
f3da05f3 428 if(GetCalorimeter()==kEMCAL)
ae2c2bc4 429 {
649b825d 430 if(GetCaloUtils()->GetEMCALChannelStatus(nModule,icol,irow)) continue;
431 }
ae2c2bc4 432 else
433 {
434 if(GetCaloUtils()->GetPHOSChannelStatus(nModule,icol,irow) ) continue;
649b825d 435 }
436 } // use bad channel map
437
438 amp = cells->GetAmplitude(iCell)*recalF;
439 time = cells->GetTime(iCell);
440 id = cells->GetCellNumber(iCell);
638916c4 441 highG = cells->GetCellHighGain(id);
649b825d 442
443 // Amplitude recalibration if set
0cea6003 444 GetCaloUtils()->RecalibrateCellAmplitude(amp, GetCalorimeter(), id);
649b825d 445
446 // Time recalibration if set
0cea6003 447 GetCaloUtils()->RecalibrateCellTime (time, GetCalorimeter(), id, GetReader()->GetInputEvent()->GetBunchCrossNumber());
649b825d 448
449 //Transform time to ns
450 time *= 1.0e9;
f15c25da 451
e6fec6f5 452 if(time < fTimeCutMin || time > fTimeCutMax)
a87e069d 453 {
2db10729 454 AliDebug(1,Form("Remove cell with Time %f",time));
649b825d 455 continue;
e6fec6f5 456 }
06f1b12a 457
458 // Remove exotic cells, defined only for EMCAL
f3da05f3 459 if(GetCalorimeter()==kEMCAL &&
06f1b12a 460 GetCaloUtils()->GetEMCALRecoUtils()->IsExoticCell(id, cells, bc)) continue;
461
649b825d 462 fhAmplitude->Fill(amp);
463 fhAmpId ->Fill(amp,id);
464 fhAmpMod ->Fill(amp,nModule);
638916c4 465 if(!highG) fhAmpIdLowGain->Fill(amp,id);
701cbf54 466 //E cross for exotic cells
467 if(amp > 0.05)
2747966a 468 {
701cbf54 469 fhCellECross->Fill(amp,1-GetECross(id,cells)/amp);
470 ecellsCut+=amp ;
471 eCellsInModule[nModule]+=amp ;
472 }
473
474 if ( amp > fCellAmpMin )
475 {
476 ncellsCut++ ;
477 nCellsInModule[nModule]++ ;
06f1b12a 478
649b825d 479 Int_t icols = icol;
480 Int_t irows = irow;
57d8227a 481
f3da05f3 482 if(GetCalorimeter()==kEMCAL)
57d8227a 483 {
649b825d 484 icols = (nModule % 2) ? icol + fNMaxCols : icol;
57d8227a 485 if(nModule < 10 )
486 irows = irow + fNMaxRows * Int_t(nModule / 2);
487 else // 1/3 SM
488 irows = irow + (fNMaxRows / 3) * Int_t(nModule / 2);
649b825d 489 }
57d8227a 490 else
491 {
06f1b12a 492 irows = irow + fNMaxRows * nModule;
649b825d 493 }
06f1b12a 494
649b825d 495 fhGridCells ->Fill(icols,irows);
496 fhGridCellsE->Fill(icols,irows,amp);
497
638916c4 498 if(!highG)
499 {
500 fhGridCellsLowGain ->Fill(icols,irows);
501 fhGridCellsELowGain->Fill(icols,irows,amp);
502 }
503
e6fec6f5 504 if(fFillAllCellTimeHisto)
a87e069d 505 {
f3da05f3 506 //printf("%s: time %g\n",GetCalorimeterString().Data(), time);
1a72f6c5 507
508 Double_t v[3] = {0,0,0}; //vertex ;
509 GetReader()->GetVertex(v);
510 if(amp > 0.5) fhTimeVz ->Fill(TMath::Abs(v[2]),time);
511
649b825d 512 fhTime ->Fill(time);
513 fhTimeId ->Fill(time,id);
514 fhTimeAmp ->Fill(amp,time);
515 fhGridCellsTime->Fill(icols,irows,time);
638916c4 516 if(!highG) fhGridCellsTimeLowGain->Fill(icols,irows,time);
649b825d 517 fhTimeMod ->Fill(time,nModule);
518 fhTimeAmpPerRCU [nModule*fNRCU+iRCU]->Fill(amp, time);
519
638916c4 520 if(!highG)
521 {
522 fhTimeIdLowGain ->Fill(time,id);
523 fhTimeAmpLowGain->Fill(amp,time);
524 }
525
649b825d 526 }
527 }
528
529 //Get Eta-Phi position of Cell
530 if(fFillAllPosHisto)
531 {
f3da05f3 532 if(GetCalorimeter()==kEMCAL && GetCaloUtils()->IsEMCALGeoMatrixSet()){
649b825d 533 Float_t celleta = 0.;
534 Float_t cellphi = 0.;
535 GetEMCALGeometry()->EtaPhiFromIndex(id, celleta, cellphi);
536
537 fhEtaPhiAmp->Fill(celleta,cellphi,amp);
538 Double_t cellpos[] = {0, 0, 0};
539 GetEMCALGeometry()->GetGlobal(id, cellpos);
540 fhXCellE->Fill(cellpos[0],amp) ;
541 fhYCellE->Fill(cellpos[1],amp) ;
542 fhZCellE->Fill(cellpos[2],amp) ;
543 Float_t rcell = TMath::Sqrt(cellpos[0]*cellpos[0]+cellpos[1]*cellpos[1]);//+cellpos[2]*cellpos[2]);
544 fhRCellE->Fill(rcell,amp) ;
545 fhXYZCell->Fill(cellpos[0],cellpos[1],cellpos[2]) ;
546 }//EMCAL Cells
f3da05f3 547 else if(GetCalorimeter()==kPHOS && GetCaloUtils()->IsPHOSGeoMatrixSet()){
649b825d 548 TVector3 xyz;
549 Int_t relId[4], module;
550 Float_t xCell, zCell;
551
552 GetPHOSGeometry()->AbsToRelNumbering(id,relId);
553 module = relId[0];
554 GetPHOSGeometry()->RelPosInModule(relId,xCell,zCell);
555 GetPHOSGeometry()->Local2Global(module,xCell,zCell,xyz);
556 Float_t rcell = TMath::Sqrt(xyz.X()*xyz.X()+xyz.Y()*xyz.Y());
557 fhXCellE ->Fill(xyz.X(),amp) ;
558 fhYCellE ->Fill(xyz.Y(),amp) ;
559 fhZCellE ->Fill(xyz.Z(),amp) ;
560 fhRCellE ->Fill(rcell ,amp) ;
561 fhXYZCell->Fill(xyz.X(),xyz.Y(),xyz.Z()) ;
562 }//PHOS cells
563 }//fill cell position histograms
564
649b825d 565 }//nmodules
566 }//cell loop
567
701cbf54 568 if( ncellsCut > 0 ) fhNCellsCutAmpMin->Fill(ncellsCut) ; //fill the cells after the cut on min amplitude and bad/exotic channels
649b825d 569
570 //Number of cells per module
45769d5b 571 for(Int_t imod = 0; imod < fNModules; imod++ )
572 {
2db10729 573 AliDebug(1,Form("Module %d calo %s cells %d", imod, GetCalorimeterString().Data(), nCellsInModule[imod]));
649b825d 574
575 fhNCellsMod->Fill(nCellsInModule[imod],imod) ;
649b825d 576 }
577
701cbf54 578 // Check energy distribution in calorimeter for selected cells
579 if(fStudyWeight)
580 {
581 for (Int_t iCell = 0; iCell < cells->GetNumberOfCells(); iCell++)
582 {
2db10729 583 AliDebug(2,Form("Cell : amp %f, absId %d", cells->GetAmplitude(iCell), cells->GetCellNumber(iCell)));
701cbf54 584
0cea6003 585 Int_t nModule = GetModuleNumberCellIndexes(cells->GetCellNumber(iCell),GetCalorimeter(), icol, irow, iRCU);
2db10729 586
587 AliDebug(2,Form("\t module %d, column %d, row %d", nModule,icol,irow));
701cbf54 588
589 if(nModule < fNModules)
590 {
591 //Check if the cell is a bad channel
592 if(GetCaloUtils()->IsBadChannelsRemovalSwitchedOn())
593 {
f3da05f3 594 if(GetCalorimeter()==kEMCAL)
701cbf54 595 {
596 if(GetCaloUtils()->GetEMCALChannelStatus(nModule,icol,irow)) continue;
597 }
598 else
599 {
600 if(GetCaloUtils()->GetPHOSChannelStatus(nModule,icol,irow) ) continue;
601 }
602 } // use bad channel map
603
604 amp = cells->GetAmplitude(iCell)*recalF;
605 time = cells->GetTime(iCell);
606 id = cells->GetCellNumber(iCell);
607
608 // Amplitude recalibration if set
0cea6003 609 GetCaloUtils()->RecalibrateCellAmplitude(amp, GetCalorimeter(), id);
701cbf54 610
611 // Time recalibration if set
0cea6003 612 GetCaloUtils()->RecalibrateCellTime (time, GetCalorimeter(), id, GetReader()->GetInputEvent()->GetBunchCrossNumber());
701cbf54 613
614 //Transform time to ns
615 time *= 1.0e9;
616
617 if(time < fTimeCutMin || time > fTimeCutMax)
618 {
2db10729 619 AliDebug(1,Form("Remove cell with Time %f",time));
701cbf54 620 continue;
621 }
622
623 // Remove exotic cells, defined only for EMCAL
f3da05f3 624 if(GetCalorimeter()==kEMCAL &&
701cbf54 625 GetCaloUtils()->GetEMCALRecoUtils()->IsExoticCell(id, cells, bc)) continue;
626
627 //E cross for exotic cells
628 if(amp > 0.05)
629 {
630 if(ecellsCut > 0)
631 {
632 Float_t ratio = amp/ecellsCut;
633 fhECellTotalRatio ->Fill(ecellsCut, ratio );
634 fhECellTotalLogRatio ->Fill(ecellsCut,TMath::Log(ratio));
635 }
636
637 if(eCellsInModule[nModule] > 0)
638 {
639 Float_t ratioMod = amp/eCellsInModule[nModule];
640 fhECellTotalRatioMod [nModule]->Fill(eCellsInModule[nModule], ratioMod );
641 fhECellTotalLogRatioMod[nModule]->Fill(eCellsInModule[nModule],TMath::Log(ratioMod));
642 }
643 }// amp > 0.5
644 }// nMod > 0 < Max
645 } // cell loop
646 } // weight studies
647
649b825d 648 delete [] nCellsInModule;
701cbf54 649 delete [] eCellsInModule;
649b825d 650
651}
652
653//__________________________________________________________________________
654void AliAnaCalorimeterQA::CellInClusterPositionHistograms(AliVCluster* clus)
655{
656 // Fill histograms releated to cell position
657
649b825d 658 Int_t nCaloCellsPerCluster = clus->GetNCells();
659 UShort_t * indexList = clus->GetCellsAbsId();
660 Float_t pos[3];
661 clus->GetPosition(pos);
662 Float_t clEnergy = clus->E();
663
664 //Loop on cluster cells
45769d5b 665 for (Int_t ipos = 0; ipos < nCaloCellsPerCluster; ipos++)
666 {
649b825d 667 // printf("Index %d\n",ipos);
668 Int_t absId = indexList[ipos];
669
670 //Get position of cell compare to cluster
671
f3da05f3 672 if(GetCalorimeter()==kEMCAL && GetCaloUtils()->IsEMCALGeoMatrixSet()){
649b825d 673
674 Double_t cellpos[] = {0, 0, 0};
675 GetEMCALGeometry()->GetGlobal(absId, cellpos);
676
677 fhDeltaCellClusterXNCells->Fill(pos[0]-cellpos[0],nCaloCellsPerCluster) ;
678 fhDeltaCellClusterYNCells->Fill(pos[1]-cellpos[1],nCaloCellsPerCluster) ;
679 fhDeltaCellClusterZNCells->Fill(pos[2]-cellpos[2],nCaloCellsPerCluster) ;
680
681 fhDeltaCellClusterXE->Fill(pos[0]-cellpos[0],clEnergy) ;
682 fhDeltaCellClusterYE->Fill(pos[1]-cellpos[1],clEnergy) ;
683 fhDeltaCellClusterZE->Fill(pos[2]-cellpos[2],clEnergy) ;
684
685 Float_t r = TMath::Sqrt(pos[0] *pos[0] + pos[1] * pos[1] );
686 Float_t rcell = TMath::Sqrt(cellpos[0]*cellpos[0] + cellpos[1]* cellpos[1]);
687
688 fhDeltaCellClusterRNCells->Fill(r-rcell, nCaloCellsPerCluster) ;
689 fhDeltaCellClusterRE ->Fill(r-rcell, clEnergy) ;
690
691 }//EMCAL and its matrices are available
f3da05f3 692 else if(GetCalorimeter()==kPHOS && GetCaloUtils()->IsPHOSGeoMatrixSet())
45769d5b 693 {
649b825d 694 TVector3 xyz;
695 Int_t relId[4], module;
696 Float_t xCell, zCell;
697
698 GetPHOSGeometry()->AbsToRelNumbering(absId,relId);
699 module = relId[0];
700 GetPHOSGeometry()->RelPosInModule(relId,xCell,zCell);
701 GetPHOSGeometry()->Local2Global(module,xCell,zCell,xyz);
702
703 fhDeltaCellClusterXNCells->Fill(pos[0]-xyz.X(),nCaloCellsPerCluster) ;
704 fhDeltaCellClusterYNCells->Fill(pos[1]-xyz.Y(),nCaloCellsPerCluster) ;
705 fhDeltaCellClusterZNCells->Fill(pos[2]-xyz.Z(),nCaloCellsPerCluster) ;
706
707 fhDeltaCellClusterXE->Fill(pos[0]-xyz.X(),clEnergy) ;
708 fhDeltaCellClusterYE->Fill(pos[1]-xyz.Y(),clEnergy) ;
709 fhDeltaCellClusterZE->Fill(pos[2]-xyz.Z(),clEnergy) ;
710
711 Float_t r = TMath::Sqrt(pos[0] * pos[0] + pos[1] * pos[1] );
712 Float_t rcell = TMath::Sqrt(xyz.X() * xyz.X() + xyz.Y() * xyz.Y());
713
714 fhDeltaCellClusterRNCells->Fill(r-rcell, nCaloCellsPerCluster) ;
715 fhDeltaCellClusterRE ->Fill(r-rcell, clEnergy) ;
716
717 }//PHOS and its matrices are available
718 }// cluster cell loop
719}
720
721//___________________________________________________________________________________________
b94e038e 722void AliAnaCalorimeterQA::ClusterAsymmetryHistograms(AliVCluster* clus, Int_t absIdMax,
723 Bool_t goodCluster)
649b825d 724{
725 // Study the shape of the cluster in cell units terms
726
727 //No use to study clusters with less than 4 cells
45769d5b 728 if( clus->GetNCells() <= 3 ) return;
649b825d 729
730 Int_t dIeta = 0;
731 Int_t dIphi = 0;
732
733 Int_t ietaMax=-1; Int_t iphiMax = 0; Int_t rcuMax = 0;
0cea6003 734 Int_t smMax = GetModuleNumberCellIndexes(absIdMax,GetCalorimeter(), ietaMax, iphiMax, rcuMax);
649b825d 735
6aba6683 736 for (Int_t ipos = 0; ipos < clus->GetNCells(); ipos++)
737 {
649b825d 738 Int_t absId = clus->GetCellsAbsId()[ipos];
739
740 Int_t ieta=-1; Int_t iphi = 0; Int_t rcu = 0;
0cea6003 741 Int_t sm = GetModuleNumberCellIndexes(absId,GetCalorimeter(), ieta, iphi, rcu);
649b825d 742
743 if(dIphi < TMath::Abs(iphi-iphiMax)) dIphi = TMath::Abs(iphi-iphiMax);
744
45769d5b 745 if(smMax==sm)
746 {
649b825d 747 if(dIeta < TMath::Abs(ieta-ietaMax)) dIeta = TMath::Abs(ieta-ietaMax);
748 }
45769d5b 749 else
750 {
649b825d 751 Int_t ietaShift = ieta;
752 Int_t ietaMaxShift = ietaMax;
753 if (ieta > ietaMax) ietaMaxShift+=48;
754 else ietaShift +=48;
755 if(dIeta < TMath::Abs(ietaShift-ietaMaxShift)) dIeta = TMath::Abs(ietaShift-ietaMaxShift);
756 }
757
649b825d 758 }// fill cell-cluster histogram loop
759
649b825d 760
1a83b960 761 Float_t dIA = 1.*(dIphi-dIeta)/(dIeta+dIphi);
762
763 if(goodCluster)
764 {
1a83b960 765 // Was cluster matched?
766 Bool_t matched = GetCaloPID()->IsTrackMatched(clus,GetCaloUtils(),GetReader()->GetInputEvent());
649b825d 767
1a83b960 768 if (clus->E() < 2 ) fhDeltaIEtaDeltaIPhiE0[matched]->Fill(dIeta,dIphi);
769 else if(clus->E() < 6 ) fhDeltaIEtaDeltaIPhiE2[matched]->Fill(dIeta,dIphi);
770 else fhDeltaIEtaDeltaIPhiE6[matched]->Fill(dIeta,dIphi);
771
772 fhDeltaIA[matched]->Fill(clus->E(),dIA);
773
6aba6683 774 if(clus->E() > 0.5)
775 {
1a83b960 776 fhDeltaIAL0[matched] ->Fill(clus->GetM02(),dIA);
777 fhDeltaIAL1[matched] ->Fill(clus->GetM20(),dIA);
778 fhDeltaIANCells[matched]->Fill(clus->GetNCells(),dIA);
649b825d 779 }
780
1a83b960 781 // Origin of clusters
782 Int_t nLabel = clus->GetNLabels();
783 Int_t* labels = clus->GetLabels();
6aba6683 784
45769d5b 785 if(IsDataMC())
786 {
0cea6003 787 Int_t tag = GetMCAnalysisUtils()->CheckOrigin(labels,nLabel, GetReader(),GetCalorimeter());
1a83b960 788 if( GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCPhoton) &&
789 !GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCPi0) &&
790 !GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCEta) &&
791 !GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCConversion) ){
792 fhDeltaIAMC[0]->Fill(clus->E(),dIA);//Pure Photon
793 }
794 else if ( GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCElectron) &&
795 !GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCConversion) ){
796 fhDeltaIAMC[1]->Fill(clus->E(),dIA);//Pure electron
797 }
798 else if ( GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCPhoton) &&
799 GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCConversion) ){
800 fhDeltaIAMC[2]->Fill(clus->E(),dIA);//Converted cluster
801 }
802 else if(!GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCPhoton)){
803 fhDeltaIAMC[3]->Fill(clus->E(),dIA);//Hadrons
804 }
805
806 } // MC
807 } // good cluster
808 else
809 {
810 if (clus->E() < 2 ) fhBadClusterDeltaIEtaDeltaIPhiE0->Fill(dIeta,dIphi);
811 else if(clus->E() < 6 ) fhBadClusterDeltaIEtaDeltaIPhiE2->Fill(dIeta,dIphi);
812 else fhBadClusterDeltaIEtaDeltaIPhiE6->Fill(dIeta,dIphi);
813
814 fhBadClusterDeltaIA->Fill(clus->E(),dIA);
1a83b960 815 }
649b825d 816}
817
b94e038e 818//__________________________________________________________________________________________________________________
819void AliAnaCalorimeterQA::ClusterHistograms(AliVCluster* clus, const TObjArray *caloClusters, AliVCaloCells * cells,
820 Int_t absIdMax, Double_t maxCellFraction, Float_t eCrossFrac,
45769d5b 821 Double_t tmax)
649b825d 822{
823 //Fill CaloCluster related histograms
824
1a83b960 825 Double_t tof = clus->GetTOF()*1.e9;
649b825d 826
1a83b960 827 fhLambda0 ->Fill(clus->E(),clus->GetM02());
828 fhLambda1 ->Fill(clus->E(),clus->GetM20());
829 fhDispersion ->Fill(clus->E(),clus->GetDispersion());
649b825d 830
a82b4462 831 fhClusterMaxCellDiff ->Fill(clus->E(),maxCellFraction);
832 fhClusterMaxCellECross->Fill(clus->E(),eCrossFrac);
833 fhClusterTimeEnergy ->Fill(clus->E(),tof);
649b825d 834
1a83b960 835 if(fStudyClustersAsymmetry) ClusterAsymmetryHistograms(clus,absIdMax,kTRUE);
836
649b825d 837 //Clusters in event time difference
45769d5b 838 for(Int_t iclus2 = 0; iclus2 < caloClusters->GetEntriesFast(); iclus2++ )
839 {
649b825d 840 AliVCluster* clus2 = (AliVCluster*) caloClusters->At(iclus2);
841
45769d5b 842 if( clus->GetID() == clus2->GetID() ) continue;
649b825d 843
45769d5b 844 if( clus->GetM02() > 0.01 && clus2->GetM02() > 0.01 )
a87e069d 845 {
a82b4462 846 Double_t tof2 = clus2->GetTOF()*1.e9;
649b825d 847 fhClusterPairDiffTimeE ->Fill(clus->E(), tof-tof2);
848 }
849 }
850
1a83b960 851 Int_t nModule = GetModuleNumber(clus);
852 Int_t nCaloCellsPerCluster = clus->GetNCells();
853
45769d5b 854 if(nCaloCellsPerCluster > 1)
855 {
649b825d 856 // check time of cells respect to max energy cell
857
e6fec6f5 858 if(fFillAllCellTimeHisto)
a87e069d 859 {
45769d5b 860 // Get some time averages
861 Double_t timeAverages[2] = {0.,0.};
862 CalculateAverageTime(clus, cells, timeAverages);
863
649b825d 864 fhClusterMaxCellDiffAverageTime ->Fill(clus->E(),tmax-timeAverages[0]);
649b825d 865 fhClusterMaxCellDiffWeightedTime ->Fill(clus->E(),tmax-timeAverages[1]);
649b825d 866 }
867
dbba06ca 868 for (Int_t ipos = 0; ipos < nCaloCellsPerCluster; ipos++)
869 {
649b825d 870 Int_t absId = clus->GetCellsAbsId()[ipos];
45769d5b 871 if( absId == absIdMax || cells->GetCellAmplitude(absIdMax) < 0.01 ) continue;
649b825d 872
873 Float_t frac = cells->GetCellAmplitude(absId)/cells->GetCellAmplitude(absIdMax);
874 fhClusterMaxCellCloseCellRatio->Fill(clus->E(),frac);
875 fhClusterMaxCellCloseCellDiff ->Fill(clus->E(),cells->GetCellAmplitude(absIdMax)-cells->GetCellAmplitude(absId));
876
e6fec6f5 877 if(fFillAllCellTimeHisto)
dbba06ca 878 {
649b825d 879 Double_t time = cells->GetCellTime(absId);
0cea6003 880 GetCaloUtils()->RecalibrateCellTime(time, GetCalorimeter(), absId,GetReader()->GetInputEvent()->GetBunchCrossNumber());
649b825d 881
882 Float_t diff = (tmax-time*1.0e9);
883 fhCellTimeSpreadRespectToCellMax->Fill(clus->E(), diff);
884 if(TMath::Abs(TMath::Abs(diff) > 100) && clus->E() > 1 ) fhCellIdCellLargeTimeSpread->Fill(absId);
885 }
886
887 }// fill cell-cluster histogram loop
888
889 }//check time and energy of cells respect to max energy cell if cluster of more than 1 cell
f4282b09 890
891 Float_t e = fClusterMomentum.E();
892 Float_t pt = fClusterMomentum.Pt();
893 Float_t eta = fClusterMomentum.Eta();
894 Float_t phi = fClusterMomentum.Phi();
649b825d 895 if(phi < 0) phi +=TMath::TwoPi();
896
2db10729 897 AliDebug(1,Form("cluster: E %2.3f, pT %2.3f, eta %2.3f, phi %2.3f",e,pt,eta,phi*TMath::RadToDeg()));
649b825d 898
899 fhE ->Fill(e);
900 if(nModule >=0 && nModule < fNModules) fhEMod->Fill(e,nModule);
45769d5b 901
902 fhPt ->Fill(pt);
903 fhPhi ->Fill(phi);
904 fhEta ->Fill(eta);
649b825d 905
906 if(fFillAllTH3)
907 fhEtaPhiE->Fill(eta,phi,e);
908
909 //Cells per cluster
910 fhNCellsPerCluster ->Fill(e, nCaloCellsPerCluster);
a82b4462 911
649b825d 912 //Position
45769d5b 913 if(fFillAllPosHisto2)
914 {
649b825d 915 Float_t pos[3] ;
916 clus->GetPosition(pos);
917
918 fhXE ->Fill(pos[0],e);
919 fhYE ->Fill(pos[1],e);
920 fhZE ->Fill(pos[2],e);
921 if(fFillAllTH3)
922 fhXYZ ->Fill(pos[0], pos[1],pos[2]);
923
924 fhXNCells->Fill(pos[0],nCaloCellsPerCluster);
925 fhYNCells->Fill(pos[1],nCaloCellsPerCluster);
926 fhZNCells->Fill(pos[2],nCaloCellsPerCluster);
927 Float_t rxyz = TMath::Sqrt(pos[0]*pos[0]+pos[1]*pos[1]);//+pos[2]*pos[2]);
928 fhRE ->Fill(rxyz,e);
929 fhRNCells->Fill(rxyz ,nCaloCellsPerCluster);
930 }
931
45769d5b 932 if( nModule >= 0 && nModule < fNModules ) fhNCellsPerClusterMod[nModule]->Fill(e, nCaloCellsPerCluster);
649b825d 933
934}
935
f3138ecf 936//____________________________________________________________________________
937void AliAnaCalorimeterQA::ClusterLoopHistograms(const TObjArray *caloClusters,
938 AliVCaloCells* cells)
649b825d 939{
940 // Fill clusters related histograms
649b825d 941 Int_t nLabel = 0 ;
942 Int_t *labels = 0x0;
943 Int_t nCaloClusters = caloClusters->GetEntriesFast() ;
944 Int_t nCaloClustersAccepted = 0 ;
945 Int_t nCaloCellsPerCluster = 0 ;
946 Bool_t matched = kFALSE;
947 Int_t nModule =-1 ;
948
949 // Get vertex for photon momentum calculation and event selection
950 Double_t v[3] = {0,0,0}; //vertex ;
1a83b960 951 //GetReader()->GetVertex(v);
95aee5e1 952
649b825d 953 Int_t *nClustersInModule = new Int_t[fNModules];
954 for(Int_t imod = 0; imod < fNModules; imod++ ) nClustersInModule[imod] = 0;
955
2db10729 956 AliDebug(1,Form("In %s there are %d clusters", GetCalorimeterString().Data(), nCaloClusters));
649b825d 957
958 // Loop over CaloClusters
45769d5b 959 for(Int_t iclus = 0; iclus < nCaloClusters; iclus++)
960 {
2db10729 961 AliDebug(1,Form("Cluster: %d/%d, data %d",iclus+1,nCaloClusters,GetReader()->GetDataType()));
649b825d 962
6aba6683 963 AliVCluster* clus = (AliVCluster*) caloClusters->At(iclus);
649b825d 964
965 // Get the fraction of the cluster energy that carries the cell with highest energy and its absId
966 Float_t maxCellFraction = 0.;
967 Int_t absIdMax = GetCaloUtils()->GetMaxEnergyCell(cells, clus,maxCellFraction);
968
969 //Cut on time of clusters
970 Double_t tof = clus->GetTOF()*1.e9;
45769d5b 971 if( tof < fTimeCutMin || tof > fTimeCutMax )
e6fec6f5 972 {
2db10729 973 AliDebug(1,Form("Remove cluster with TOF %f",tof));
649b825d 974 continue;
975 }
976
977 // Get cluster kinematics
f4282b09 978 clus->GetMomentum(fClusterMomentum,v);
649b825d 979
980 // Check only certain regions
981 Bool_t in = kTRUE;
f4282b09 982 if(IsFiducialCutOn()) in = GetFiducialCut()->IsInFiducialCut(fClusterMomentum.Eta(),fClusterMomentum.Phi(),GetCalorimeter()) ;
649b825d 983 if(!in) continue;
984
1a83b960 985 // MC labels
649b825d 986 nLabel = clus->GetNLabels();
987 labels = clus->GetLabels();
988
1a83b960 989 // SuperModule number of cluster
990 nModule = GetModuleNumber(clus);
991
649b825d 992 // Cells per cluster
993 nCaloCellsPerCluster = clus->GetNCells();
994
995 // Cluster mathed with track?
49b5c49b 996 matched = GetCaloPID()->IsTrackMatched(clus,GetCaloUtils(), GetReader()->GetInputEvent());
649b825d 997
649b825d 998 //Get time of max cell
999 Double_t tmax = cells->GetCellTime(absIdMax);
0cea6003 1000 GetCaloUtils()->RecalibrateCellTime(tmax, GetCalorimeter(), absIdMax,GetReader()->GetInputEvent()->GetBunchCrossNumber());
649b825d 1001 tmax*=1.e9;
1002
1a83b960 1003 // Fill histograms related to single cluster
1004
1a83b960 1005 // Fill some histograms before applying the exotic cell / bad map cut
1006 fhNCellsPerClusterNoCut ->Fill(clus->E(), nCaloCellsPerCluster);
1007 if(nModule >=0 && nModule < fNModules) fhNCellsPerClusterModNoCut[nModule]->Fill(clus->E(), nCaloCellsPerCluster);
1008
1009 fhClusterMaxCellDiffNoCut->Fill(clus->E(),maxCellFraction);
1010
1011 Float_t ampMax = cells->GetCellAmplitude(absIdMax);
0cea6003 1012 GetCaloUtils()->RecalibrateCellAmplitude(ampMax,GetCalorimeter(), absIdMax);
1a83b960 1013
f1538a5f 1014 if(fStudyExotic) ExoticHistograms(absIdMax, ampMax, clus, cells);
1015
649b825d 1016 //Check bad clusters if requested and rejection was not on
a82b4462 1017 Bool_t goodCluster = IsGoodCluster(absIdMax, cells);
f1538a5f 1018
2747966a 1019 Float_t eCrossFrac = 0;
1020 if(ampMax > 0.01) eCrossFrac = 1-GetECross(absIdMax,cells)/ampMax;
649b825d 1021
649b825d 1022 if(!goodCluster)
1a83b960 1023 {
649b825d 1024 BadClusterHistograms(clus, caloClusters, cells, absIdMax,
45769d5b 1025 maxCellFraction, eCrossFrac, tmax);
1a83b960 1026 continue;
1027 }
649b825d 1028
1029 ClusterHistograms(clus, caloClusters, cells, absIdMax,
45769d5b 1030 maxCellFraction, eCrossFrac, tmax);
649b825d 1031
1032 nCaloClustersAccepted++;
49214ef9 1033 nModule = GetModuleNumber(clus);
f4282b09 1034 if(nModule >=0 && nModule < fNModules && fClusterMomentum.E() > 2*fCellAmpMin)
701cbf54 1035 nClustersInModule[nModule]++;
1a83b960 1036
649b825d 1037 // Cluster weights
1038 if(fStudyWeight) WeightHistograms(clus, cells);
1039
1040 // Cells in cluster position
1041 if(fFillAllPosHisto) CellInClusterPositionHistograms(clus);
1042
1043 // Fill histograms related to single cluster, mc vs data
1044 Int_t mcOK = kFALSE;
1045 Int_t pdg = -1;
1046 if(IsDataMC() && nLabel > 0 && labels)
f4282b09 1047 mcOK = ClusterMCHistograms(matched, labels, nLabel, pdg);
008693e5 1048
649b825d 1049 // Matched clusters with tracks, also do some MC comparison, needs input from ClusterMCHistograms
1050 if( matched && fFillAllTMHisto)
f4282b09 1051 ClusterMatchedWithTrackHistograms(clus,mcOK,pdg);
649b825d 1052
1053 // Invariant mass
d07278cf 1054 // Try to reduce background with a mild shower shape cut and no more than 1 maxima
1055 // in cluster and remove low energy clusters
1056 if(fFillAllPi0Histo && nCaloClusters > 1 && nCaloCellsPerCluster > 1 &&
1057 GetCaloUtils()->GetNumberOfLocalMaxima(clus,cells) == 1 &&
07e4c878 1058 clus->GetM02() < 0.5 && clus->E() > fMinInvMassECut)
f4282b09 1059 InvariantMassHistograms(iclus, nModule, caloClusters,cells);
649b825d 1060
1061 }//cluster loop
1062
1063 // Number of clusters histograms
1064 if(nCaloClustersAccepted > 0) fhNClusters->Fill(nCaloClustersAccepted);
1065
1066 // Number of clusters per module
2747966a 1067 for(Int_t imod = 0; imod < fNModules; imod++ )
1068 {
2db10729 1069 AliDebug(1,Form("Module %d calo %s clusters %d", imod, GetCalorimeterString().Data(), nClustersInModule[imod]));
649b825d 1070 fhNClustersMod->Fill(nClustersInModule[imod],imod);
1071 }
1072
1073 delete [] nClustersInModule;
1074
1075}
1076
f4282b09 1077//__________________________________________________________________________________
1078Bool_t AliAnaCalorimeterQA::ClusterMCHistograms(Bool_t matched,const Int_t * labels,
1079 Int_t nLabels, Int_t & pdg )
649b825d 1080{
1081
1082 //Fill histograms only possible when simulation
1083
2747966a 1084 if(!labels || nLabels<=0)
1085 {
2db10729 1086 AliWarning(Form("Strange, labels array %p, n labels %d", labels,nLabels));
42d47cb7 1087 return kFALSE;
1088 }
1089
2db10729 1090 AliDebug(1,Form("Primaries: nlabels %d",nLabels));
649b825d 1091
f4282b09 1092 Float_t e = fClusterMomentum.E();
1093 Float_t eta = fClusterMomentum.Eta();
1094 Float_t phi = fClusterMomentum.Phi();
649b825d 1095 if(phi < 0) phi +=TMath::TwoPi();
1096
1097 AliAODMCParticle * aodprimary = 0x0;
1098 TParticle * primary = 0x0;
1099
1100 //Play with the MC stack if available
1101 Int_t label = labels[0];
1102
2747966a 1103 if(label < 0)
1104 {
2db10729 1105 AliDebug(1,Form(" *** bad label ***: label %d", label));
649b825d 1106 return kFALSE;
1107 }
1108
45769d5b 1109 Int_t pdg0 =-1; Int_t status = -1; Int_t iMother = -1; Int_t iParent = -1;
1110 Float_t vxMC = 0; Float_t vyMC = 0;
1111 Float_t eMC = 0; //Float_t ptMC= 0;
1112 Float_t phiMC = 0; Float_t etaMC = 0;
1113 Int_t charge = 0;
649b825d 1114
1115 //Check the origin.
0cea6003 1116 Int_t tag = GetMCAnalysisUtils()->CheckOrigin(labels,nLabels, GetReader(),GetCalorimeter());
649b825d 1117
d07278cf 1118 if ( GetReader()->ReadStack() &&
1119 !GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCUnknown))
1120 { //it MC stack and known tag
649b825d 1121
d07278cf 1122 if( label >= GetMCStack()->GetNtrack())
1123 {
2db10729 1124 AliDebug(1,Form("*** large label ***: label %d, n tracks %d", label, GetMCStack()->GetNtrack()));
649b825d 1125 return kFALSE;
1126 }
1127
1128 primary = GetMCStack()->Particle(label);
1129 iMother = label;
1130 pdg0 = TMath::Abs(primary->GetPdgCode());
1131 pdg = pdg0;
1132 status = primary->GetStatusCode();
1133 vxMC = primary->Vx();
1134 vyMC = primary->Vy();
1135 iParent = primary->GetFirstMother();
1136
2db10729 1137 AliDebug(1,"Cluster most contributing mother:");
1138 AliDebug(1,Form("\t Mother label %d, pdg %d, %s, status %d, parent %d",iMother, pdg0, primary->GetName(),status, iParent));
1139
649b825d 1140
1141 //Get final particle, no conversion products
d07278cf 1142 if(GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCConversion))
1143 {
649b825d 1144 //Get the parent
1145 primary = GetMCStack()->Particle(iParent);
1146 pdg = TMath::Abs(primary->GetPdgCode());
008693e5 1147
2db10729 1148 AliDebug(2,"Converted cluster!. Find before conversion:");
008693e5 1149
d07278cf 1150 while((pdg == 22 || pdg == 11) && status != 1)
1151 {
008693e5 1152 Int_t iMotherOrg = iMother;
649b825d 1153 iMother = iParent;
1154 primary = GetMCStack()->Particle(iMother);
1155 status = primary->GetStatusCode();
649b825d 1156 pdg = TMath::Abs(primary->GetPdgCode());
008693e5 1157 iParent = primary->GetFirstMother();
1158
1159 // If gone too back and non stable, assign the decay photon/electron
1160 // there are other possible decays, ignore them for the moment
1161 if(pdg==111 || pdg==221)
1162 {
1163 primary = GetMCStack()->Particle(iMotherOrg);
1164 break;
1165 }
1166
1167 if( iParent < 0 )
1168 {
1169 iParent = iMother;
008693e5 1170 break;
1171 }
1172
2db10729 1173 AliDebug(2,Form("\t pdg %d, index %d, %s, status %d",pdg, iMother, primary->GetName(),status));
649b825d 1174 }
008693e5 1175
2db10729 1176 AliDebug(1,"Converted Cluster mother before conversion:");
1177 AliDebug(1,Form("\t Mother label %d, pdg %d, %s, status %d, parent %d",iMother, pdg, primary->GetName(), status, iParent));
649b825d 1178
1179 }
1180
1181 //Overlapped pi0 (or eta, there will be very few), get the meson
1182 if(GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCPi0) ||
d07278cf 1183 GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCEta))
1184 {
2db10729 1185 AliDebug(2,"Overlapped Meson decay!, Find it:");
008693e5 1186
d07278cf 1187 while(pdg != 111 && pdg != 221)
008693e5 1188 {
1189 //printf("iMother %d, pdg %d, iParent %d, pdg %d\n",iMother,pdg,iParent,GetMCStack()->Particle(iParent)->GetPdgCode());
649b825d 1190 iMother = iParent;
1191 primary = GetMCStack()->Particle(iMother);
1192 status = primary->GetStatusCode();
649b825d 1193 pdg = TMath::Abs(primary->GetPdgCode());
008693e5 1194 iParent = primary->GetFirstMother();
1195
2db10729 1196 if( iParent < 0 ) break;
008693e5 1197
2db10729 1198 AliDebug(2,Form("\t pdg %d, %s, index %d",pdg, primary->GetName(),iMother));
008693e5 1199
d07278cf 1200 if(iMother==-1)
1201 {
2db10729 1202 AliWarning("Tagged as Overlapped photon but meson not found, why?");
649b825d 1203 //break;
1204 }
1205 }
008693e5 1206
2db10729 1207 AliDebug(2,Form("Overlapped %s decay, label %d",primary->GetName(),iMother));
649b825d 1208 }
1209
1210 eMC = primary->Energy();
ecdde216 1211 //ptMC = primary->Pt();
649b825d 1212 phiMC = primary->Phi();
1213 etaMC = primary->Eta();
1214 pdg = TMath::Abs(primary->GetPdgCode());
1215 charge = (Int_t) TDatabasePDG::Instance()->GetParticle(pdg)->Charge();
1216
1217 }
d07278cf 1218 else if( GetReader()->ReadAODMCParticles() &&
1219 !GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCUnknown))
1220 {//it MC AOD and known tag
649b825d 1221 //Get the list of MC particles
2644ead9 1222 if(!GetReader()->GetAODMCParticles())
649b825d 1223 AliFatal("MCParticles not available!");
1224
2644ead9 1225 aodprimary = (AliAODMCParticle*) (GetReader()->GetAODMCParticles())->At(label);
649b825d 1226 iMother = label;
1227 pdg0 = TMath::Abs(aodprimary->GetPdgCode());
1228 pdg = pdg0;
1229 status = aodprimary->IsPrimary();
1230 vxMC = aodprimary->Xv();
1231 vyMC = aodprimary->Yv();
1232 iParent = aodprimary->GetMother();
1233
2db10729 1234 AliDebug(1,"Cluster most contributing mother:");
1235 AliDebug(1,Form("\t Mother label %d, pdg %d, Primary? %d, Physical Primary? %d, parent %d",
1236 iMother, pdg0, aodprimary->IsPrimary(), aodprimary->IsPhysicalPrimary(), iParent));
649b825d 1237
1238 //Get final particle, no conversion products
45769d5b 1239 if( GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCConversion) )
008693e5 1240 {
2db10729 1241 AliDebug(2,"Converted cluster!. Find before conversion:");
649b825d 1242 //Get the parent
2644ead9 1243 aodprimary = (AliAODMCParticle*)(GetReader()->GetAODMCParticles())->At(iParent);
649b825d 1244 pdg = TMath::Abs(aodprimary->GetPdgCode());
d07278cf 1245 while ((pdg == 22 || pdg == 11) && !aodprimary->IsPhysicalPrimary())
1246 {
008693e5 1247 Int_t iMotherOrg = iMother;
649b825d 1248 iMother = iParent;
2644ead9 1249 aodprimary = (AliAODMCParticle*)(GetReader()->GetAODMCParticles())->At(iMother);
649b825d 1250 status = aodprimary->IsPrimary();
1251 iParent = aodprimary->GetMother();
1252 pdg = TMath::Abs(aodprimary->GetPdgCode());
008693e5 1253
1254 // If gone too back and non stable, assign the decay photon/electron
1255 // there are other possible decays, ignore them for the moment
45769d5b 1256 if( pdg == 111 || pdg == 221 )
008693e5 1257 {
2644ead9 1258 aodprimary = (AliAODMCParticle*)(GetReader()->GetAODMCParticles())->At(iMotherOrg);
008693e5 1259 break;
1260 }
1261
45769d5b 1262 if( iParent < 0 )
008693e5 1263 {
1264 iParent = iMother;
1265 break;
1266 }
1267
2db10729 1268 AliDebug(2,Form("\t pdg %d, index %d, Primary? %d, Physical Primary? %d",
1269 pdg, iMother, aodprimary->IsPrimary(), aodprimary->IsPhysicalPrimary()));
649b825d 1270 }
1271
2db10729 1272 AliDebug(1,"Converted Cluster mother before conversion:");
1273 AliDebug(1,Form("\t Mother label %d, pdg %d, parent %d, Primary? %d, Physical Primary? %d",
1274 iMother, pdg, iParent, aodprimary->IsPrimary(), aodprimary->IsPhysicalPrimary()));
1275
649b825d 1276 }
1277
1278 //Overlapped pi0 (or eta, there will be very few), get the meson
1279 if(GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCPi0) ||
2747966a 1280 GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCEta))
1281 {
2db10729 1282 AliDebug(2,Form("Overlapped Meson decay!, Find it: PDG %d, mom %d",pdg, iMother));
45769d5b 1283
008693e5 1284 while(pdg != 111 && pdg != 221)
1285 {
649b825d 1286 iMother = iParent;
2644ead9 1287 aodprimary = (AliAODMCParticle*)(GetReader()->GetAODMCParticles())->At(iMother);
649b825d 1288 status = aodprimary->IsPrimary();
1289 iParent = aodprimary->GetMother();
1290 pdg = TMath::Abs(aodprimary->GetPdgCode());
008693e5 1291
45769d5b 1292 if( iParent < 0 ) break;
649b825d 1293
2db10729 1294 AliDebug(2,Form("\t pdg %d, index %d",pdg, iMother));
649b825d 1295
2747966a 1296 if(iMother==-1)
1297 {
2db10729 1298 AliWarning("Tagged as Overlapped photon but meson not found, why?");
649b825d 1299 //break;
1300 }
1301 }
1302
2db10729 1303 AliDebug(2,Form("Overlapped %s decay, label %d",aodprimary->GetName(),iMother));
649b825d 1304 }
1305
1306 status = aodprimary->IsPrimary();
1307 eMC = aodprimary->E();
ecdde216 1308 //ptMC = aodprimary->Pt();
649b825d 1309 phiMC = aodprimary->Phi();
1310 etaMC = aodprimary->Eta();
1311 pdg = TMath::Abs(aodprimary->GetPdgCode());
1312 charge = aodprimary->Charge();
1313
1314 }
1315
1316 //Float_t vz = primary->Vz();
1317 Float_t rVMC = TMath::Sqrt(vxMC*vxMC + vyMC*vyMC);
45769d5b 1318 if( ( pdg == 22 || TMath::Abs(pdg) == 11 ) && status != 1 )
2747966a 1319 {
649b825d 1320 fhEMVxyz ->Fill(vxMC,vyMC);//,vz);
1321 fhEMR ->Fill(e,rVMC);
1322 }
1323
1324 //printf("reco e %f, pt %f, phi %f, eta %f \n", e, pt, phi, eta);
1325 //printf("prim e %f, pt %f, phi %f, eta %f \n", eMC,ptMC,phiMC ,etaMC );
1326 //printf("vertex: vx %f, vy %f, vz %f, r %f \n", vxMC, vyMC, vz, r);
1327
1328 //Overlapped pi0 (or eta, there will be very few)
6aba6683 1329 Int_t mcIndex = -1;
1330 if ( GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCPi0 ) )
2747966a 1331 {
6aba6683 1332 mcIndex = kmcPi0;
649b825d 1333 }//Overlapped pizero decay
6aba6683 1334 else if( GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCEta ) )
2747966a 1335 {
6aba6683 1336 mcIndex = kmcEta;
649b825d 1337 }//Overlapped eta decay
6aba6683 1338 else if( GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCPhoton ) )
2747966a 1339 {
6aba6683 1340 if( GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCConversion))
1341 mcIndex = kmcPhotonConv ;
1342 else
1343 mcIndex = kmcPhoton ;
649b825d 1344 }//photon
6aba6683 1345 else if( GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCElectron) )
2747966a 1346 {
6aba6683 1347 mcIndex = kmcElectron;
1348 fhEMVxyz ->Fill(vxMC,vyMC);//,vz);
1349 fhEMR ->Fill(e,rVMC);
649b825d 1350 }
2747966a 1351 else if(charge == 0)
1352 {
6aba6683 1353 mcIndex = kmcNeHadron;
649b825d 1354 fhHaVxyz ->Fill(vxMC,vyMC);//,vz);
1355 fhHaR ->Fill(e,rVMC);
1356 }
2747966a 1357 else if(charge!=0)
1358 {
6aba6683 1359 mcIndex = kmcChHadron;
1360 fhHaVxyz ->Fill(vxMC,vyMC);//,vz);
1361 fhHaR ->Fill(e,rVMC);
1362 }
1363
1364 //printf("mc index %d\n",mcIndex);
1365
1366 if( mcIndex >= 0 && mcIndex < 7 )
1367 {
1368 fhRecoMCE [mcIndex][(matched)] ->Fill(e,eMC);
1369 if(e > 0.5 && eMC > 0.5) fhRecoMCEta[mcIndex][(matched)]->Fill(eta,etaMC);
1370 if(e > 0.5 && eMC > 0.5) fhRecoMCPhi[mcIndex][(matched)]->Fill(phi,phiMC);
1371 if(eMC > 0) fhRecoMCRatioE [mcIndex][(matched)]->Fill(e,e/eMC);
1372 fhRecoMCDeltaE [mcIndex][(matched)]->Fill(e,eMC-e);
1373 fhRecoMCDeltaPhi[mcIndex][(matched)]->Fill(e,phiMC-phi);
1374 fhRecoMCDeltaEta[mcIndex][(matched)]->Fill(e,etaMC-eta);
649b825d 1375 }
1376
45769d5b 1377 if( primary || aodprimary ) return kTRUE ;
1378 else return kFALSE;
649b825d 1379
1380}
1381
f4282b09 1382//_________________________________________________________________________________________________________
1383void AliAnaCalorimeterQA::ClusterMatchedWithTrackHistograms(AliVCluster *clus, Bool_t okPrimary, Int_t pdg)
649b825d 1384{
1385 //Histograms for clusters matched with tracks
42d47cb7 1386
f4282b09 1387 Float_t e = fClusterMomentum.E();
1388 Float_t pt = fClusterMomentum.Pt();
1389 Float_t eta = fClusterMomentum.Eta();
1390 Float_t phi = fClusterMomentum.Phi();
649b825d 1391 if(phi < 0) phi +=TMath::TwoPi();
45769d5b 1392
1393 fhECharged ->Fill(e);
1394 fhPtCharged ->Fill(pt);
1395 fhPhiCharged ->Fill(phi);
1396 fhEtaCharged ->Fill(eta);
a82b4462 1397
42d47cb7 1398 //Study the track and matched cluster if track exists.
1399
4bfeae64 1400 AliVTrack *track = GetCaloUtils()->GetMatchedTrack(clus, GetReader()->GetInputEvent());
649b825d 1401
1402 if(!track) return ;
649b825d 1403
a87e069d 1404 Double_t tpt = track->Pt();
1405 Double_t tmom = track->P();
1406 Double_t dedx = track->GetTPCsignal();
1407 Int_t nITS = track->GetNcls(0);
1408 Int_t nTPC = track->GetNcls(1);
653aed3c 1409 Bool_t positive = kFALSE;
1410 if(track) positive = (track->Charge()>0);
a87e069d 1411
1412 // Residuals
1413 Float_t deta = clus->GetTrackDz();
1414 Float_t dphi = clus->GetTrackDx();
1415 Double_t dR = TMath::Sqrt(dphi*dphi + deta*deta);
1416
653aed3c 1417 if(TMath::Abs(dphi) < 999)
1418 {
1419 fhTrackMatchedDEta->Fill(e,deta);
1420 fhTrackMatchedDPhi->Fill(e,dphi);
1421 if(e > 0.5) fhTrackMatchedDEtaDPhi->Fill(deta,dphi);
1422
1423 if(track && positive)
1424 {
653aed3c 1425 fhTrackMatchedDEtaPos->Fill(e,deta);
1426 fhTrackMatchedDPhiPos->Fill(e,dphi);
1427 if(e > 0.5) fhTrackMatchedDEtaDPhiPos->Fill(deta,dphi);
1428 }
1429 }
a87e069d 1430
653aed3c 1431 Double_t eOverP = e/tmom;
d55bb5e1 1432 fh1EOverP->Fill(tpt, eOverP);
653aed3c 1433 if(dR < 0.02)
1434 {
a054a582 1435 fh1EOverPR02->Fill(tpt,eOverP);
1436 if(dedx > 60 && dedx < 100) fh1EleEOverP->Fill(tpt,eOverP);
1437 }
a87e069d 1438
1439 fh2dR->Fill(e,dR);
1440 fh2MatchdEdx->Fill(tmom,dedx);
1441
1442 if(IsDataMC() && okPrimary)
1443 {
1444 Double_t charge = TDatabasePDG::Instance()->GetParticle(pdg)->Charge();
649b825d 1445
a87e069d 1446 if(TMath::Abs(pdg) == 11)
1447 {
d55bb5e1 1448 fhMCEle1EOverP->Fill(tpt,eOverP);
a87e069d 1449 fhMCEle1dR->Fill(dR);
1450 fhMCEle2MatchdEdx->Fill(tmom,dedx);
653aed3c 1451 if(dR < 0.02)
1452 {
1453 fhMCEle1EOverPR02->Fill(tpt,eOverP);
1454 if(dedx > 60 && dedx < 100) fhMCEle1EleEOverP->Fill(tpt,eOverP);
a054a582 1455 }
a87e069d 1456 }
1457 else if(charge!=0)
1458 {
d55bb5e1 1459 fhMCChHad1EOverP->Fill(tpt,eOverP);
a87e069d 1460 fhMCChHad1dR->Fill(dR);
1461 fhMCChHad2MatchdEdx->Fill(tmom,dedx);
653aed3c 1462 if(dR < 0.02)
1463 {
1464 fhMCChHad1EOverPR02->Fill(tpt,eOverP);
1465 if(dedx > 60 && dedx < 100) fhMCChHad1EleEOverP->Fill(tpt,eOverP);
a054a582 1466 }
a87e069d 1467 }
1468 else if(charge == 0)
1469 {
d55bb5e1 1470 fhMCNeutral1EOverP->Fill(tpt,eOverP);
a87e069d 1471 fhMCNeutral1dR->Fill(dR);
1472 fhMCNeutral2MatchdEdx->Fill(tmom,dedx);
653aed3c 1473 if(dR < 0.02)
1474 {
1475 fhMCNeutral1EOverPR02->Fill(tpt,eOverP);
1476 if(dedx > 60 && dedx < 100) fhMCNeutral1EleEOverP->Fill(tpt,eOverP);
a054a582 1477 }
649b825d 1478 }
a87e069d 1479 }//DataMC
1480
45769d5b 1481 if(dR < 0.02 && eOverP > 0.6 && eOverP < 1.2
a87e069d 1482 && clus->GetNCells() > 1 && nITS > 3 && nTPC > 20)
1483 {
1484 fh2EledEdx->Fill(tmom,dedx);
649b825d 1485 }
a87e069d 1486
649b825d 1487}
1488
1489//___________________________________
1490void AliAnaCalorimeterQA::Correlate()
1491{
1492 // Correlate information from PHOS and EMCAL and with V0 and track multiplicity
1493
95aee5e1 1494 //Clusters arrays
649b825d 1495 TObjArray * caloClustersEMCAL = GetEMCALClusters();
1496 TObjArray * caloClustersPHOS = GetPHOSClusters();
1497
95aee5e1 1498 if(!caloClustersEMCAL || !caloClustersPHOS)
1499 {
2db10729 1500 AliDebug(1,Form("PHOS (%p) or EMCAL (%p) clusters array not available, do not correlate",caloClustersPHOS,caloClustersEMCAL));
95aee5e1 1501 return ;
1502 }
1503
1504 //Cells arrays
1505 AliVCaloCells * cellsEMCAL = GetEMCALCells();
1506 AliVCaloCells * cellsPHOS = GetPHOSCells();
1507
1508 if(!cellsEMCAL || !cellsPHOS)
1509 {
2db10729 1510 AliDebug(1,Form("PHOS (%p) or EMCAL (%p) cells array ot available, do not correlate",cellsPHOS,cellsEMCAL));
95aee5e1 1511 return ;
1512 }
1513
1514 // Clusters parameters
649b825d 1515 Int_t nclEMCAL = caloClustersEMCAL->GetEntriesFast();
1516 Int_t nclPHOS = caloClustersPHOS ->GetEntriesFast();
1517
653aed3c 1518 Float_t cen = GetEventCentrality();
1519 Float_t ep = GetEventPlaneAngle();
1520
649b825d 1521 Float_t sumClusterEnergyEMCAL = 0;
1522 Float_t sumClusterEnergyPHOS = 0;
1523 Int_t iclus = 0;
1524 for(iclus = 0 ; iclus < caloClustersEMCAL->GetEntriesFast() ; iclus++)
1525 sumClusterEnergyEMCAL += ((AliVCluster*)caloClustersEMCAL->At(iclus))->E();
1526 for(iclus = 0 ; iclus < caloClustersPHOS->GetEntriesFast(); iclus++)
1527 sumClusterEnergyPHOS += ((AliVCluster*)caloClustersPHOS->At(iclus))->E();
1528
95aee5e1 1529 //Cells parameters
649b825d 1530 Int_t ncellsEMCAL = cellsEMCAL->GetNumberOfCells();
1531 Int_t ncellsPHOS = cellsPHOS ->GetNumberOfCells();
1532
1533 Float_t sumCellEnergyEMCAL = 0;
1534 Float_t sumCellEnergyPHOS = 0;
1535 Int_t icell = 0;
1536 for(icell = 0 ; icell < cellsEMCAL->GetNumberOfCells() ; icell++)
1537 sumCellEnergyEMCAL += cellsEMCAL->GetAmplitude(icell);
1538 for(icell = 0 ; icell < cellsPHOS->GetNumberOfCells(); icell++)
1539 sumCellEnergyPHOS += cellsPHOS->GetAmplitude(icell);
1540
1541
1542 //Fill Histograms
1543 fhCaloCorrNClusters->Fill(nclEMCAL,nclPHOS);
1544 fhCaloCorrEClusters->Fill(sumClusterEnergyEMCAL,sumClusterEnergyPHOS);
1545 fhCaloCorrNCells ->Fill(ncellsEMCAL,ncellsPHOS);
1546 fhCaloCorrECells ->Fill(sumCellEnergyEMCAL,sumCellEnergyPHOS);
1547
1548 Int_t v0S = GetV0Signal(0)+GetV0Signal(1);
1549 Int_t v0M = GetV0Multiplicity(0)+GetV0Multiplicity(1);
1550 Int_t trM = GetTrackMultiplicity();
f3da05f3 1551 if(GetCalorimeter()==kPHOS)
653aed3c 1552 {
649b825d 1553 fhCaloV0MCorrNClusters ->Fill(v0M,nclPHOS);
1554 fhCaloV0MCorrEClusters ->Fill(v0M,sumClusterEnergyPHOS);
1555 fhCaloV0MCorrNCells ->Fill(v0M,ncellsPHOS);
1556 fhCaloV0MCorrECells ->Fill(v0M,sumCellEnergyPHOS);
1557
1558 fhCaloV0SCorrNClusters ->Fill(v0S,nclPHOS);
1559 fhCaloV0SCorrEClusters ->Fill(v0S,sumClusterEnergyPHOS);
1560 fhCaloV0SCorrNCells ->Fill(v0S,ncellsPHOS);
1561 fhCaloV0SCorrECells ->Fill(v0S,sumCellEnergyPHOS);
1562
1563 fhCaloTrackMCorrNClusters->Fill(trM,nclPHOS);
1564 fhCaloTrackMCorrEClusters->Fill(trM,sumClusterEnergyPHOS);
1565 fhCaloTrackMCorrNCells ->Fill(trM,ncellsPHOS);
1566 fhCaloTrackMCorrECells ->Fill(trM,sumCellEnergyPHOS);
653aed3c 1567
1568 fhCaloCenNClusters ->Fill(cen,nclPHOS);
1569 fhCaloCenEClusters ->Fill(cen,sumClusterEnergyPHOS);
1570 fhCaloCenNCells ->Fill(cen,ncellsPHOS);
1571 fhCaloCenECells ->Fill(cen,sumCellEnergyPHOS);
1572
1573 fhCaloEvPNClusters ->Fill(ep ,nclPHOS);
1574 fhCaloEvPEClusters ->Fill(ep ,sumClusterEnergyPHOS);
1575 fhCaloEvPNCells ->Fill(ep ,ncellsPHOS);
1576 fhCaloEvPECells ->Fill(ep ,sumCellEnergyPHOS);
649b825d 1577 }
653aed3c 1578 else
1579 {
649b825d 1580 fhCaloV0MCorrNClusters ->Fill(v0M,nclEMCAL);
1581 fhCaloV0MCorrEClusters ->Fill(v0M,sumClusterEnergyEMCAL);
1582 fhCaloV0MCorrNCells ->Fill(v0M,ncellsEMCAL);
1583 fhCaloV0MCorrECells ->Fill(v0M,sumCellEnergyEMCAL);
1584
1585 fhCaloV0SCorrNClusters ->Fill(v0S,nclEMCAL);
1586 fhCaloV0SCorrEClusters ->Fill(v0S,sumClusterEnergyEMCAL);
1587 fhCaloV0SCorrNCells ->Fill(v0S,ncellsEMCAL);
1588 fhCaloV0SCorrECells ->Fill(v0S,sumCellEnergyEMCAL);
1589
1590 fhCaloTrackMCorrNClusters->Fill(trM,nclEMCAL);
1591 fhCaloTrackMCorrEClusters->Fill(trM,sumClusterEnergyEMCAL);
1592 fhCaloTrackMCorrNCells ->Fill(trM,ncellsEMCAL);
1593 fhCaloTrackMCorrECells ->Fill(trM,sumCellEnergyEMCAL);
653aed3c 1594
1595 fhCaloCenNClusters ->Fill(cen,nclEMCAL);
1596 fhCaloCenEClusters ->Fill(cen,sumClusterEnergyEMCAL);
1597 fhCaloCenNCells ->Fill(cen,ncellsEMCAL);
1598 fhCaloCenECells ->Fill(cen,sumCellEnergyEMCAL);
1599
1600 fhCaloEvPNClusters ->Fill(ep ,nclEMCAL);
1601 fhCaloEvPEClusters ->Fill(ep ,sumClusterEnergyEMCAL);
1602 fhCaloEvPNCells ->Fill(ep ,ncellsEMCAL);
1603 fhCaloEvPECells ->Fill(ep ,sumCellEnergyEMCAL);
649b825d 1604 }
1605
2db10729 1606 AliDebug(1,"Correlate():");
1607 AliDebug(1,Form("\t EMCAL: N cells %d, N clusters %d, summed E cells %f, summed E clusters %f",
1608 ncellsEMCAL,nclEMCAL, sumCellEnergyEMCAL,sumClusterEnergyEMCAL));
1609 AliDebug(1,Form("\t PHOS : N cells %d, N clusters %d, summed E cells %f, summed E clusters %f",
1610 ncellsPHOS,nclPHOS,sumCellEnergyPHOS,sumClusterEnergyPHOS));
1611 AliDebug(1,Form("\t V0 : Signal %d, Multiplicity %d, Track Multiplicity %d", v0S,v0M,trM));
1612 AliDebug(1,Form("\t centrality : %f, Event plane angle %f", cen,ep));
1613
649b825d 1614}
1615
1616//__________________________________________________
1617TObjString * AliAnaCalorimeterQA::GetAnalysisCuts()
1618{
1619 //Save parameters used for analysis
1620 TString parList ; //this will be list of parameters used for this analysis.
1621 const Int_t buffersize = 255;
1622 char onePar[buffersize] ;
1623
2db10729 1624 snprintf(onePar,buffersize,"--- AliAnaCalorimeterQA ---:") ;
649b825d 1625 parList+=onePar ;
2db10729 1626 snprintf(onePar,buffersize,"Calorimeter: %s;",GetCalorimeterString().Data()) ;
649b825d 1627 parList+=onePar ;
2db10729 1628 snprintf(onePar,buffersize,"Time Cut : %2.2f < T < %2.2f ns;",fTimeCutMin, fTimeCutMax) ;
649b825d 1629 parList+=onePar ;
2db10729 1630 snprintf(onePar,buffersize,"PHOS Cell Amplitude > %2.2f GeV, EMCAL Cell Amplitude > %2.2f GeV;",fPHOSCellAmpMin, fEMCALCellAmpMin) ;
649b825d 1631 parList+=onePar ;
2db10729 1632 snprintf(onePar,buffersize,"Inv. Mass E1, E2 > %2.2f GeV;",fMinInvMassECut) ;
07e4c878 1633 parList+=onePar ;
1634
649b825d 1635 //Get parameters set in base class.
1636 //parList += GetBaseParametersList() ;
1637
1638 //Get parameters set in FiducialCut class (not available yet)
1639 //parlist += GetFidCut()->GetFidCutParametersList()
1640
1641 return new TObjString(parList) ;
1642}
1643
b94e038e 1644//_________________________________________________________________________________
1645void AliAnaCalorimeterQA::ExoticHistograms(Int_t absIdMax, Float_t ampMax,
f1538a5f 1646 AliVCluster *clus, AliVCaloCells* cells)
1647{
1648 // Calculate weights
1649
1650 if(ampMax < 0.01)
1651 {
2db10729 1652 AliDebug(1,Form("Low amplitude energy %f",ampMax));
f1538a5f 1653 return;
1654 }
1655
cc11121e 1656 Float_t l0 = clus->GetM02();
765206a5 1657 Float_t l1 = clus->GetM20();
cc11121e 1658 Float_t en = clus->E();
1659 Int_t nc = clus->GetNCells();
765206a5 1660 Double_t tmax = clus->GetTOF()*1.e9; // recalibrated elsewhere
1661
1662 Float_t eCrossFrac = 1-GetECross(absIdMax,cells, 10000000)/ampMax;
1663
1664 if(en > 5)
1665 {
1666 fhExoL0ECross->Fill(eCrossFrac,l0);
1667 fhExoL1ECross->Fill(eCrossFrac,l1);
1668 }
f1538a5f 1669
1670 for(Int_t ie = 0; ie < fExoNECrossCuts; ie++)
1671 {
1672 for(Int_t idt = 0; idt < fExoNDTimeCuts; idt++)
1673 {
765206a5 1674 eCrossFrac = 1-GetECross(absIdMax,cells, fExoDTimeCuts[idt])/ampMax;
f1538a5f 1675
1676 if(eCrossFrac > fExoECrossCuts[ie])
1677 {
1678 //Exotic
1679 fhExoL0 [ie][idt]->Fill(en,l0 );
765206a5 1680 fhExoL1 [ie][idt]->Fill(en,l1 );
1681 fhExoTime [ie][idt]->Fill(en,tmax);
1682
1683 if(en > 5)
1684 {
1685 fhExoL0NCell[ie][idt]->Fill(nc,l0);
1686 fhExoL1NCell[ie][idt]->Fill(nc,l1);
1687 }
f1538a5f 1688
1689 // Diff time, do for one cut in e cross
1690 if(ie == 0)
1691 {
1692 for (Int_t icell = 0; icell < clus->GetNCells(); icell++)
1693 {
1694 Int_t absId = clus->GetCellsAbsId()[icell];
1695 Double_t time = cells->GetCellTime(absId);
0cea6003 1696 GetCaloUtils()->RecalibrateCellTime(time, GetCalorimeter(), absId,GetReader()->GetInputEvent()->GetBunchCrossNumber());
f1538a5f 1697
1698 Float_t diff = (tmax-time)*1e9;
1699 fhExoDTime[idt]->Fill(en, diff);
1700 }
1701 }
1702 }
1703 else
1704 {
1705 fhExoECross[ie][idt]->Fill(en,eCrossFrac);
cc11121e 1706 fhExoNCell [ie][idt]->Fill(en,nc);
f1538a5f 1707 }
1708 } // D time cut loop
1709 } // e cross cut loop
1710}
1711
649b825d 1712//____________________________________________________
1713TList * AliAnaCalorimeterQA::GetCreateOutputObjects()
1714{
1715 // Create histograms to be saved in output file and
1716 // store them in outputContainer
1717
1718 TList * outputContainer = new TList() ;
1719 outputContainer->SetName("QAHistos") ;
1720
9f48b3f0 1721 // Init the number of modules, set in the class AliCalorimeterUtils
1722 fNModules = GetCaloUtils()->GetNumberOfSuperModulesUsed();
f3da05f3 1723 if(GetCalorimeter()==kPHOS && fNModules > 4) fNModules = 4;
9f48b3f0 1724
649b825d 1725 //Histograms
745913ae 1726 Int_t nptbins = GetHistogramRanges()->GetHistoPtBins(); Float_t ptmax = GetHistogramRanges()->GetHistoPtMax(); Float_t ptmin = GetHistogramRanges()->GetHistoPtMin();
1727 Int_t nfineptbins = GetHistogramRanges()->GetHistoFinePtBins(); Float_t ptfinemax = GetHistogramRanges()->GetHistoFinePtMax(); Float_t ptfinemin = GetHistogramRanges()->GetHistoFinePtMin();
1728 Int_t nphibins = GetHistogramRanges()->GetHistoPhiBins(); Float_t phimax = GetHistogramRanges()->GetHistoPhiMax(); Float_t phimin = GetHistogramRanges()->GetHistoPhiMin();
1729 Int_t netabins = GetHistogramRanges()->GetHistoEtaBins(); Float_t etamax = GetHistogramRanges()->GetHistoEtaMax(); Float_t etamin = GetHistogramRanges()->GetHistoEtaMin();
1730 Int_t nmassbins = GetHistogramRanges()->GetHistoMassBins(); Float_t massmax = GetHistogramRanges()->GetHistoMassMax(); Float_t massmin = GetHistogramRanges()->GetHistoMassMin();
1731 Int_t nasymbins = GetHistogramRanges()->GetHistoAsymmetryBins(); Float_t asymmax = GetHistogramRanges()->GetHistoAsymmetryMax(); Float_t asymmin = GetHistogramRanges()->GetHistoAsymmetryMin();
ec58c056 1732 Int_t nPoverEbins = GetHistogramRanges()->GetHistoPOverEBins(); Float_t eOverPmax = GetHistogramRanges()->GetHistoPOverEMax(); Float_t eOverPmin = GetHistogramRanges()->GetHistoPOverEMin();
745913ae 1733 Int_t ndedxbins = GetHistogramRanges()->GetHistodEdxBins(); Float_t dedxmax = GetHistogramRanges()->GetHistodEdxMax(); Float_t dedxmin = GetHistogramRanges()->GetHistodEdxMin();
1734 Int_t ndRbins = GetHistogramRanges()->GetHistodRBins(); Float_t dRmax = GetHistogramRanges()->GetHistodRMax(); Float_t dRmin = GetHistogramRanges()->GetHistodRMin();
1735 Int_t ntimebins = GetHistogramRanges()->GetHistoTimeBins(); Float_t timemax = GetHistogramRanges()->GetHistoTimeMax(); Float_t timemin = GetHistogramRanges()->GetHistoTimeMin();
1736 Int_t nclbins = GetHistogramRanges()->GetHistoNClustersBins(); Int_t nclmax = GetHistogramRanges()->GetHistoNClustersMax(); Int_t nclmin = GetHistogramRanges()->GetHistoNClustersMin();
1737 Int_t ncebins = GetHistogramRanges()->GetHistoNCellsBins(); Int_t ncemax = GetHistogramRanges()->GetHistoNCellsMax(); Int_t ncemin = GetHistogramRanges()->GetHistoNCellsMin();
1738 Int_t nceclbins = GetHistogramRanges()->GetHistoNClusterCellBins(); Int_t nceclmax = GetHistogramRanges()->GetHistoNClusterCellMax(); Int_t nceclmin = GetHistogramRanges()->GetHistoNClusterCellMin();
1739 Int_t nvdistbins = GetHistogramRanges()->GetHistoVertexDistBins(); Float_t vdistmax = GetHistogramRanges()->GetHistoVertexDistMax(); Float_t vdistmin = GetHistogramRanges()->GetHistoVertexDistMin();
1740 Int_t rbins = GetHistogramRanges()->GetHistoRBins(); Float_t rmax = GetHistogramRanges()->GetHistoRMax(); Float_t rmin = GetHistogramRanges()->GetHistoRMin();
1741 Int_t xbins = GetHistogramRanges()->GetHistoXBins(); Float_t xmax = GetHistogramRanges()->GetHistoXMax(); Float_t xmin = GetHistogramRanges()->GetHistoXMin();
1742 Int_t ybins = GetHistogramRanges()->GetHistoYBins(); Float_t ymax = GetHistogramRanges()->GetHistoYMax(); Float_t ymin = GetHistogramRanges()->GetHistoYMin();
1743 Int_t zbins = GetHistogramRanges()->GetHistoZBins(); Float_t zmax = GetHistogramRanges()->GetHistoZMax(); Float_t zmin = GetHistogramRanges()->GetHistoZMin();
1744 Int_t ssbins = GetHistogramRanges()->GetHistoShowerShapeBins(); Float_t ssmax = GetHistogramRanges()->GetHistoShowerShapeMax(); Float_t ssmin = GetHistogramRanges()->GetHistoShowerShapeMin();
1745 Int_t tdbins = GetHistogramRanges()->GetHistoDiffTimeBins() ; Float_t tdmax = GetHistogramRanges()->GetHistoDiffTimeMax(); Float_t tdmin = GetHistogramRanges()->GetHistoDiffTimeMin();
1746
1747 Int_t nv0sbins = GetHistogramRanges()->GetHistoV0SignalBins(); Int_t nv0smax = GetHistogramRanges()->GetHistoV0SignalMax(); Int_t nv0smin = GetHistogramRanges()->GetHistoV0SignalMin();
1748 Int_t nv0mbins = GetHistogramRanges()->GetHistoV0MultiplicityBins(); Int_t nv0mmax = GetHistogramRanges()->GetHistoV0MultiplicityMax(); Int_t nv0mmin = GetHistogramRanges()->GetHistoV0MultiplicityMin();
1749 Int_t ntrmbins = GetHistogramRanges()->GetHistoTrackMultiplicityBins(); Int_t ntrmmax = GetHistogramRanges()->GetHistoTrackMultiplicityMax(); Int_t ntrmmin = GetHistogramRanges()->GetHistoTrackMultiplicityMin();
649b825d 1750
1751 //EMCAL
1752 fNMaxCols = 48;
1753 fNMaxRows = 24;
1754 fNRCU = 2 ;
1755 //PHOS
f3da05f3 1756 if(GetCalorimeter()==kPHOS)
f1538a5f 1757 {
649b825d 1758 fNMaxCols = 56;
1759 fNMaxRows = 64;
1760 fNRCU = 4 ;
1761 }
1762
78451bcd 1763 fhE = new TH1F ("hE","#it{E} reconstructed clusters ", nptbins*5,ptmin,ptmax*5);
1764 fhE->SetXTitle("#it{E} (GeV)");
649b825d 1765 outputContainer->Add(fhE);
1766
78451bcd 1767 fhPt = new TH1F ("hPt","#it{p}_{T} reconstructed clusters", nptbins,ptmin,ptmax);
1768 fhPt->SetXTitle("#it{p}_{T} (GeV/#it{c})");
45769d5b 1769 outputContainer->Add(fhPt);
1770
1771 fhPhi = new TH1F ("hPhi","#phi reconstructed clusters ",nphibins,phimin,phimax);
1772 fhPhi->SetXTitle("#phi (rad)");
1773 outputContainer->Add(fhPhi);
1774
1775 fhEta = new TH1F ("hEta","#eta reconstructed clusters ",netabins,etamin,etamax);
1776 fhEta->SetXTitle("#eta ");
1777 outputContainer->Add(fhEta);
1778
649b825d 1779
f1538a5f 1780 if(fFillAllTH3)
1781 {
649b825d 1782 fhEtaPhiE = new TH3F ("hEtaPhiE","#eta vs #phi vs energy, reconstructed clusters",
1783 netabins,etamin,etamax,nphibins,phimin,phimax,nptbins,ptmin,ptmax);
1784 fhEtaPhiE->SetXTitle("#eta ");
1785 fhEtaPhiE->SetYTitle("#phi (rad)");
78451bcd 1786 fhEtaPhiE->SetZTitle("#it{E} (GeV) ");
649b825d 1787 outputContainer->Add(fhEtaPhiE);
1788 }
1789
1790 fhClusterTimeEnergy = new TH2F ("hClusterTimeEnergy","energy vs TOF, reconstructed clusters",
1791 nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
78451bcd 1792 fhClusterTimeEnergy->SetXTitle("#it{E} (GeV) ");
649b825d 1793 fhClusterTimeEnergy->SetYTitle("TOF (ns)");
1794 outputContainer->Add(fhClusterTimeEnergy);
1795
1796 fhClusterPairDiffTimeE = new TH2F("hClusterPairDiffTimeE","cluster pair time difference vs E, only good clusters",
1797 nptbins,ptmin,ptmax, tdbins,tdmin,tdmax);
78451bcd 1798 fhClusterPairDiffTimeE->SetXTitle("#it{E}_{cluster} (GeV)");
1799 fhClusterPairDiffTimeE->SetYTitle("#Delta #it{t} (ns)");
649b825d 1800 outputContainer->Add(fhClusterPairDiffTimeE);
1801
1802 fhLambda0 = new TH2F ("hLambda0","shower shape, #lambda^{2}_{0} vs E",
1803 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
78451bcd 1804 fhLambda0->SetXTitle("#it{E}_{cluster}");
649b825d 1805 fhLambda0->SetYTitle("#lambda^{2}_{0}");
1806 outputContainer->Add(fhLambda0);
1807
1808 fhLambda1 = new TH2F ("hLambda1","shower shape, #lambda^{2}_{1} vs E for bad cluster ",
1809 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
78451bcd 1810 fhLambda1->SetXTitle("#it{E}_{cluster}");
649b825d 1811 fhLambda1->SetYTitle("#lambda^{2}_{1}");
1812 outputContainer->Add(fhLambda1);
1813
1814 fhDispersion = new TH2F ("hDispersion","shower shape, Dispersion^{2} vs E for bad cluster ",
1815 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
78451bcd 1816 fhDispersion->SetXTitle("#it{E}_{cluster}");
649b825d 1817 fhDispersion->SetYTitle("Dispersion");
1818 outputContainer->Add(fhDispersion);
1819
1820 fhClusterMaxCellCloseCellRatio = new TH2F ("hClusterMaxCellCloseCellRatio","energy vs ratio of max cell / neighbour cell, reconstructed clusters",
1821 nptbins,ptmin,ptmax, 100,0,1.);
78451bcd 1822 fhClusterMaxCellCloseCellRatio->SetXTitle("#it{E}_{cluster} (GeV) ");
1823 fhClusterMaxCellCloseCellRatio->SetYTitle("#it{E}_{cell i}/#it{E}_{cell max}");
649b825d 1824 outputContainer->Add(fhClusterMaxCellCloseCellRatio);
1825
1826 fhClusterMaxCellCloseCellDiff = new TH2F ("hClusterMaxCellCloseCellDiff","energy vs ratio of max cell / neighbour cell, reconstructed clusters",
1827 nptbins,ptmin,ptmax, 500,0,100.);
78451bcd 1828 fhClusterMaxCellCloseCellDiff->SetXTitle("#it{E}_{cluster} (GeV) ");
1829 fhClusterMaxCellCloseCellDiff->SetYTitle("#it{E}_{cell max}-#it{E}_{cell i} (GeV)");
649b825d 1830 outputContainer->Add(fhClusterMaxCellCloseCellDiff);
1831
1832 fhClusterMaxCellDiff = new TH2F ("hClusterMaxCellDiff","energy vs difference of cluster energy - max cell energy / cluster energy, good clusters",
1833 nptbins,ptmin,ptmax, 500,0,1.);
78451bcd 1834 fhClusterMaxCellDiff->SetXTitle("#it{E}_{cluster} (GeV) ");
1835 fhClusterMaxCellDiff->SetYTitle("(#it{E}_{cluster} - #it{E}_{cell max})/ #it{E}_{cluster}");
649b825d 1836 outputContainer->Add(fhClusterMaxCellDiff);
1837
1838 fhClusterMaxCellDiffNoCut = new TH2F ("hClusterMaxCellDiffNoCut","energy vs difference of cluster energy - max cell energy / cluster energy",
1839 nptbins,ptmin,ptmax, 500,0,1.);
78451bcd 1840 fhClusterMaxCellDiffNoCut->SetXTitle("#it{E}_{cluster} (GeV) ");
1841 fhClusterMaxCellDiffNoCut->SetYTitle("(#it{E}_{cluster} - #it{E}_{cell max})/ #it{E}_{cluster}");
649b825d 1842 outputContainer->Add(fhClusterMaxCellDiffNoCut);
1843
1a72f6c5 1844 fhClusterMaxCellECross = new TH2F ("hClusterMaxCellECross","1 - Energy in cross around max energy cell / max energy cell vs cluster energy, good clusters",
1845 nptbins,ptmin,ptmax, 400,-1,1.);
78451bcd 1846 fhClusterMaxCellECross->SetXTitle("#it{E}_{cluster} (GeV) ");
1847 fhClusterMaxCellECross->SetYTitle("1- #it{E}_{cross}/#it{E}_{cell max}");
1a72f6c5 1848 outputContainer->Add(fhClusterMaxCellECross);
1849
f1538a5f 1850 fhNCellsPerClusterNoCut = new TH2F ("hNCellsPerClusterNoCut","# cells per cluster vs energy, no bad clusters cut",
1851 nptbins,ptmin,ptmax, nceclbins,nceclmin,nceclmax);
78451bcd 1852 fhNCellsPerClusterNoCut->SetXTitle("#it{E} (GeV)");
1853 fhNCellsPerClusterNoCut->SetYTitle("#it{n}_{cells}");
f1538a5f 1854 outputContainer->Add(fhNCellsPerClusterNoCut);
1855
1856 fhNCellsPerCluster = new TH2F ("hNCellsPerCluster","# cells per cluster vs energy",nptbins,ptmin,ptmax, nceclbins,nceclmin,nceclmax);
78451bcd 1857 fhNCellsPerCluster->SetXTitle("#it{E} (GeV)");
1858 fhNCellsPerCluster->SetYTitle("#it{n}_{cells}");
f1538a5f 1859 outputContainer->Add(fhNCellsPerCluster);
1860
1861 fhNClusters = new TH1F ("hNClusters","# clusters", nclbins,nclmin,nclmax);
78451bcd 1862 fhNClusters->SetXTitle("#it{n}_{clusters}");
f1538a5f 1863 outputContainer->Add(fhNClusters);
1864
1a83b960 1865 if(fStudyBadClusters)
1866 {
649b825d 1867 fhBadClusterEnergy = new TH1F ("hBadClusterEnergy","Bad cluster energy", nptbins,ptmin,ptmax);
78451bcd 1868 fhBadClusterEnergy->SetXTitle("#it{E}_{cluster} (GeV) ");
649b825d 1869 outputContainer->Add(fhBadClusterEnergy);
1870
1871 fhBadClusterMaxCellCloseCellRatio = new TH2F ("hBadClusterMaxCellCloseCellRatio","energy vs ratio of max cell / neighbour cell constributing cell, reconstructed bad clusters",
1872 nptbins,ptmin,ptmax, 100,0,1.);
78451bcd 1873 fhBadClusterMaxCellCloseCellRatio->SetXTitle("#it{E}_{cluster} (GeV) ");
649b825d 1874 fhBadClusterMaxCellCloseCellRatio->SetYTitle("ratio");
1875 outputContainer->Add(fhBadClusterMaxCellCloseCellRatio);
1876
1877 fhBadClusterMaxCellCloseCellDiff = new TH2F ("hBadClusterMaxCellCloseCellDiff","energy vs ratio of max cell - neighbour cell constributing cell, reconstructed bad clusters",
1878 nptbins,ptmin,ptmax, 500,0,100);
78451bcd 1879 fhBadClusterMaxCellCloseCellDiff->SetXTitle("#it{E}_{cluster} (GeV) ");
1880 fhBadClusterMaxCellCloseCellDiff->SetYTitle("#it{E}_{cell max} - #it{E}_{cell i} (GeV)");
649b825d 1881 outputContainer->Add(fhBadClusterMaxCellCloseCellDiff);
1882
1883 fhBadClusterMaxCellDiff = new TH2F ("hBadClusterMaxCellDiff","energy vs difference of cluster energy - max cell energy / cluster energy for bad clusters",
1884 nptbins,ptmin,ptmax, 500,0,1.);
78451bcd 1885 fhBadClusterMaxCellDiff->SetXTitle("#it{E}_{cluster} (GeV) ");
1886 fhBadClusterMaxCellDiff->SetYTitle("(#it{E}_{cluster} - #it{E}_{cell max}) / #it{E}_{cluster}");
649b825d 1887 outputContainer->Add(fhBadClusterMaxCellDiff);
1888
1889 fhBadClusterTimeEnergy = new TH2F ("hBadClusterTimeEnergy","energy vs TOF of reconstructed bad clusters",
1890 nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
78451bcd 1891 fhBadClusterTimeEnergy->SetXTitle("#it{E}_{cluster} (GeV) ");
1892 fhBadClusterTimeEnergy->SetYTitle("#it{t} (ns)");
649b825d 1893 outputContainer->Add(fhBadClusterTimeEnergy);
1894
1895 fhBadClusterPairDiffTimeE = new TH2F("hBadClusterPairDiffTimeE","cluster pair time difference (bad - good) vs E from bad cluster",nptbins,ptmin,ptmax, tdbins,tdmin,tdmax);
78451bcd 1896 fhBadClusterPairDiffTimeE->SetXTitle("#it{E}_{bad cluster} (GeV)");
1897 fhBadClusterPairDiffTimeE->SetYTitle("#Delta #it{t} (ns)");
649b825d 1898 outputContainer->Add(fhBadClusterPairDiffTimeE);
1899
78451bcd 1900 fhBadClusterMaxCellECross = new TH2F ("hBadClusterMaxCellECross","1 - #it{E}_{+} around max energy cell / max energy cell vs cluster energy, bad clusters",
1a72f6c5 1901 nptbins,ptmin,ptmax, 400,-1,1.);
78451bcd 1902 fhBadClusterMaxCellECross->SetXTitle("#it{E}_{cluster} (GeV) ");
1903 fhBadClusterMaxCellECross->SetYTitle("1- #it{E}_{cross}/#it{E}_{cell max}");
1a72f6c5 1904 outputContainer->Add(fhBadClusterMaxCellECross);
1905
e6fec6f5 1906 if(fFillAllCellTimeHisto)
1907 {
78451bcd 1908 fhBadCellTimeSpreadRespectToCellMax = new TH2F ("hBadCellTimeSpreadRespectToCellMax","#it{t}_{cell max}-#it{t}_{cell i} from bad cluster", nptbins,ptmin,ptmax, tdbins,tdmin,tdmax);
1909 fhBadCellTimeSpreadRespectToCellMax->SetXTitle("#it{E} (GeV)");
1910 fhBadCellTimeSpreadRespectToCellMax->SetYTitle("#Delta #it{t}_{cell max - i} (ns)");
649b825d 1911 outputContainer->Add(fhBadCellTimeSpreadRespectToCellMax);
1912
78451bcd 1913 fhBadClusterMaxCellDiffAverageTime = new TH2F ("hBadClusterMaxCellDiffAverageTime","#it{t}_{cell max}-#it{t}_{average} from bad cluster", nptbins,ptmin,ptmax, tdbins,tdmin,tdmax);
1914 fhBadClusterMaxCellDiffAverageTime->SetXTitle("#it{E} (GeV)");
1915 fhBadClusterMaxCellDiffAverageTime->SetYTitle("#Delta #it{t}_{cell max - average} (ns)");
649b825d 1916 outputContainer->Add(fhBadClusterMaxCellDiffAverageTime);
a82b4462 1917
78451bcd 1918 fhBadClusterMaxCellDiffWeightedTime = new TH2F ("hBadClusterMaxCellDiffWeightedTime","#it{t}_{cell max}-#it{t}_{weighted} from bad cluster", nptbins,ptmin,ptmax, tdbins,tdmin,tdmax);
1919 fhBadClusterMaxCellDiffWeightedTime->SetXTitle("#it{E} (GeV)");
1920 fhBadClusterMaxCellDiffWeightedTime->SetYTitle("#Delta #it{t}_{cell max - weighted} (ns)");
649b825d 1921 outputContainer->Add(fhBadClusterMaxCellDiffWeightedTime);
1922
649b825d 1923 }
1924
1925 }
1926
f1538a5f 1927 if(fStudyExotic)
1928 {
765206a5 1929 fhExoL0ECross = new TH2F("hExoL0_ECross",
78451bcd 1930 "#lambda^{2}_{0} vs 1-#it{E}_{+}/#it{E}_{max} for E > 5 GeV",
e49010f3 1931 400,0,1,ssbins,ssmin,ssmax);
78451bcd 1932 fhExoL0ECross ->SetXTitle("1-#it{E}_{+}/#it{E}_{cell max}");
cdacd584 1933 fhExoL0ECross ->SetYTitle("#lambda^{2}_{0}");
765206a5 1934 outputContainer->Add(fhExoL0ECross) ;
1935
1936 fhExoL1ECross = new TH2F("hExoL1_ECross",
78451bcd 1937 "#lambda^{2}_{1} vs 1-#it{E}_{+}/#it{E}_{max} for E > 5 GeV",
e49010f3 1938 400,0,1,ssbins,ssmin,ssmax);
78451bcd 1939 fhExoL1ECross ->SetXTitle("1-#it{E}_{+}/#it{E}_{cell max}");
cdacd584 1940 fhExoL1ECross ->SetYTitle("#lambda^{2}_{1}");
765206a5 1941 outputContainer->Add(fhExoL1ECross) ;
1942
f1538a5f 1943 for(Int_t ie = 0; ie <fExoNECrossCuts; ie++)
1944 {
1945
1946 fhExoDTime[ie] = new TH2F(Form("hExoDTime_ECross%d",ie),
78451bcd 1947 Form("#Delta time = t_{max}-t_{cells} vs #it{E}_{cluster} for exotic, 1-#it{E}_{+}/#it{E}_{max} < %2.2f",fExoECrossCuts[ie]),
f1538a5f 1948 nptbins,ptmin,ptmax,tdbins,tdmin,tdmax);
78451bcd 1949 fhExoDTime[ie] ->SetYTitle("#Delta #it{t} (ns)");
1950 fhExoDTime[ie] ->SetXTitle("#it{E} (GeV)");
f1538a5f 1951 outputContainer->Add(fhExoDTime[ie]) ;
1952
1953 for(Int_t idt = 0; idt < fExoNDTimeCuts; idt++)
1954 {
1955 fhExoNCell[ie][idt] = new TH2F(Form("hExoNCell_ECross%d_DT%d",ie,idt),
78451bcd 1956 Form("N cells per cluster vs E cluster, 1-#it{E}_{+}/#it{E}_{max} < %2.2f, #Delta t < %2.0f",fExoECrossCuts[ie],fExoDTimeCuts[idt]),
f1538a5f 1957 nptbins,ptmin,ptmax,nceclbins,nceclmin,nceclmax);
78451bcd 1958 fhExoNCell[ie][idt] ->SetYTitle("#it{n}_cells");
1959 fhExoNCell[ie][idt] ->SetXTitle("#it{E} (GeV)");
f1538a5f 1960 outputContainer->Add(fhExoNCell[ie][idt]) ;
1961
1962 fhExoL0 [ie][idt] = new TH2F(Form("hExoL0_ECross%d_DT%d",ie,idt),
78451bcd 1963 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]),
f1538a5f 1964 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
1965 fhExoL0 [ie][idt] ->SetYTitle("#lambda^{2}_{0}");
78451bcd 1966 fhExoL0 [ie][idt] ->SetXTitle("#it{E} (GeV)");
f1538a5f 1967 outputContainer->Add(fhExoL0[ie][idt]) ;
765206a5 1968
1969 fhExoL1 [ie][idt] = new TH2F(Form("hExoL1_ECross%d_DT%d",ie,idt),
78451bcd 1970 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]),
765206a5 1971 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
1972 fhExoL1 [ie][idt] ->SetYTitle("#lambda^{2}_{1}");
78451bcd 1973 fhExoL1 [ie][idt] ->SetXTitle("#it{E} (GeV)");
765206a5 1974 outputContainer->Add(fhExoL1[ie][idt]) ;
f1538a5f 1975
1976 fhExoECross[ie][idt] = new TH2F(Form("hExoECross_ECross%d_DT%d",ie,idt),
78451bcd 1977 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]),
f1538a5f 1978 nptbins,ptmin,ptmax,400,0,1);
78451bcd 1979 fhExoECross[ie][idt] ->SetYTitle("1-#it{E}_{+}/#it{E}_{cell max}");
1980 fhExoECross[ie][idt] ->SetXTitle("#it{E}_{cell} (GeV)");
f1538a5f 1981 outputContainer->Add(fhExoECross[ie][idt]) ;
1982
1983 fhExoTime [ie][idt] = new TH2F(Form("hExoTime_ECross%d_DT%d",ie,idt),
78451bcd 1984 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]),
f1538a5f 1985 nptbins,ptmin,ptmax,ntimebins,timemin,timemax);
78451bcd 1986 fhExoTime [ie][idt] ->SetYTitle("#it{t}_{max} (ns)");
1987 fhExoTime [ie][idt] ->SetXTitle("#it{E} (GeV)");
f1538a5f 1988 outputContainer->Add(fhExoTime[ie][idt]) ;
1989
765206a5 1990 fhExoL0NCell[ie][idt] = new TH2F(Form("hExoL0_NCell%d_DT%d",ie,idt),
78451bcd 1991 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]),
765206a5 1992 nptbins,ptmin,ptmax,ntimebins,timemin,timemax);
78451bcd 1993 fhExoL0NCell[ie][idt] ->SetYTitle("#it{n}_{cells}");
765206a5 1994 fhExoL0NCell[ie][idt] ->SetXTitle("#lambda^{2}_{0}");
1995 outputContainer->Add(fhExoL0NCell[ie][idt]) ;
1996
1997 fhExoL1NCell[ie][idt] = new TH2F(Form("hExoL1_NCell%d_DT%d",ie,idt),
78451bcd 1998 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]),
765206a5 1999 nptbins,ptmin,ptmax,ntimebins,timemin,timemax);
78451bcd 2000 fhExoL1NCell[ie][idt] ->SetYTitle("#it{n}_{cells}");
765206a5 2001 fhExoL1NCell[ie][idt] ->SetXTitle("#lambda^{2}_{1}");
2002 outputContainer->Add(fhExoL1NCell[ie][idt]) ;
2003
f1538a5f 2004 }
2005 }
2006 }
2007
649b825d 2008 // Cluster size in terms of cells
f1538a5f 2009 if(fStudyClustersAsymmetry)
2010 {
78451bcd 2011 fhDeltaIEtaDeltaIPhiE0[0] = new TH2F ("hDeltaIEtaDeltaIPhiE0"," Cluster size in columns vs rows for E < 2 GeV, #it{n}_{cells} > 3",
649b825d 2012 50,0,50,50,0,50);
2013 fhDeltaIEtaDeltaIPhiE0[0]->SetXTitle("#Delta Column");
2014 fhDeltaIEtaDeltaIPhiE0[0]->SetYTitle("#Delta Row");
2015 outputContainer->Add(fhDeltaIEtaDeltaIPhiE0[0]);
2016
78451bcd 2017 fhDeltaIEtaDeltaIPhiE2[0] = new TH2F ("hDeltaIEtaDeltaIPhiE2"," Cluster size in columns vs rows for 2 <E < 6 GeV, #it{n}_{cells} > 3",
649b825d 2018 50,0,50,50,0,50);
2019 fhDeltaIEtaDeltaIPhiE2[0]->SetXTitle("#Delta Column");
2020 fhDeltaIEtaDeltaIPhiE2[0]->SetYTitle("#Delta Row");
2021 outputContainer->Add(fhDeltaIEtaDeltaIPhiE2[0]);
2022
78451bcd 2023 fhDeltaIEtaDeltaIPhiE6[0] = new TH2F ("hDeltaIEtaDeltaIPhiE6"," Cluster size in columns vs rows for E > 6 GeV, #it{n}_{cells} > 3",
649b825d 2024 50,0,50,50,0,50);
2025 fhDeltaIEtaDeltaIPhiE6[0]->SetXTitle("#Delta Column");
2026 fhDeltaIEtaDeltaIPhiE6[0]->SetYTitle("#Delta Row");
2027 outputContainer->Add(fhDeltaIEtaDeltaIPhiE6[0]);
2028
2029 fhDeltaIA[0] = new TH2F ("hDeltaIA"," Cluster *asymmetry* in cell units vs E",
2030 nptbins,ptmin,ptmax,21,-1.05,1.05);
78451bcd 2031 fhDeltaIA[0]->SetXTitle("#it{E}_{cluster}");
2032 fhDeltaIA[0]->SetYTitle("#it{A}_{cell in cluster}");
649b825d 2033 outputContainer->Add(fhDeltaIA[0]);
2034
2035 fhDeltaIAL0[0] = new TH2F ("hDeltaIAL0"," Cluster *asymmetry* in cell units vs #lambda^{2}_{0}",
2036 ssbins,ssmin,ssmax,21,-1.05,1.05);
2037 fhDeltaIAL0[0]->SetXTitle("#lambda^{2}_{0}");
78451bcd 2038 fhDeltaIAL0[0]->SetYTitle("#it{A}_{cell in cluster}");
649b825d 2039 outputContainer->Add(fhDeltaIAL0[0]);
2040
2041 fhDeltaIAL1[0] = new TH2F ("hDeltaIAL1"," Cluster *asymmetry* in cell units vs #lambda^{2}_{1}",
2042 ssbins,ssmin,ssmax,21,-1.05,1.05);
2043 fhDeltaIAL1[0]->SetXTitle("#lambda^{2}_{1}");
78451bcd 2044 fhDeltaIAL1[0]->SetYTitle("#it{A}_{cell in cluster}");
649b825d 2045 outputContainer->Add(fhDeltaIAL1[0]);
2046
2047 fhDeltaIANCells[0] = new TH2F ("hDeltaIANCells"," Cluster *asymmetry* in cell units vs N cells in cluster",
2048 nceclbins,nceclmin,nceclmax,21,-1.05,1.05);
78451bcd 2049 fhDeltaIANCells[0]->SetXTitle("#it{n}_{cell in cluster}");
2050 fhDeltaIANCells[0]->SetYTitle("#it{A}_{cell in cluster}");
649b825d 2051 outputContainer->Add(fhDeltaIANCells[0]);
2052
2053
78451bcd 2054 fhDeltaIEtaDeltaIPhiE0[1] = new TH2F ("hDeltaIEtaDeltaIPhiE0Charged"," Cluster size in columns vs rows for E < 2 GeV, #it{n}_{cells} > 3, matched with track",
649b825d 2055 50,0,50,50,0,50);
2056 fhDeltaIEtaDeltaIPhiE0[1]->SetXTitle("#Delta Column");
2057 fhDeltaIEtaDeltaIPhiE0[1]->SetYTitle("#Delta Row");
2058 outputContainer->Add(fhDeltaIEtaDeltaIPhiE0[1]);
2059
78451bcd 2060 fhDeltaIEtaDeltaIPhiE2[1] = new TH2F ("hDeltaIEtaDeltaIPhiE2Charged"," Cluster size in columns vs rows for 2 <E < 6 GeV, #it{n}_{cells} > 3, matched with track",
649b825d 2061 50,0,50,50,0,50);
2062 fhDeltaIEtaDeltaIPhiE2[1]->SetXTitle("#Delta Column");
2063 fhDeltaIEtaDeltaIPhiE2[1]->SetYTitle("#Delta Row");
2064 outputContainer->Add(fhDeltaIEtaDeltaIPhiE2[1]);
2065
78451bcd 2066 fhDeltaIEtaDeltaIPhiE6[1] = new TH2F ("hDeltaIEtaDeltaIPhiE6Charged"," Cluster size in columns vs rows for E > 6 GeV, #it{n}_{cells} > 3, matched with track",
649b825d 2067 50,0,50,50,0,50);
2068 fhDeltaIEtaDeltaIPhiE6[1]->SetXTitle("#Delta Column");
2069 fhDeltaIEtaDeltaIPhiE6[1]->SetYTitle("#Delta Row");
2070 outputContainer->Add(fhDeltaIEtaDeltaIPhiE6[1]);
2071
2072 fhDeltaIA[1] = new TH2F ("hDeltaIACharged"," Cluster *asymmetry* in cell units vs E, matched with track",
2073 nptbins,ptmin,ptmax,21,-1.05,1.05);
78451bcd 2074 fhDeltaIA[1]->SetXTitle("#it{E}_{cluster}");
2075 fhDeltaIA[1]->SetYTitle("#it{A}_{cell in cluster}");
649b825d 2076 outputContainer->Add(fhDeltaIA[1]);
2077
2078 fhDeltaIAL0[1] = new TH2F ("hDeltaIAL0Charged"," Cluster *asymmetry* in cell units vs #lambda^{2}_{0}, matched with track",
2079 ssbins,ssmin,ssmax,21,-1.05,1.05);
2080 fhDeltaIAL0[1]->SetXTitle("#lambda^{2}_{0}");
78451bcd 2081 fhDeltaIAL0[1]->SetYTitle("#it{A}_{cell in cluster}");
649b825d 2082 outputContainer->Add(fhDeltaIAL0[1]);
2083
2084 fhDeltaIAL1[1] = new TH2F ("hDeltaIAL1Charged"," Cluster *asymmetry* in cell units vs #lambda^{2}_{1}, matched with track",
2085 ssbins,ssmin,ssmax,21,-1.05,1.05);
2086 fhDeltaIAL1[1]->SetXTitle("#lambda^{2}_{1}");
78451bcd 2087 fhDeltaIAL1[1]->SetYTitle("#it{A}_{cell in cluster}");
649b825d 2088 outputContainer->Add(fhDeltaIAL1[1]);
2089
2090 fhDeltaIANCells[1] = new TH2F ("hDeltaIANCellsCharged"," Cluster *asymmetry* in cell units vs N cells in cluster, matched with track",
2091 nceclbins,nceclmin,nceclmax,21,-1.05,1.05);
78451bcd 2092 fhDeltaIANCells[1]->SetXTitle("#it{n}_{cell in cluster}");
2093 fhDeltaIANCells[1]->SetYTitle("#it{A}_{cell in cluster}");
649b825d 2094 outputContainer->Add(fhDeltaIANCells[1]);
2095
2096 if(IsDataMC()){
2097 TString particle[]={"Photon","Electron","Conversion","Hadron"};
2098 for (Int_t iPart = 0; iPart < 4; iPart++) {
2099
2100 fhDeltaIAMC[iPart] = new TH2F (Form("hDeltaIA_MC%s",particle[iPart].Data()),Form(" Cluster *asymmetry* in cell units vs E, from %s",particle[iPart].Data()),
2101 nptbins,ptmin,ptmax,21,-1.05,1.05);
78451bcd 2102 fhDeltaIAMC[iPart]->SetXTitle("#it{E}_{cluster}");
2103 fhDeltaIAMC[iPart]->SetYTitle("#it{A}_{cell in cluster}");
649b825d 2104 outputContainer->Add(fhDeltaIAMC[iPart]);
2105 }
2106 }
f1538a5f 2107
1a83b960 2108 if(fStudyBadClusters)
2109 {
78451bcd 2110 fhBadClusterDeltaIEtaDeltaIPhiE0 = new TH2F ("hBadClusterDeltaIEtaDeltaIPhiE0"," Cluster size in columns vs rows for E < 2 GeV, #it{n}_{cells} > 3",
1a83b960 2111 50,0,50,50,0,50);
2112 fhBadClusterDeltaIEtaDeltaIPhiE0->SetXTitle("#Delta Column");
2113 fhBadClusterDeltaIEtaDeltaIPhiE0->SetYTitle("#Delta Row");
2114 outputContainer->Add(fhBadClusterDeltaIEtaDeltaIPhiE0);
2115
78451bcd 2116 fhBadClusterDeltaIEtaDeltaIPhiE2 = new TH2F ("hBadClusterDeltaIEtaDeltaIPhiE2"," Cluster size in columns vs rows for 2 <E < 6 GeV, #it{n}_{cells} > 3",
1a83b960 2117 50,0,50,50,0,50);
2118 fhBadClusterDeltaIEtaDeltaIPhiE2->SetXTitle("#Delta Column");
2119 fhBadClusterDeltaIEtaDeltaIPhiE2->SetYTitle("#Delta Row");
2120 outputContainer->Add(fhBadClusterDeltaIEtaDeltaIPhiE2);
2121
78451bcd 2122 fhBadClusterDeltaIEtaDeltaIPhiE6 = new TH2F ("hBadClusterDeltaIEtaDeltaIPhiE6"," Cluster size in columns vs rows for E > 6 GeV, #it{n}_{cells} > 3",
1a83b960 2123 50,0,50,50,0,50);
2124 fhBadClusterDeltaIEtaDeltaIPhiE6->SetXTitle("#Delta Column");
2125 fhBadClusterDeltaIEtaDeltaIPhiE6->SetYTitle("#Delta Row");
2126 outputContainer->Add(fhBadClusterDeltaIEtaDeltaIPhiE6);
2127
2128 fhBadClusterDeltaIA = new TH2F ("hBadClusterDeltaIA"," Cluster *asymmetry* in cell units vs E",
2129 nptbins,ptmin,ptmax,21,-1.05,1.05);
78451bcd 2130 fhBadClusterDeltaIA->SetXTitle("#it{E}_{cluster}");
2131 fhBadClusterDeltaIA->SetYTitle("#it{A}_{cell in cluster}");
1a83b960 2132 outputContainer->Add(fhBadClusterDeltaIA);
2133 }
649b825d 2134 }
2135
f1538a5f 2136 if(fStudyWeight)
2137 {
649b825d 2138 fhECellClusterRatio = new TH2F ("hECellClusterRatio"," cell energy / cluster energy vs cluster energy",
2139 nptbins,ptmin,ptmax, 100,0,1.);
78451bcd 2140 fhECellClusterRatio->SetXTitle("#it{E}_{cluster} (GeV) ");
2141 fhECellClusterRatio->SetYTitle("#it{E}_{cell i}/#it{E}_{cluster}");
649b825d 2142 outputContainer->Add(fhECellClusterRatio);
2143
2144 fhECellClusterLogRatio = new TH2F ("hECellClusterLogRatio"," Log(cell energy / cluster energy) vs cluster energy",
1a72f6c5 2145 nptbins,ptmin,ptmax, 100,-10,0);
78451bcd 2146 fhECellClusterLogRatio->SetXTitle("#it{E}_{cluster} (GeV) ");
2147 fhECellClusterLogRatio->SetYTitle("Log(#it{E}_{cell i}/#it{E}_{cluster})");
649b825d 2148 outputContainer->Add(fhECellClusterLogRatio);
2149
2150 fhEMaxCellClusterRatio = new TH2F ("hEMaxCellClusterRatio"," max cell energy / cluster energy vs cluster energy",
2151 nptbins,ptmin,ptmax, 100,0,1.);
78451bcd 2152 fhEMaxCellClusterRatio->SetXTitle("#it{E}_{cluster} (GeV) ");
2153 fhEMaxCellClusterRatio->SetYTitle("#it{E}_{max cell}/#it{E}_{cluster}");
649b825d 2154 outputContainer->Add(fhEMaxCellClusterRatio);
2155
2156 fhEMaxCellClusterLogRatio = new TH2F ("hEMaxCellClusterLogRatio"," Log(max cell energy / cluster energy) vs cluster energy",
1a72f6c5 2157 nptbins,ptmin,ptmax, 100,-10,0);
78451bcd 2158 fhEMaxCellClusterLogRatio->SetXTitle("#it{E}_{cluster} (GeV) ");
2159 fhEMaxCellClusterLogRatio->SetYTitle("Log (#it{E}_{max cell}/#it{E}_{cluster})");
649b825d 2160 outputContainer->Add(fhEMaxCellClusterLogRatio);
2161
701cbf54 2162 fhECellTotalRatio = new TH2F ("hECellTotalRatio"," cell energy / sum all energy vs all energy",
2163 nptbins*2,ptmin,ptmax*2, 100,0,1.);
78451bcd 2164 fhECellTotalRatio->SetXTitle("#it{E}_{total} (GeV) ");
2165 fhECellTotalRatio->SetYTitle("#it{E}_{cell i}/#it{E}_{total}");
701cbf54 2166 outputContainer->Add(fhECellTotalRatio);
2167
2168 fhECellTotalLogRatio = new TH2F ("hECellTotalLogRatio"," Log(cell energy / sum all energy) vs all energy",
2169 nptbins*2,ptmin,ptmax*2, 100,-10,0);
78451bcd 2170 fhECellTotalLogRatio->SetXTitle("#it{E}_{total} (GeV) ");
2171 fhECellTotalLogRatio->SetYTitle("Log(#it{E}_{cell i}/#it{E}_{total})");
701cbf54 2172 outputContainer->Add(fhECellTotalLogRatio);
2173
2174 fhECellTotalRatioMod = new TH2F*[fNModules];
2175 fhECellTotalLogRatioMod = new TH2F*[fNModules];
2176
2177 for(Int_t imod = 0; imod < fNModules; imod++)
2178 {
2179 fhECellTotalRatioMod[imod] = new TH2F (Form("hECellTotalRatio_Mod%d",imod),
2180 Form("#cell energy / sum all energy vs all energy in Module %d",imod),
2181 nptbins*2,ptmin,ptmax*2, 100,0,1.);
78451bcd 2182 fhECellTotalRatioMod[imod]->SetXTitle("#it{E} (GeV)");
2183 fhECellTotalRatioMod[imod]->SetYTitle("#it{n}_{cells}");
701cbf54 2184 outputContainer->Add(fhECellTotalRatioMod[imod]);
2185
2186 fhECellTotalLogRatioMod[imod] = new TH2F (Form("hECellTotalLogRatio_Mod%d",imod),
2187 Form("Log(cell energy / sum all energy) vs all energy in Module %d",imod),
2188 nptbins*2,ptmin,ptmax*2, 100,-10,0);
78451bcd 2189 fhECellTotalLogRatioMod[imod]->SetXTitle("#it{E} (GeV)");
2190 fhECellTotalLogRatioMod[imod]->SetYTitle("#it{n}_{cells}");
701cbf54 2191 outputContainer->Add(fhECellTotalLogRatioMod[imod]);
2192
2193 }
2194
2195 for(Int_t iw = 0; iw < 12; iw++)
2196 {
2197 Float_t w0 = 3+0.25*iw;
2198 fhLambda0ForW0[iw] = new TH2F (Form("hLambda0ForW0%d",iw),Form("shower shape, #lambda^{2}_{0} vs E, w0 = %1.1f",w0),
649b825d 2199 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
78451bcd 2200 fhLambda0ForW0[iw]->SetXTitle("#it{E}_{cluster}");
649b825d 2201 fhLambda0ForW0[iw]->SetYTitle("#lambda^{2}_{0}");
2202 outputContainer->Add(fhLambda0ForW0[iw]);
2203
701cbf54 2204// fhLambda1ForW0[iw] = new TH2F (Form("hLambda1ForW0%d",iw),Form("shower shape, #lambda^{2}_{1} vs E, w0 = %1.1f",w0),
1a72f6c5 2205// nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
78451bcd 2206// fhLambda1ForW0[iw]->SetXTitle("#it{E}_{cluster}");
1a72f6c5 2207// fhLambda1ForW0[iw]->SetYTitle("#lambda^{2}_{1}");
2208// outputContainer->Add(fhLambda1ForW0[iw]);
649b825d 2209
2210 if(IsDataMC()){
2211 TString mcnames[] = {"Photon", "Electron","Conversion","Pi0","Hadron"};
2212 for(Int_t imc = 0; imc < 5; imc++){
2213 fhLambda0ForW0MC[iw][imc] = new TH2F (Form("hLambda0ForW0%d_MC%s",iw,mcnames[imc].Data()),
701cbf54 2214 Form("shower shape, #lambda^{2}_{0} vs E, w0 = %1.1f, for MC %s",w0,mcnames[imc].Data()),
649b825d 2215 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
78451bcd 2216 fhLambda0ForW0MC[iw][imc]->SetXTitle("#it{E}_{cluster}");
649b825d 2217 fhLambda0ForW0MC[iw][imc]->SetYTitle("#lambda^{2}_{0}");
2218 outputContainer->Add(fhLambda0ForW0MC[iw][imc]);
2219
1a72f6c5 2220// fhLambda1ForW0MC[iw][imc] = new TH2F (Form("hLambda1ForW0%d_MC%s",iw,mcnames[imc].Data()),
701cbf54 2221// Form("shower shape, #lambda^{2}_{1} vs E, w0 = %1.1f, for MC %s",w0,mcnames[imc].Data()),
1a72f6c5 2222// nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
78451bcd 2223// fhLambda1ForW0MC[iw][imc]->SetXTitle("#it{E}_{cluster}");
1a72f6c5 2224// fhLambda1ForW0MC[iw][imc]->SetYTitle("#lambda^{2}_{1}");
2225// outputContainer->Add(fhLambda1ForW0MC[iw][imc]);
649b825d 2226 }
2227 }
f1538a5f 2228 }
649b825d 2229 }
2230
2231 //Track Matching
f1538a5f 2232 if(fFillAllTMHisto)
2233 {
653aed3c 2234 Int_t nresetabins = GetHistogramRanges()->GetHistoTrackResidualEtaBins();
2235 Float_t resetamax = GetHistogramRanges()->GetHistoTrackResidualEtaMax();
2236 Float_t resetamin = GetHistogramRanges()->GetHistoTrackResidualEtaMin();
2237 Int_t nresphibins = GetHistogramRanges()->GetHistoTrackResidualPhiBins();
2238 Float_t resphimax = GetHistogramRanges()->GetHistoTrackResidualPhiMax();
2239 Float_t resphimin = GetHistogramRanges()->GetHistoTrackResidualPhiMin();
2240
2241 fhTrackMatchedDEta = new TH2F("hTrackMatchedDEta","d#eta of cluster-track vs cluster energy",
2242 nptbins,ptmin,ptmax,nresetabins,resetamin,resetamax);
2243 fhTrackMatchedDEta->SetYTitle("d#eta");
78451bcd 2244 fhTrackMatchedDEta->SetXTitle("#it{E}_{cluster} (GeV)");
653aed3c 2245
2246 fhTrackMatchedDPhi = new TH2F("hTrackMatchedDPhi","d#phi of cluster-track vs cluster energy",
2247 nptbins,ptmin,ptmax,nresphibins,resphimin,resphimax);
2248 fhTrackMatchedDPhi->SetYTitle("d#phi (rad)");
78451bcd 2249 fhTrackMatchedDPhi->SetXTitle("#it{E}_{cluster} (GeV)");
653aed3c 2250
2251 fhTrackMatchedDEtaDPhi = new TH2F("hTrackMatchedDEtaDPhi","d#eta vs d#phi of cluster-track vs cluster energy",
2252 nresetabins,resetamin,resetamax,nresphibins,resphimin,resphimax);
2253 fhTrackMatchedDEtaDPhi->SetYTitle("d#phi (rad)");
2254 fhTrackMatchedDEtaDPhi->SetXTitle("d#eta");
2255
2256 fhTrackMatchedDEtaPos = new TH2F("hTrackMatchedDEtaPos","d#eta of cluster-track vs cluster energy",
2257 nptbins,ptmin,ptmax,nresetabins,resetamin,resetamax);
2258 fhTrackMatchedDEtaPos->SetYTitle("d#eta");
78451bcd 2259 fhTrackMatchedDEtaPos->SetXTitle("#it{E}_{cluster} (GeV)");
653aed3c 2260
2261 fhTrackMatchedDPhiPos = new TH2F("hTrackMatchedDPhiPos","d#phi of cluster-track vs cluster energy",
2262 nptbins,ptmin,ptmax,nresphibins,resphimin,resphimax);
2263 fhTrackMatchedDPhiPos->SetYTitle("d#phi (rad)");
78451bcd 2264 fhTrackMatchedDPhiPos->SetXTitle("#it{E}_{cluster} (GeV)");
653aed3c 2265
2266 fhTrackMatchedDEtaDPhiPos = new TH2F("hTrackMatchedDEtaDPhiPos","d#eta vs d#phi of cluster-track vs cluster energy",
2267 nresetabins,resetamin,resetamax,nresphibins,resphimin,resphimax);
2268 fhTrackMatchedDEtaDPhiPos->SetYTitle("d#phi (rad)");
2269 fhTrackMatchedDEtaDPhiPos->SetXTitle("d#eta");
2270
2271 outputContainer->Add(fhTrackMatchedDEta) ;
2272 outputContainer->Add(fhTrackMatchedDPhi) ;
2273 outputContainer->Add(fhTrackMatchedDEtaDPhi) ;
2274 outputContainer->Add(fhTrackMatchedDEtaPos) ;
2275 outputContainer->Add(fhTrackMatchedDPhiPos) ;
2276 outputContainer->Add(fhTrackMatchedDEtaDPhiPos) ;
2277
78451bcd 2278 fhECharged = new TH1F ("hECharged","#it{E} reconstructed clusters, matched with track", nptbins,ptmin,ptmax);
2279 fhECharged->SetXTitle("#it{E} (GeV)");
45769d5b 2280 outputContainer->Add(fhECharged);
2281
78451bcd 2282 fhPtCharged = new TH1F ("hPtCharged","#it{p}_{T} reconstructed clusters, matched with track", nptbins,ptmin,ptmax);
2283 fhPtCharged->SetXTitle("#it{p}_{T} (GeV/#it{c})");
45769d5b 2284 outputContainer->Add(fhPtCharged);
2285
2286 fhPhiCharged = new TH1F ("hPhiCharged","#phi reconstructed clusters, matched with track",nphibins,phimin,phimax);
2287 fhPhiCharged->SetXTitle("#phi (rad)");
2288 outputContainer->Add(fhPhiCharged);
2289
2290 fhEtaCharged = new TH1F ("hEtaCharged","#eta reconstructed clusters, matched with track",netabins,etamin,etamax);
2291 fhEtaCharged->SetXTitle("#eta ");
2292 outputContainer->Add(fhEtaCharged);
2293
2294 if(fFillAllTH3)
f1538a5f 2295 {
3748ffb5 2296 fhEtaPhiECharged = new TH3F ("hEtaPhiECharged","#eta vs #phi, reconstructed clusters, matched with track",
521636d2 2297 netabins,etamin,etamax,nphibins,phimin,phimax,nptbins,ptmin,ptmax);
3748ffb5 2298 fhEtaPhiECharged->SetXTitle("#eta ");
2299 fhEtaPhiECharged->SetYTitle("#phi ");
78451bcd 2300 fhEtaPhiECharged->SetZTitle("#it{E} (GeV) ");
3748ffb5 2301 outputContainer->Add(fhEtaPhiECharged);
2302 }
55c05f8c 2303
78451bcd 2304 fh1EOverP = new TH2F("h1EOverP","TRACK matches #it{E}/#it{p}",nptbins,ptmin,ptmax, nPoverEbins,eOverPmin,eOverPmax);
2305 fh1EOverP->SetYTitle("#it{E}/#it{p}");
2306 fh1EOverP->SetXTitle("#it{p}_{T} (GeV/#it{c})");
d55bb5e1 2307 outputContainer->Add(fh1EOverP);
3bfc4732 2308
78451bcd 2309 fh2dR = new TH2F("h2dR","TRACK matches #Delta #it{R}",nptbins,ptmin,ptmax,ndRbins,dRmin,dRmax);
2310 fh2dR->SetXTitle("#Delta #it{R} (rad)");
2311 fh2dR->SetXTitle("#it{E} cluster (GeV)");
a87e069d 2312 outputContainer->Add(fh2dR) ;
3bfc4732 2313
78451bcd 2314 fh2MatchdEdx = new TH2F("h2MatchdEdx","#it{dE/dx} vs. #it{p} for all matches",nptbins,ptmin,ptmax,ndedxbins,dedxmin,dedxmax);
2315 fh2MatchdEdx->SetXTitle("p (GeV/#it{c})");
2316 fh2MatchdEdx->SetYTitle("#it{dE/dx}>");
3bfc4732 2317 outputContainer->Add(fh2MatchdEdx);
2318
78451bcd 2319 fh2EledEdx = new TH2F("h2EledEdx","#it{dE/dx} vs. #it{p} for electrons",nptbins,ptmin,ptmax,ndedxbins,dedxmin,dedxmax);
2320 fh2EledEdx->SetXTitle("p (GeV/#it{c})");
2321 fh2EledEdx->SetYTitle("<#it{dE/dx}>");
3bfc4732 2322 outputContainer->Add(fh2EledEdx) ;
2323
78451bcd 2324 fh1EOverPR02 = new TH2F("h1EOverPR02","TRACK matches #it{E}/#it{p}, all",nptbins,ptmin,ptmax, nPoverEbins,eOverPmin,eOverPmax);
2325 fh1EOverPR02->SetYTitle("#it{E}/#it{p}");
2326 fh1EOverPR02->SetXTitle("#it{p}_{T} (GeV/#it{c})");
d55bb5e1 2327 outputContainer->Add(fh1EOverPR02);
a054a582 2328
78451bcd 2329 fh1EleEOverP = new TH2F("h1EleEOverP","Electron candidates #it{E}/#it{p} (60<#it{dE/dx}<100)",nptbins,ptmin,ptmax, nPoverEbins,eOverPmin,eOverPmax);
2330 fh1EleEOverP->SetYTitle("#it{E}/#it{p}");
2331 fh1EleEOverP->SetXTitle("#it{p}_{T} (GeV/#it{c})");
a054a582 2332 outputContainer->Add(fh1EleEOverP);
55c05f8c 2333 }
2334
f1538a5f 2335 if(fFillAllPi0Histo)
2336 {
9e9f04cb 2337 fhIM = new TH2F ("hIM","Cluster pairs Invariant mass vs reconstructed pair energy, ncell > 1",nptbins,ptmin,ptmax,nmassbins,massmin,massmax);
78451bcd 2338 fhIM->SetXTitle("#it{p}_{T, cluster pairs} (GeV) ");
2339 fhIM->SetYTitle("M_{cluster pairs} (GeV/#it{c}^{2})");
55c05f8c 2340 outputContainer->Add(fhIM);
a6f26052 2341
55c05f8c 2342 fhAsym = new TH2F ("hAssym","Cluster pairs Asymmetry vs reconstructed pair energy",nptbins,ptmin,ptmax,nasymbins,asymmin,asymmax);
78451bcd 2343 fhAsym->SetXTitle("#it{p}_{T, cluster pairs} (GeV) ");
2344 fhAsym->SetYTitle("#it{Asymmetry}");
55c05f8c 2345 outputContainer->Add(fhAsym);
a6f26052 2346 }
649b825d 2347
c8fe2783 2348
f1538a5f 2349 if(fFillAllPosHisto2)
2350 {
2351 if(fFillAllTH3)
2352 {
78451bcd 2353 fhXYZ = new TH3F ("hXYZ","Cluster: #it{x} vs #it{y} vs #it{z}",xbins,xmin,xmax,ybins,ymin,ymax,zbins,zmin,zmax);
2354 fhXYZ->SetXTitle("#it{x} (cm)");
2355 fhXYZ->SetYTitle("#it{y} (cm)");
2356 fhXYZ->SetZTitle("#it{z} (cm) ");
55c05f8c 2357 outputContainer->Add(fhXYZ);
2358 }
2359
35c71d5c 2360 fhXNCells = new TH2F ("hXNCells","Cluster X position vs N Cells per Cluster",xbins,xmin,xmax,nceclbins,nceclmin,nceclmax);
78451bcd 2361 fhXNCells->SetXTitle("#it{x} (cm)");
55c05f8c 2362 fhXNCells->SetYTitle("N cells per cluster");
2363 outputContainer->Add(fhXNCells);
2364
35c71d5c 2365 fhZNCells = new TH2F ("hZNCells","Cluster Z position vs N Cells per Cluster",zbins,zmin,zmax,nceclbins,nceclmin,nceclmax);
78451bcd 2366 fhZNCells->SetXTitle("#it{z} (cm)");
55c05f8c 2367 fhZNCells->SetYTitle("N cells per cluster");
2368 outputContainer->Add(fhZNCells);
2369
2370 fhXE = new TH2F ("hXE","Cluster X position vs cluster energy",xbins,xmin,xmax,nptbins,ptmin,ptmax);
78451bcd 2371 fhXE->SetXTitle("#it{x} (cm)");
2372 fhXE->SetYTitle("#it{E} (GeV)");
55c05f8c 2373 outputContainer->Add(fhXE);
2374
2375 fhZE = new TH2F ("hZE","Cluster Z position vs cluster energy",zbins,zmin,zmax,nptbins,ptmin,ptmax);
78451bcd 2376 fhZE->SetXTitle("#it{z} (cm)");
2377 fhZE->SetYTitle("#it{E} (GeV)");
55c05f8c 2378 outputContainer->Add(fhZE);
2379
35c71d5c 2380 fhRNCells = new TH2F ("hRNCells","Cluster R position vs N Cells per Cluster",rbins,rmin,rmax,nceclbins,nceclmin,nceclmax);
55c05f8c 2381 fhRNCells->SetXTitle("r = #sqrt{x^{2}+y^{2}} (cm)");
2382 fhRNCells->SetYTitle("N cells per cluster");
2383 outputContainer->Add(fhRNCells);
2384
2385
35c71d5c 2386 fhYNCells = new TH2F ("hYNCells","Cluster Y position vs N Cells per Cluster",ybins,ymin,ymax,nceclbins,nceclmin,nceclmax);
78451bcd 2387 fhYNCells->SetXTitle("#it{y} (cm)");
55c05f8c 2388 fhYNCells->SetYTitle("N cells per cluster");
2389 outputContainer->Add(fhYNCells);
2390
2391 fhRE = new TH2F ("hRE","Cluster R position vs cluster energy",rbins,rmin,rmax,nptbins,ptmin,ptmax);
2392 fhRE->SetXTitle("r = #sqrt{x^{2}+y^{2}} (cm)");
78451bcd 2393 fhRE->SetYTitle("#it{E} (GeV)");
55c05f8c 2394 outputContainer->Add(fhRE);
2395
2396 fhYE = new TH2F ("hYE","Cluster Y position vs cluster energy",ybins,ymin,ymax,nptbins,ptmin,ptmax);
78451bcd 2397 fhYE->SetXTitle("#it{y} (cm)");
2398 fhYE->SetYTitle("#it{E} (GeV)");
55c05f8c 2399 outputContainer->Add(fhYE);
2400 }
f1538a5f 2401
2402 if(fFillAllPosHisto)
2403 {
a6f26052 2404 fhRCellE = new TH2F ("hRCellE","Cell R position vs cell energy",rbins,rmin,rmax,nptbins,ptmin,ptmax);
2405 fhRCellE->SetXTitle("r = #sqrt{x^{2}+y^{2}} (cm)");
78451bcd 2406 fhRCellE->SetYTitle("#it{E} (GeV)");
a6f26052 2407 outputContainer->Add(fhRCellE);
2408
2409 fhXCellE = new TH2F ("hXCellE","Cell X position vs cell energy",xbins,xmin,xmax,nptbins,ptmin,ptmax);
78451bcd 2410 fhXCellE->SetXTitle("#it{x} (cm)");
2411 fhXCellE->SetYTitle("#it{E} (GeV)");
a6f26052 2412 outputContainer->Add(fhXCellE);
2413
2414 fhYCellE = new TH2F ("hYCellE","Cell Y position vs cell energy",ybins,ymin,ymax,nptbins,ptmin,ptmax);
78451bcd 2415 fhYCellE->SetXTitle("#it{y} (cm)");
2416 fhYCellE->SetYTitle("#it{E} (GeV)");
a6f26052 2417 outputContainer->Add(fhYCellE);
2418
2419 fhZCellE = new TH2F ("hZCellE","Cell Z position vs cell energy",zbins,zmin,zmax,nptbins,ptmin,ptmax);
78451bcd 2420 fhZCellE->SetXTitle("#it{z} (cm)");
2421 fhZCellE->SetYTitle("#it{E} (GeV)");
a6f26052 2422 outputContainer->Add(fhZCellE);
2423
78451bcd 2424 fhXYZCell = new TH3F ("hXYZCell","Cell : #it{x} vs #it{y} vs #it{z}",xbins,xmin,xmax,ybins,ymin,ymax,zbins,zmin,zmax);
2425 fhXYZCell->SetXTitle("#it{x} (cm)");
2426 fhXYZCell->SetYTitle("#it{y} (cm)");
2427 fhXYZCell->SetZTitle("#it{z} (cm)");
a6f26052 2428 outputContainer->Add(fhXYZCell);
2429
2430
2431 Float_t dx = TMath::Abs(xmin)+TMath::Abs(xmax);
2432 Float_t dy = TMath::Abs(ymin)+TMath::Abs(ymax);
2433 Float_t dz = TMath::Abs(zmin)+TMath::Abs(zmax);
2434 Float_t dr = TMath::Abs(rmin)+TMath::Abs(rmax);
2435
35c71d5c 2436 fhDeltaCellClusterRNCells = new TH2F ("hDeltaCellClusterRNCells","Cluster-Cell R position vs N Cells per Cluster",rbins*2,-dr,dr,nceclbins,nceclmin,nceclmax);
a6f26052 2437 fhDeltaCellClusterRNCells->SetXTitle("r = #sqrt{x^{2}+y^{2}} (cm)");
78451bcd 2438 fhDeltaCellClusterRNCells->SetYTitle("#it{n}_{cells per cluster}");
a6f26052 2439 outputContainer->Add(fhDeltaCellClusterRNCells);
2440
35c71d5c 2441 fhDeltaCellClusterXNCells = new TH2F ("hDeltaCellClusterXNCells","Cluster-Cell X position vs N Cells per Cluster",xbins*2,-dx,dx,nceclbins,nceclmin,nceclmax);
78451bcd 2442 fhDeltaCellClusterXNCells->SetXTitle("#it{x} (cm)");
2443 fhDeltaCellClusterXNCells->SetYTitle("#it{n}_{cells per cluster}");
a6f26052 2444 outputContainer->Add(fhDeltaCellClusterXNCells);
2445
35c71d5c 2446 fhDeltaCellClusterYNCells = new TH2F ("hDeltaCellClusterYNCells","Cluster-Cell Y position vs N Cells per Cluster",ybins*2,-dy,dy,nceclbins,nceclmin,nceclmax);
78451bcd 2447 fhDeltaCellClusterYNCells->SetXTitle("#it{y} (cm)");
a6f26052 2448 fhDeltaCellClusterYNCells->SetYTitle("N cells per cluster");
2449 outputContainer->Add(fhDeltaCellClusterYNCells);
2450
35c71d5c 2451 fhDeltaCellClusterZNCells = new TH2F ("hDeltaCellClusterZNCells","Cluster-Cell Z position vs N Cells per Cluster",zbins*2,-dz,dz,nceclbins,nceclmin,nceclmax);
78451bcd 2452 fhDeltaCellClusterZNCells->SetXTitle("#it{z} (cm)");
2453 fhDeltaCellClusterZNCells->SetYTitle("#it{n}_{cells per cluster}");
a6f26052 2454 outputContainer->Add(fhDeltaCellClusterZNCells);
2455
2456 fhDeltaCellClusterRE = new TH2F ("hDeltaCellClusterRE","Cluster-Cell R position vs cluster energy",rbins*2,-dr,dr,nptbins,ptmin,ptmax);
2457 fhDeltaCellClusterRE->SetXTitle("r = #sqrt{x^{2}+y^{2}} (cm)");
78451bcd 2458 fhDeltaCellClusterRE->SetYTitle("#it{E} (GeV)");
a6f26052 2459 outputContainer->Add(fhDeltaCellClusterRE);
2460
2461 fhDeltaCellClusterXE = new TH2F ("hDeltaCellClusterXE","Cluster-Cell X position vs cluster energy",xbins*2,-dx,dx,nptbins,ptmin,ptmax);
78451bcd 2462 fhDeltaCellClusterXE->SetXTitle("#it{x} (cm)");
2463 fhDeltaCellClusterXE->SetYTitle("#it{E} (GeV)");
a6f26052 2464 outputContainer->Add(fhDeltaCellClusterXE);
2465
2466 fhDeltaCellClusterYE = new TH2F ("hDeltaCellClusterYE","Cluster-Cell Y position vs cluster energy",ybins*2,-dy,dy,nptbins,ptmin,ptmax);
78451bcd 2467 fhDeltaCellClusterYE->SetXTitle("#it{y} (cm)");
2468 fhDeltaCellClusterYE->SetYTitle("#it{E} (GeV)");
a6f26052 2469 outputContainer->Add(fhDeltaCellClusterYE);
2470
2471 fhDeltaCellClusterZE = new TH2F ("hDeltaCellClusterZE","Cluster-Cell Z position vs cluster energy",zbins*2,-dz,dz,nptbins,ptmin,ptmax);
78451bcd 2472 fhDeltaCellClusterZE->SetXTitle("#it{z} (cm)");
2473 fhDeltaCellClusterZE->SetYTitle("#it{E} (GeV)");
a6f26052 2474 outputContainer->Add(fhDeltaCellClusterZE);
2475
2476 fhEtaPhiAmp = new TH3F ("hEtaPhiAmp","Cell #eta vs cell #phi vs cell energy",netabins,etamin,etamax,nphibins,phimin,phimax,nptbins,ptmin,ptmax);
2477 fhEtaPhiAmp->SetXTitle("#eta ");
2478 fhEtaPhiAmp->SetYTitle("#phi (rad)");
78451bcd 2479 fhEtaPhiAmp->SetZTitle("#it{E} (GeV) ");
a6f26052 2480 outputContainer->Add(fhEtaPhiAmp);
2481
2302a644 2482 }
c8fe2783 2483
a6f26052 2484 //Calo cells
e3300002 2485 fhNCells = new TH1F ("hNCells","# cells", ncebins,ncemin+0.5,ncemax);
78451bcd 2486 fhNCells->SetXTitle("#it{n}_{cells}");
2302a644 2487 outputContainer->Add(fhNCells);
701cbf54 2488
2489 fhNCellsCutAmpMin = new TH1F ("hNCellsCutAmpMin",Form("# cells amp > %1.2f-%1.2f",fEMCALCellAmpMin,fPHOSCellAmpMin), ncebins,ncemin+0.5,ncemax);
78451bcd 2490 fhNCellsCutAmpMin->SetXTitle("#it{n}_{cells}");
701cbf54 2491 outputContainer->Add(fhNCellsCutAmpMin);
2302a644 2492
78451bcd 2493 fhAmplitude = new TH1F ("hAmplitude","#it{E}_{cell}", nptbins*2,ptmin,ptmax);
2494 fhAmplitude->SetXTitle("#it{E}_{cell} (GeV)");
2302a644 2495 outputContainer->Add(fhAmplitude);
2496
78451bcd 2497 fhAmpId = new TH2F ("hAmpId","#it{E}_{cell}", nfineptbins,ptfinemin,ptfinemax,fNMaxRows*fNMaxCols*fNModules,0,fNMaxRows*fNMaxCols*fNModules);
2498 fhAmpId->SetXTitle("#it{E}_{cell} (GeV)");
2302a644 2499 outputContainer->Add(fhAmpId);
638916c4 2500
2501 fhAmpIdLowGain = new TH2F ("hAmpIdLG","Low gain: #it{E}_{cell}", nfineptbins,ptfinemin,ptfinemax,fNMaxRows*fNMaxCols*fNModules,0,fNMaxRows*fNMaxCols*fNModules);
2502 fhAmpIdLowGain->SetXTitle("#it{E}_{cell} (GeV)");
2503 outputContainer->Add(fhAmpIdLowGain);
2302a644 2504
638916c4 2505 if(fFillAllCellTimeHisto)
e6fec6f5 2506 {
649b825d 2507 fhCellTimeSpreadRespectToCellMax = new TH2F ("hCellTimeSpreadRespectToCellMax","t_{cell max}-t_{cell i} per cluster", nptbins,ptmin,ptmax,tdbins,tdmin,tdmax);
78451bcd 2508 fhCellTimeSpreadRespectToCellMax->SetXTitle("#it{E} (GeV)");
2509 fhCellTimeSpreadRespectToCellMax->SetYTitle("#Delta #it{t}_{cell max-i} (ns)");
2302a644 2510 outputContainer->Add(fhCellTimeSpreadRespectToCellMax);
2511
649b825d 2512 fhClusterMaxCellDiffAverageTime = new TH2F ("hClusterMaxCellDiffAverageTime","t_{cell max}-t_{average} per cluster", nptbins,ptmin,ptmax, tdbins,tdmin,tdmax);
78451bcd 2513 fhClusterMaxCellDiffAverageTime->SetXTitle("#it{E} (GeV)");
2514 fhClusterMaxCellDiffAverageTime->SetYTitle("#Delta #it{t}_{cell max - average} (ns)");
9e9f04cb 2515 outputContainer->Add(fhClusterMaxCellDiffAverageTime);
a82b4462 2516
649b825d 2517 fhClusterMaxCellDiffWeightedTime = new TH2F ("hClusterMaxCellDiffWeightedTime","t_{cell max}-t_{weighted} per cluster", nptbins,ptmin,ptmax, tdbins,tdmin,tdmax);
78451bcd 2518 fhClusterMaxCellDiffWeightedTime->SetXTitle("#it{E} (GeV)");
2519 fhClusterMaxCellDiffWeightedTime->SetYTitle("#Delta #it{t}_{cell max - weighted} (ns)");
649b825d 2520 outputContainer->Add(fhClusterMaxCellDiffWeightedTime);
35c71d5c 2521
924e319f 2522 fhCellIdCellLargeTimeSpread= new TH1F ("hCellIdCellLargeTimeSpread","Cells with time 100 ns larger than cell max in cluster ",
35c71d5c 2523 fNMaxCols*fNMaxRows*fNModules,0,fNMaxCols*fNMaxRows*fNModules);
2302a644 2524 fhCellIdCellLargeTimeSpread->SetXTitle("Absolute Cell Id");
2525 outputContainer->Add(fhCellIdCellLargeTimeSpread);
521636d2 2526
78451bcd 2527 fhTime = new TH1F ("hTime","#it{t}_{cell}",ntimebins,timemin,timemax);
2528 fhTime->SetXTitle("#it{t}_{cell} (ns)");
2302a644 2529 outputContainer->Add(fhTime);
2530
78451bcd 2531 fhTimeVz = new TH2F ("hTimeVz","#it{t}_{cell} vs vertex, amplitude > 0.5 GeV",100, 0, 50,ntimebins,timemin,timemax);
1a72f6c5 2532 fhTimeVz->SetXTitle("|v_{z}| (cm)");
78451bcd 2533 fhTimeVz->SetYTitle("#it{t}_{cell} (ns)");
1a72f6c5 2534 outputContainer->Add(fhTimeVz);
2535
78451bcd 2536 fhTimeId = new TH2F ("hTimeId","#it{t}_{cell} vs Absolute Id",
35c71d5c 2537 ntimebins,timemin,timemax,fNMaxRows*fNMaxCols*fNModules,0,fNMaxRows*fNMaxCols*fNModules);
78451bcd 2538 fhTimeId->SetXTitle("#it{t}_{cell} (ns)");
2302a644 2539 fhTimeId->SetYTitle("Cell Absolute Id");
2540 outputContainer->Add(fhTimeId);
2541
78451bcd 2542 fhTimeAmp = new TH2F ("hTimeAmp","#it{t}_{cell} vs #it{E}_{cell}",nptbins*2,ptmin,ptmax,ntimebins,timemin,timemax);
2543 fhTimeAmp->SetYTitle("#it{t}_{cell} (ns)");
2544 fhTimeAmp->SetXTitle("#it{E}_{cell} (GeV)");
2302a644 2545 outputContainer->Add(fhTimeAmp);
2546
638916c4 2547 fhTimeIdLowGain = new TH2F ("hTimeIdLG","Low gain: #it{t}_{cell} vs Absolute Id",
2548 ntimebins,timemin,timemax,fNMaxRows*fNMaxCols*fNModules,0,fNMaxRows*fNMaxCols*fNModules);
2549 fhTimeIdLowGain->SetXTitle("#it{t}_{cell} (ns)");
2550 fhTimeIdLowGain->SetYTitle("Cell Absolute Id");
2551 outputContainer->Add(fhTimeIdLowGain);
2552
2553 fhTimeAmpLowGain = new TH2F ("hTimeAmpLG","Low gain: #it{t}_{cell} vs #it{E}_{cell}",nptbins*2,ptmin,ptmax,ntimebins,timemin,timemax);
2554 fhTimeAmpLowGain->SetYTitle("#it{t}_{cell} (ns)");
2555 fhTimeAmpLowGain->SetXTitle("#it{E}_{cell} (GeV)");
2556 outputContainer->Add(fhTimeAmpLowGain);
2557
2302a644 2558 }
06f1b12a 2559
2560 fhCellECross = new TH2F ("hCellECross","1 - Energy in cross around cell / cell energy",
2561 nptbins,ptmin,ptmax, 400,-1,1.);
78451bcd 2562 fhCellECross->SetXTitle("#it{E}_{cell} (GeV) ");
2563 fhCellECross->SetYTitle("1- #it{E}_{cross}/#it{E}_{cell}");
06f1b12a 2564 outputContainer->Add(fhCellECross);
2565
1a72f6c5 2566
f1538a5f 2567 if(fCorrelate)
2568 {
798a9b04 2569 //PHOS vs EMCAL
35c71d5c 2570 fhCaloCorrNClusters = new TH2F ("hCaloCorrNClusters","# clusters in EMCAL vs PHOS", nclbins,nclmin,nclmax,nclbins,nclmin,nclmax);
2302a644 2571 fhCaloCorrNClusters->SetXTitle("number of clusters in EMCAL");
2572 fhCaloCorrNClusters->SetYTitle("number of clusters in PHOS");
2573 outputContainer->Add(fhCaloCorrNClusters);
2574
653aed3c 2575 fhCaloCorrEClusters = new TH2F ("hCaloCorrEClusters","summed energy of clusters in EMCAL vs PHOS", nptbins,ptmin,ptmax*2,nptbins,ptmin,ptmax*2);
78451bcd 2576 fhCaloCorrEClusters->SetXTitle("#Sigma #it{E} of clusters in EMCAL (GeV)");
2577 fhCaloCorrEClusters->SetYTitle("#Sigma #it{E} of clusters in PHOS (GeV)");
2302a644 2578 outputContainer->Add(fhCaloCorrEClusters);
2579
35c71d5c 2580 fhCaloCorrNCells = new TH2F ("hCaloCorrNCells","# Cells in EMCAL vs PHOS", ncebins,ncemin,ncemax, ncebins,ncemin,ncemax);
649b825d 2581 fhCaloCorrNCells->SetXTitle("number of Cells in EMCAL");
2582 fhCaloCorrNCells->SetYTitle("number of Cells in PHOS");
2583 outputContainer->Add(fhCaloCorrNCells);
2302a644 2584
653aed3c 2585 fhCaloCorrECells = new TH2F ("hCaloCorrECells","summed energy of Cells in EMCAL vs PHOS", nptbins*2,ptmin,ptmax*4,nptbins*2,ptmin,ptmax*4);
78451bcd 2586 fhCaloCorrECells->SetXTitle("#Sigma #it{E} of Cells in EMCAL (GeV)");
2587 fhCaloCorrECells->SetYTitle("#Sigma #it{E} of Cells in PHOS (GeV)");
649b825d 2588 outputContainer->Add(fhCaloCorrECells);
2302a644 2589
649b825d 2590 //Calorimeter VS V0 signal
f3da05f3 2591 fhCaloV0SCorrNClusters = new TH2F ("hCaloV0SNClusters",Form("# clusters in %s vs V0 signal",GetCalorimeterString().Data()), nv0sbins,nv0smin,nv0smax,nclbins,nclmin,nclmax);
649b825d 2592 fhCaloV0SCorrNClusters->SetXTitle("V0 signal");
f3da05f3 2593 fhCaloV0SCorrNClusters->SetYTitle(Form("number of clusters in %s",GetCalorimeterString().Data()));
649b825d 2594 outputContainer->Add(fhCaloV0SCorrNClusters);
2302a644 2595
f3da05f3 2596 fhCaloV0SCorrEClusters = new TH2F ("hCaloV0SEClusters",Form("summed energy of clusters in %s vs V0 signal",GetCalorimeterString().Data()), nv0sbins,nv0smin,nv0smax,nptbins,ptmin,ptmax*2);
649b825d 2597 fhCaloV0SCorrEClusters->SetXTitle("V0 signal");
f3da05f3 2598 fhCaloV0SCorrEClusters->SetYTitle(Form("#Sigma #it{E} of clusters in %s (GeV)",GetCalorimeterString().Data()));
649b825d 2599 outputContainer->Add(fhCaloV0SCorrEClusters);
2302a644 2600
f3da05f3 2601 fhCaloV0SCorrNCells = new TH2F ("hCaloV0SNCells",Form("# Cells in %s vs V0 signal",GetCalorimeterString().Data()), nv0sbins,nv0smin,nv0smax, ncebins,ncemin,ncemax);
649b825d 2602 fhCaloV0SCorrNCells->SetXTitle("V0 signal");
f3da05f3 2603 fhCaloV0SCorrNCells->SetYTitle(Form("number of Cells in %s",GetCalorimeterString().Data()));
649b825d 2604 outputContainer->Add(fhCaloV0SCorrNCells);
3bfc4732 2605
f3da05f3 2606 fhCaloV0SCorrECells = new TH2F ("hCaloV0SECells",Form("summed energy of Cells in %s vs V0 signal",GetCalorimeterString().Data()), nv0sbins,nv0smin,nv0smax,nptbins,ptmin,ptmax*2);
649b825d 2607 fhCaloV0SCorrECells->SetXTitle("V0 signal");
f3da05f3 2608 fhCaloV0SCorrECells->SetYTitle(Form("#Sigma #it{E} of Cells in %s (GeV)",GetCalorimeterString().Data()));
649b825d 2609 outputContainer->Add(fhCaloV0SCorrECells);
3bfc4732 2610
649b825d 2611 //Calorimeter VS V0 multiplicity
f3da05f3 2612 fhCaloV0MCorrNClusters = new TH2F ("hCaloV0MNClusters",Form("# clusters in %s vs V0 signal",GetCalorimeterString().Data()), nv0mbins,nv0mmin,nv0mmax,nclbins,nclmin,nclmax);
649b825d 2613 fhCaloV0MCorrNClusters->SetXTitle("V0 signal");
f3da05f3 2614 fhCaloV0MCorrNClusters->SetYTitle(Form("number of clusters in %s",GetCalorimeterString().Data()));
649b825d 2615 outputContainer->Add(fhCaloV0MCorrNClusters);
3bfc4732 2616
f3da05f3 2617 fhCaloV0MCorrEClusters = new TH2F ("hCaloV0MEClusters",Form("summed energy of clusters in %s vs V0 signal",GetCalorimeterString().Data()), nv0mbins,nv0mmin,nv0mmax,nptbins,ptmin,ptmax*2);
649b825d 2618 fhCaloV0MCorrEClusters->SetXTitle("V0 signal");
f3da05f3 2619 fhCaloV0MCorrEClusters->SetYTitle(Form("#Sigma #it{E} of clusters in %s (GeV)",GetCalorimeterString().Data()));
649b825d 2620 outputContainer->Add(fhCaloV0MCorrEClusters);
3bfc4732 2621
f3da05f3 2622 fhCaloV0MCorrNCells = new TH2F ("hCaloV0MNCells",Form("# Cells in %s vs V0 signal",GetCalorimeterString().Data()), nv0mbins,nv0mmin,nv0mmax, ncebins,ncemin,ncemax);
649b825d 2623 fhCaloV0MCorrNCells->SetXTitle("V0 signal");
f3da05f3 2624 fhCaloV0MCorrNCells->SetYTitle(Form("number of Cells in %s",GetCalorimeterString().Data()));
649b825d 2625 outputContainer->Add(fhCaloV0MCorrNCells);
3bfc4732 2626
f3da05f3 2627 fhCaloV0MCorrECells = new TH2F ("hCaloV0MECells",Form("summed energy of Cells in %s vs V0 signal",GetCalorimeterString().Data()), nv0mbins,nv0mmin,nv0mmax,nptbins,ptmin,ptmax*2);
649b825d 2628 fhCaloV0MCorrECells->SetXTitle("V0 signal");
f3da05f3 2629 fhCaloV0MCorrECells->SetYTitle(Form("#Sigma #it{E} of Cells in %s (GeV)",GetCalorimeterString().Data()));
649b825d 2630 outputContainer->Add(fhCaloV0MCorrECells);
3bfc4732 2631
649b825d 2632 //Calorimeter VS Track multiplicity
f3da05f3 2633 fhCaloTrackMCorrNClusters = new TH2F ("hCaloTrackMNClusters",Form("# clusters in %s vs # tracks",GetCalorimeterString().Data()), ntrmbins,ntrmmin,ntrmmax,nclbins,nclmin,nclmax);
649b825d 2634 fhCaloTrackMCorrNClusters->SetXTitle("# tracks");
f3da05f3 2635 fhCaloTrackMCorrNClusters->SetYTitle(Form("number of clusters in %s",GetCalorimeterString().Data()));
649b825d 2636 outputContainer->Add(fhCaloTrackMCorrNClusters);
3bfc4732 2637
f3da05f3 2638 fhCaloTrackMCorrEClusters = new TH2F ("hCaloTrackMEClusters",Form("summed energy of clusters in %s vs # tracks",GetCalorimeterString().Data()), ntrmbins,ntrmmin,ntrmmax,nptbins,ptmin,ptmax*2);
649b825d 2639 fhCaloTrackMCorrEClusters->SetXTitle("# tracks");
f3da05f3 2640 fhCaloTrackMCorrEClusters->SetYTitle(Form("#Sigma #it{E} of clusters in %s (GeV)",GetCalorimeterString().Data()));
649b825d 2641 outputContainer->Add(fhCaloTrackMCorrEClusters);
3bfc4732 2642
f3da05f3 2643 fhCaloTrackMCorrNCells = new TH2F ("hCaloTrackMNCells",Form("# Cells in %s vs # tracks",GetCalorimeterString().Data()), ntrmbins,ntrmmin,ntrmmax, ncebins,ncemin,ncemax);
649b825d 2644 fhCaloTrackMCorrNCells->SetXTitle("# tracks");
f3da05f3 2645 fhCaloTrackMCorrNCells->SetYTitle(Form("number of Cells in %s",GetCalorimeterString().Data()));
649b825d 2646 outputContainer->Add(fhCaloTrackMCorrNCells);
3bfc4732 2647
f3da05f3 2648 fhCaloTrackMCorrECells = new TH2F ("hCaloTrackMECells",Form("summed energy of Cells in %s vs # tracks",GetCalorimeterString().Data()), ntrmbins,ntrmmin,ntrmmax,nptbins,ptmin,ptmax*2);
649b825d 2649 fhCaloTrackMCorrECells->SetXTitle("# tracks");
f3da05f3 2650 fhCaloTrackMCorrECells->SetYTitle(Form("#Sigma #it{E} of Cells in %s (GeV)",GetCalorimeterString().Data()));
649b825d 2651 outputContainer->Add(fhCaloTrackMCorrECells);
3bfc4732 2652
653aed3c 2653 fhCaloCenNClusters = new TH2F ("hCaloCenNClusters","# clusters in calorimeter vs centrality",100,0,100,nclbins,nclmin,nclmax);
2654 fhCaloCenNClusters->SetYTitle("number of clusters in calorimeter");
2655 fhCaloCenNClusters->SetXTitle("Centrality");
2656 outputContainer->Add(fhCaloCenNClusters);
2657
2658 fhCaloCenEClusters = new TH2F ("hCaloCenEClusters","summed energy of clusters in calorimeter vs centrality",100,0,100,nptbins,ptmin,ptmax*2);
78451bcd 2659 fhCaloCenEClusters->SetYTitle("#Sigma #it{E} of clusters in calorimeter (GeV)");
653aed3c 2660 fhCaloCenEClusters->SetXTitle("Centrality");
2661 outputContainer->Add(fhCaloCenEClusters);
2662
2663 fhCaloCenNCells = new TH2F ("hCaloCenNCells","# Cells in calorimeter vs centrality",100,0,100,ncebins,ncemin,ncemax);
2664 fhCaloCenNCells->SetYTitle("number of Cells in calorimeter");
2665 fhCaloCenNCells->SetXTitle("Centrality");
2666 outputContainer->Add(fhCaloCenNCells);
2667
2668 fhCaloCenECells = new TH2F ("hCaloCenECells","summed energy of Cells in calorimeter vs centrality",100,0,100,nptbins*2,ptmin,ptmax*4);
78451bcd 2669 fhCaloCenECells->SetYTitle("#Sigma #it{E} of Cells in calorimeter (GeV)");
653aed3c 2670 fhCaloCenECells->SetXTitle("Centrality");
2671 outputContainer->Add(fhCaloCenECells);
2672
2673 fhCaloEvPNClusters = new TH2F ("hCaloEvPNClusters","# clusters in calorimeter vs event plane angle",100,0,TMath::Pi(),nclbins,nclmin,nclmax);
2674 fhCaloEvPNClusters->SetYTitle("number of clusters in calorimeter");
2675 fhCaloEvPNClusters->SetXTitle("Event plane angle (rad)");
2676 outputContainer->Add(fhCaloEvPNClusters);
2677
2678 fhCaloEvPEClusters = new TH2F ("hCaloEvPEClusters","summed energy of clusters in calorimeter vs event plane angle",100,0,TMath::Pi(),nptbins,ptmin,ptmax*2);
78451bcd 2679 fhCaloEvPEClusters->SetYTitle("#Sigma #it{E} of clusters in calorimeter (GeV)");
653aed3c 2680 fhCaloEvPEClusters->SetXTitle("Event plane angle (rad)");
2681 outputContainer->Add(fhCaloEvPEClusters);
2682
2683 fhCaloEvPNCells = new TH2F ("hCaloEvPNCells","# Cells in calorimeter vs event plane angle",100,0,TMath::Pi(),ncebins,ncemin,ncemax);
2684 fhCaloEvPNCells->SetYTitle("number of Cells in calorimeter");
2685 fhCaloEvPNCells->SetXTitle("Event plane angle (rad)");
2686 outputContainer->Add(fhCaloEvPNCells);
2687
2688 fhCaloEvPECells = new TH2F ("hCaloEvPECells","summed energy of Cells in calorimeter vs event plane angle",100,0,TMath::Pi(),nptbins*2,ptmin,ptmax*4);
78451bcd 2689 fhCaloEvPECells->SetYTitle("#Sigma #it{E} of Cells in calorimeter (GeV)");
653aed3c 2690 fhCaloEvPECells->SetXTitle("Event plane angle (rad)");
2691 outputContainer->Add(fhCaloEvPECells);
2692
3bfc4732 2693
649b825d 2694 }//correlate calorimeters
c8fe2783 2695
649b825d 2696 //Module histograms
9725fd2a 2697
649b825d 2698 fhEMod = new TH2F ("hE_Mod","Cluster reconstructed Energy in each present Module",nptbins,ptmin,ptmax,fNModules,0,fNModules);
78451bcd 2699 fhEMod->SetXTitle("#it{E} (GeV)");
649b825d 2700 fhEMod->SetYTitle("Module");
2701 outputContainer->Add(fhEMod);
c8fe2783 2702
649b825d 2703 fhAmpMod = new TH2F ("hAmp_Mod","Cell energy in each present Module",nptbins,ptmin,ptmax,fNModules,0,fNModules);
78451bcd 2704 fhAmpMod->SetXTitle("#it{E} (GeV)");
649b825d 2705 fhAmpMod->SetYTitle("Module");
2706 outputContainer->Add(fhAmpMod);
3bfc4732 2707
e6fec6f5 2708 if(fFillAllCellTimeHisto)
2709 {
0fb69ade 2710 fhTimeMod = new TH2F ("hTime_Mod","Cell time in each present Module",ntimebins,timemin,timemax,fNModules,0,fNModules);
2711 fhTimeMod->SetXTitle("t (ns)");
2712 fhTimeMod->SetYTitle("Module");
2713 outputContainer->Add(fhTimeMod);
2714 }
3bfc4732 2715
e3300002 2716 fhNClustersMod = new TH2F ("hNClusters_Mod","# clusters vs Module", nclbins,nclmin+0.5,nclmax,fNModules,0,fNModules);
649b825d 2717 fhNClustersMod->SetXTitle("number of clusters");
2718 fhNClustersMod->SetYTitle("Module");
2719 outputContainer->Add(fhNClustersMod);
2302a644 2720
e3300002 2721 fhNCellsMod = new TH2F ("hNCells_Mod","# cells vs Module", ncebins,ncemin+0.5,ncemax,fNModules,0,fNModules);
78451bcd 2722 fhNCellsMod->SetXTitle("#it{n}_{cells}");
649b825d 2723 fhNCellsMod->SetYTitle("Module");
2724 outputContainer->Add(fhNCellsMod);
17708df9 2725
649b825d 2726 Int_t colmaxs = fNMaxCols;
2727 Int_t rowmaxs = fNMaxRows;
f3da05f3 2728 if(GetCalorimeter()==kEMCAL)
e6fec6f5 2729 {
649b825d 2730 colmaxs=2*fNMaxCols;
2731 rowmaxs=Int_t(fNModules/2)*fNMaxRows;
2732 }
e6fec6f5 2733 else
2734 {
649b825d 2735 rowmaxs=fNModules*fNMaxRows;
2302a644 2736 }
3bfc4732 2737
649b825d 2738 fhGridCells = new TH2F ("hGridCells",Form("Entries in grid of cells"),
2739 colmaxs+2,-1.5,colmaxs+0.5, rowmaxs+2,-1.5,rowmaxs+0.5);
2740 fhGridCells->SetYTitle("row (phi direction)");
2741 fhGridCells->SetXTitle("column (eta direction)");
2742 outputContainer->Add(fhGridCells);
3bfc4732 2743
649b825d 2744 fhGridCellsE = new TH2F ("hGridCellsE","Accumulated energy in grid of cells",
2745 colmaxs+2,-1.5,colmaxs+0.5, rowmaxs+2,-1.5,rowmaxs+0.5);
2746 fhGridCellsE->SetYTitle("row (phi direction)");
2747 fhGridCellsE->SetXTitle("column (eta direction)");
2748 outputContainer->Add(fhGridCellsE);
3bfc4732 2749
638916c4 2750 fhGridCellsLowGain = new TH2F ("hGridCellsLG",Form("Low gain: Entries in grid of cells"),
2751 colmaxs+2,-1.5,colmaxs+0.5, rowmaxs+2,-1.5,rowmaxs+0.5);
2752 fhGridCellsLowGain->SetYTitle("row (phi direction)");
2753 fhGridCellsLowGain->SetXTitle("column (eta direction)");
2754 outputContainer->Add(fhGridCellsLowGain);
2755
2756 fhGridCellsELowGain = new TH2F ("hGridCellsELG","Low gain: Accumulated energy in grid of cells",
2757 colmaxs+2,-1.5,colmaxs+0.5, rowmaxs+2,-1.5,rowmaxs+0.5);
2758 fhGridCellsELowGain->SetYTitle("row (phi direction)");
2759 fhGridCellsELowGain->SetXTitle("column (eta direction)");
2760 outputContainer->Add(fhGridCellsELowGain);
2761
2762
e6fec6f5 2763 if(fFillAllCellTimeHisto)
2764 {
638916c4 2765 fhGridCellsTime = new TH2F ("hGridCellsTime","Accumulated time in grid of cells",
2766 colmaxs+2,-1.5,colmaxs+0.5, rowmaxs+2,-1.5,rowmaxs+0.5);
0fb69ade 2767 fhGridCellsTime->SetYTitle("row (phi direction)");
2768 fhGridCellsTime->SetXTitle("column (eta direction)");
638916c4 2769 outputContainer->Add(fhGridCellsTime);
2770
2771 fhGridCellsTimeLowGain = new TH2F ("hGridCellsTimeLG","Low gain: Accumulated time in grid of cells",
2772 colmaxs+2,-1.5,colmaxs+0.5, rowmaxs+2,-1.5,rowmaxs+0.5);
2773 fhGridCellsTimeLowGain->SetYTitle("row (phi direction)");
2774 fhGridCellsTimeLowGain->SetXTitle("column (eta direction)");
2775 outputContainer->Add(fhGridCellsTimeLowGain);
0fb69ade 2776 }
3bfc4732 2777
e6fec6f5 2778 fhNCellsPerClusterMod = new TH2F*[fNModules];
649b825d 2779 fhNCellsPerClusterModNoCut = new TH2F*[fNModules];
e6fec6f5 2780 fhIMMod = new TH2F*[fNModules];
2781 if(fFillAllCellTimeHisto) fhTimeAmpPerRCU = new TH2F*[fNModules*fNRCU];
0fb69ade 2782
e6fec6f5 2783 for(Int_t imod = 0; imod < fNModules; imod++)
2784 {
649b825d 2785 fhNCellsPerClusterMod[imod] = new TH2F (Form("hNCellsPerCluster_Mod%d",imod),
2786 Form("# cells per cluster vs cluster energy in Module %d",imod),
2787 nptbins,ptmin,ptmax, nceclbins,nceclmin,nceclmax);
78451bcd 2788 fhNCellsPerClusterMod[imod]->SetXTitle("#it{E} (GeV)");
2789 fhNCellsPerClusterMod[imod]->SetYTitle("#it{n}_{cells}");
649b825d 2790 outputContainer->Add(fhNCellsPerClusterMod[imod]);
3bfc4732 2791
649b825d 2792 fhNCellsPerClusterModNoCut[imod] = new TH2F (Form("hNCellsPerClusterNoCut_Mod%d",imod),
2793 Form("# cells per cluster vs cluster energy in Module %d, no cut",imod),
2794 nptbins,ptmin,ptmax, nceclbins,nceclmin,nceclmax);
78451bcd 2795 fhNCellsPerClusterModNoCut[imod]->SetXTitle("#it{E} (GeV)");
2796 fhNCellsPerClusterModNoCut[imod]->SetYTitle("#it{n}_{cells}");
649b825d 2797 outputContainer->Add(fhNCellsPerClusterModNoCut[imod]);
3bfc4732 2798
e6fec6f5 2799 if(fFillAllCellTimeHisto)
2800 {
2801 for(Int_t ircu = 0; ircu < fNRCU; ircu++)
2802 {
649b825d 2803 fhTimeAmpPerRCU[imod*fNRCU+ircu] = new TH2F (Form("hTimeAmp_Mod%d_RCU%d",imod,ircu),
78451bcd 2804 Form("#it{E}_{cell} vs #it{t}_{cell} in Module %d, RCU %d ",imod,ircu),
649b825d 2805 nptbins,ptmin,ptmax,ntimebins,timemin,timemax);
78451bcd 2806 fhTimeAmpPerRCU[imod*fNRCU+ircu]->SetXTitle("#it{E} (GeV)");
2807 fhTimeAmpPerRCU[imod*fNRCU+ircu]->SetYTitle("#it{t} (ns)");
649b825d 2808 outputContainer->Add(fhTimeAmpPerRCU[imod*fNRCU+ircu]);
3bfc4732 2809
649b825d 2810 }
3bfc4732 2811 }
2812
e6fec6f5 2813 if(fFillAllPi0Histo)
2814 {
649b825d 2815 fhIMMod[imod] = new TH2F (Form("hIM_Mod%d",imod),
2816 Form("Cluster pairs Invariant mass vs reconstructed pair energy in Module %d, n cell > 1",imod),
2817 nptbins,ptmin,ptmax,nmassbins,massmin,massmax);
78451bcd 2818 fhIMMod[imod]->SetXTitle("#it{p}_{T, cluster pairs} (GeV) ");
2819 fhIMMod[imod]->SetYTitle("#it{M}_{cluster pairs} (GeV/#it{c}^{2})");
649b825d 2820 outputContainer->Add(fhIMMod[imod]);
3bfc4732 2821
649b825d 2822 }
2823 }
2824
a82b4462 2825 // Monte Carlo Histograms
649b825d 2826
6aba6683 2827 TString particleName[] = { "Photon","PhotonConv","Pi0", "Eta", "Electron", "NeutralHadron", "ChargedHadron" };
649b825d 2828
f1538a5f 2829 if(IsDataMC())
2830 {
6aba6683 2831 for(Int_t iPart = 0; iPart < 7; iPart++)
f1538a5f 2832 {
2833 for(Int_t iCh = 0; iCh < 2; iCh++)
2834 {
649b825d 2835 fhRecoMCRatioE[iPart][iCh] = new TH2F (Form("hRecoMCRatioE_%s_Match%d",particleName[iPart].Data(),iCh),
3fa37dff 2836 Form("Reconstructed/Generated E, %s, Matched %d",particleName[iPart].Data(),iCh),
649b825d 2837 nptbins, ptmin, ptmax, 200,0,2);
78451bcd 2838 fhRecoMCRatioE[iPart][iCh]->SetYTitle("#it{E}_{reconstructed}/#it{E}_{generated}");
2839 fhRecoMCRatioE[iPart][iCh]->SetXTitle("#it{E}_{reconstructed} (GeV)");
649b825d 2840 outputContainer->Add(fhRecoMCRatioE[iPart][iCh]);
3bfc4732 2841
3bfc4732 2842
649b825d 2843 fhRecoMCDeltaE[iPart][iCh] = new TH2F (Form("hRecoMCDeltaE_%s_Match%d",particleName[iPart].Data(),iCh),
3fa37dff 2844 Form("Generated - Reconstructed E, %s, Matched %d",particleName[iPart].Data(),iCh),
649b825d 2845 nptbins, ptmin, ptmax, nptbins*2,-ptmax,ptmax);
78451bcd 2846 fhRecoMCDeltaE[iPart][iCh]->SetYTitle("#Delta #it{E} (GeV)");
2847 fhRecoMCDeltaE[iPart][iCh]->SetXTitle("#it{E}_{reconstructed} (GeV)");
649b825d 2848 outputContainer->Add(fhRecoMCDeltaE[iPart][iCh]);
3bfc4732 2849
649b825d 2850 fhRecoMCDeltaPhi[iPart][iCh] = new TH2F (Form("hRecoMCDeltaPhi_%s_Match%d",particleName[iPart].Data(),iCh),
3fa37dff 2851 Form("Generated - Reconstructed #phi, %s, Matched %d",particleName[iPart].Data(),iCh),
649b825d 2852 nptbins, ptmin, ptmax, nphibins*2,-phimax,phimax);
3fa37dff 2853 fhRecoMCDeltaPhi[iPart][iCh]->SetYTitle("#Delta #phi (rad)");
78451bcd 2854 fhRecoMCDeltaPhi[iPart][iCh]->SetXTitle("#it{E}_{reconstructed} (GeV)");
649b825d 2855 outputContainer->Add(fhRecoMCDeltaPhi[iPart][iCh]);
3bfc4732 2856
649b825d 2857 fhRecoMCDeltaEta[iPart][iCh] = new TH2F (Form("hRecoMCDeltaEta_%s_Match%d",particleName[iPart].Data(),iCh),
3fa37dff 2858 Form("Generated - Reconstructed #eta, %s, Matched %d",particleName[iPart].Data(),iCh),
649b825d 2859 nptbins, ptmin, ptmax,netabins*2,-etamax,etamax);
3fa37dff 2860 fhRecoMCDeltaEta[iPart][iCh]->SetYTitle("#Delta #eta ");
78451bcd 2861 fhRecoMCDeltaEta[iPart][iCh]->SetXTitle("#it{E}_{reconstructed} (GeV)");
649b825d 2862 outputContainer->Add(fhRecoMCDeltaEta[iPart][iCh]);
e1e62b89 2863
649b825d 2864 fhRecoMCE[iPart][iCh] = new TH2F (Form("hRecoMCE_%s_Match%d",particleName[iPart].Data(),iCh),
78451bcd 2865 Form("#it{E} distribution, reconstructed vs generated, %s, Matched %d",particleName[iPart].Data(),iCh),
649b825d 2866 nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
78451bcd 2867 fhRecoMCE[iPart][iCh]->SetXTitle("#it{E}_{rec} (GeV)");
2868 fhRecoMCE[iPart][iCh]->SetYTitle("#it{E}_{gen} (GeV)");
649b825d 2869 outputContainer->Add(fhRecoMCE[iPart][iCh]);
3bfc4732 2870
649b825d 2871 fhRecoMCPhi[iPart][iCh] = new TH2F (Form("hRecoMCPhi_%s_Match%d",particleName[iPart].Data(),iCh),
2872 Form("#phi distribution, reconstructed vs generated, %s, Matched %d",particleName[iPart].Data(),iCh),
2873 nphibins,phimin,phimax, nphibins,phimin,phimax);
3fa37dff 2874 fhRecoMCPhi[iPart][iCh]->SetXTitle("#phi_{reconstructed} (rad)");
2875 fhRecoMCPhi[iPart][iCh]->SetYTitle("#phi_{generated} (rad)");
649b825d 2876 outputContainer->Add(fhRecoMCPhi[iPart][iCh]);
e1e62b89 2877
649b825d 2878 fhRecoMCEta[iPart][iCh] = new TH2F (Form("hRecoMCEta_%s_Match%d",particleName[iPart].Data(),iCh),
2879 Form("#eta distribution, reconstructed vs generated, %s, Matched %d",particleName[iPart].Data(),iCh),
2880 netabins,etamin,etamax,netabins,etamin,etamax);
3fa37dff 2881 fhRecoMCEta[iPart][iCh]->SetXTitle("#eta_{reconstructed} ");
2882 fhRecoMCEta[iPart][iCh]->SetYTitle("#eta_{generated} ");
649b825d 2883 outputContainer->Add(fhRecoMCEta[iPart][iCh]);
3bfc4732 2884 }
649b825d 2885 }
c8fe2783 2886
649b825d 2887 //Pure MC
f1538a5f 2888 for(Int_t iPart = 0; iPart < 4; iPart++)
2889 {
95aee5e1 2890 fhGenMCE [iPart] = new TH1F(Form("hGenMCE_%s",particleName[iPart].Data()) ,
2891 Form("#it{E} of generated %s",particleName[iPart].Data()),
2892 nptbins,ptmin,ptmax);
2893
2894 fhGenMCPt[iPart] = new TH1F(Form("hGenMCPt_%s",particleName[iPart].Data()) ,
78451bcd 2895 Form("#it{p}_{T} of generated %s",particleName[iPart].Data()),
649b825d 2896 nptbins,ptmin,ptmax);
95aee5e1 2897
649b825d 2898 fhGenMCEtaPhi[iPart] = new TH2F(Form("hGenMCEtaPhi_%s",particleName[iPart].Data()),
2899 Form("Y vs #phi of generated %s",particleName[iPart].Data()),
95aee5e1 2900 200,-1,1,360,0,TMath::TwoPi());
649b825d 2901
95aee5e1 2902 fhGenMCE [iPart] ->SetXTitle("#it{E} (GeV)");
2903 fhGenMCPt[iPart] ->SetXTitle("#it{p}_{T} (GeV/#it{c})");
649b825d 2904 fhGenMCEtaPhi[iPart]->SetXTitle("#eta");
2905 fhGenMCEtaPhi[iPart]->SetYTitle("#phi (rad)");
95aee5e1 2906
2907 outputContainer->Add(fhGenMCE [iPart]);
2908 outputContainer->Add(fhGenMCPt [iPart]);
649b825d 2909 outputContainer->Add(fhGenMCEtaPhi[iPart]);
521636d2 2910
521636d2 2911
95aee5e1 2912 fhGenMCAccE [iPart] = new TH1F(Form("hGenMCAccE_%s",particleName[iPart].Data()) ,
78451bcd 2913 Form("#it{E} of generated %s",particleName[iPart].Data()),
649b825d 2914 nptbins,ptmin,ptmax);
95aee5e1 2915 fhGenMCAccPt[iPart] = new TH1F(Form("hGenMCAccPt_%s",particleName[iPart].Data()) ,
2916 Form("#it{p}_{T} of generated %s",particleName[iPart].Data()),
2917 nptbins,ptmin,ptmax);
649b825d 2918 fhGenMCAccEtaPhi[iPart] = new TH2F(Form("hGenMCAccEtaPhi_%s",particleName[iPart].Data()),
2919 Form("Y vs #phi of generated %s",particleName[iPart].Data()),
2920 netabins,etamin,etamax,nphibins,phimin,phimax);
2921
95aee5e1 2922 fhGenMCAccE [iPart] ->SetXTitle("#it{E} (GeV)");
2923 fhGenMCAccPt[iPart] ->SetXTitle("#it{p}_{T} (GeV/#it{c})");
649b825d 2924 fhGenMCAccEtaPhi[iPart]->SetXTitle("#eta");
2925 fhGenMCAccEtaPhi[iPart]->SetYTitle("#phi (rad)");
2926
95aee5e1 2927 outputContainer->Add(fhGenMCAccE [iPart]);
2928 outputContainer->Add(fhGenMCAccPt [iPart]);
649b825d 2929 outputContainer->Add(fhGenMCAccEtaPhi[iPart]);
2930
2931 }
f5036bcb 2932
649b825d 2933 //Vertex of generated particles
2934
2935 fhEMVxyz = new TH2F ("hEMVxyz","Production vertex of reconstructed ElectroMagnetic particles",nvdistbins,vdistmin,vdistmax,nvdistbins,vdistmin,vdistmax);//,100,0,500);
78451bcd 2936 fhEMVxyz->SetXTitle("#it{v}_{x}");
2937 fhEMVxyz->SetYTitle("#it{v}_{y}");
649b825d 2938 //fhEMVxyz->SetZTitle("v_{z}");
2939 outputContainer->Add(fhEMVxyz);
2940
2941 fhHaVxyz = new TH2F ("hHaVxyz","Production vertex of reconstructed hadrons",nvdistbins,vdistmin,vdistmax,nvdistbins,vdistmin,vdistmax);//,100,0,500);
78451bcd 2942 fhHaVxyz->SetXTitle("#it{v}_{x}");
2943 fhHaVxyz->SetYTitle("#it{v}_{y}");
649b825d 2944 //fhHaVxyz->SetZTitle("v_{z}");
2945 outputContainer->Add(fhHaVxyz);
2946
2947 fhEMR = new TH2F ("hEMR","Distance to production vertex of reconstructed ElectroMagnetic particles vs E rec",nptbins,ptmin,ptmax,nvdistbins,vdistmin,vdistmax);
78451bcd 2948 fhEMR->SetXTitle("#it{E} (GeV)");
649b825d 2949 fhEMR->SetYTitle("TMath::Sqrt(v_{x}^{2}+v_{y}^{2})");
2950 outputContainer->Add(fhEMR);
2951
2952 fhHaR = new TH2F ("hHaR","Distance to production vertex of reconstructed Hadrons vs E rec",nptbins,ptmin,ptmax,nvdistbins,vdistmin,vdistmax);
78451bcd 2953 fhHaR->SetXTitle("#it{E} (GeV)");
649b825d 2954 fhHaR->SetYTitle("TMath::Sqrt(v_{x}^{2}+v_{y}^{2})");
2955 outputContainer->Add(fhHaR);
2956
2957
2958 //Track Matching
2959
78451bcd 2960 fhMCEle1EOverP = new TH2F("hMCEle1EOverP","TRACK matches #it{E}/#it{p}, MC electrons",nptbins,ptmin,ptmax, nPoverEbins,eOverPmin,eOverPmax);
2961 fhMCEle1EOverP->SetYTitle("#it{E}/#it{p}");
2962 fhMCEle1EOverP->SetXTitle("#it{p}_{T} (GeV/#it{c})");
d55bb5e1 2963 outputContainer->Add(fhMCEle1EOverP);
649b825d 2964
2965 fhMCEle1dR = new TH1F("hMCEle1dR","TRACK matches dR, MC electrons",ndRbins,dRmin,dRmax);
78451bcd 2966 fhMCEle1dR->SetXTitle("#Delta #it{R} (rad)");
649b825d 2967 outputContainer->Add(fhMCEle1dR) ;
2968
78451bcd 2969 fhMCEle2MatchdEdx = new TH2F("hMCEle2MatchdEdx","#it{dE/dx} vs. #it{p} for all matches, MC electrons",nptbins,ptmin,ptmax,ndedxbins,dedxmin,dedxmax);
2970 fhMCEle2MatchdEdx->SetXTitle("#it{p} (GeV/#it{c})");
2971 fhMCEle2MatchdEdx->SetYTitle("<#it{dE/dx}>");
649b825d 2972 outputContainer->Add(fhMCEle2MatchdEdx);
2973
78451bcd 2974 fhMCChHad1EOverP = new TH2F("hMCChHad1EOverP","TRACK matches #it{E}/#it{p}, MC charged hadrons",nptbins,ptmin,ptmax, nPoverEbins,eOverPmin,eOverPmax);
2975 fhMCChHad1EOverP->SetYTitle("#it{E}/#it{p}");
2976 fhMCChHad1EOverP->SetXTitle("#it{p}_{T} (GeV/#it{c})");
d55bb5e1 2977 outputContainer->Add(fhMCChHad1EOverP);
649b825d 2978
2979 fhMCChHad1dR = new TH1F("hMCChHad1dR","TRACK matches dR, MC charged hadrons",ndRbins,dRmin,dRmax);
2980 fhMCChHad1dR->SetXTitle("#Delta R (rad)");
2981 outputContainer->Add(fhMCChHad1dR) ;
2982
78451bcd 2983 fhMCChHad2MatchdEdx = new TH2F("hMCChHad2MatchdEdx","#it{dE/dx} vs. #it{p} for all matches, MC charged hadrons",nptbins,ptmin,ptmax,ndedxbins,dedxmin,dedxmax);
2984 fhMCChHad2MatchdEdx->SetXTitle("#it{p} (GeV/#it{c})");
2985 fhMCChHad2MatchdEdx->SetYTitle("#it{dE/dx}>");
649b825d 2986 outputContainer->Add(fhMCChHad2MatchdEdx);
2987
78451bcd 2988 fhMCNeutral1EOverP = new TH2F("hMCNeutral1EOverP","TRACK matches #it{E}/#it{p}, MC neutrals",nptbins,ptmin,ptmax, nPoverEbins,eOverPmin,eOverPmax);
2989 fhMCNeutral1EOverP->SetYTitle("#it{E}/#it{p}");
2990 fhMCNeutral1EOverP->SetXTitle("#it{p}_{T} (GeV/#it{c})");
d55bb5e1 2991 outputContainer->Add(fhMCNeutral1EOverP);
649b825d 2992
2993 fhMCNeutral1dR = new TH1F("hMCNeutral1dR","TRACK matches dR, MC neutrals",ndRbins,dRmin,dRmax);
78451bcd 2994 fhMCNeutral1dR->SetXTitle("#Delta #it{R} (rad)");
649b825d 2995 outputContainer->Add(fhMCNeutral1dR) ;
2996
78451bcd 2997 fhMCNeutral2MatchdEdx = new TH2F("hMCNeutral2MatchdEdx","#it{dE/dx} vs. #it{p} for all matches, MC neutrals",nptbins,ptmin,ptmax,ndedxbins,dedxmin,dedxmax);
2998 fhMCNeutral2MatchdEdx->SetXTitle("#it{p} (GeV/#it{c})");
2999 fhMCNeutral2MatchdEdx->SetYTitle("#it{dE/dx}>");
649b825d 3000 outputContainer->Add(fhMCNeutral2MatchdEdx);
3001
78451bcd 3002 fhMCEle1EOverPR02 = new TH2F("hMCEle1EOverPR02","TRACK matches #it{E}/#it{p}, MC electrons",nptbins,ptmin,ptmax, nPoverEbins,eOverPmin,eOverPmax);
3003 fhMCEle1EOverPR02->SetYTitle("#it{E}/#it{p}");
3004 fhMCEle1EOverPR02->SetXTitle("#it{p}_{T} (GeV/#it{c})");
d55bb5e1 3005 outputContainer->Add(fhMCEle1EOverPR02);
649b825d 3006
78451bcd 3007 fhMCChHad1EOverPR02 = new TH2F("hMCChHad1EOverPR02","TRACK matches #it{E}/#it{p}, MC charged hadrons",nptbins,ptmin,ptmax, nPoverEbins,eOverPmin,eOverPmax);
3008 fhMCChHad1EOverPR02->SetYTitle("#it{E}/#it{p}");
3009 fhMCChHad1EOverPR02->SetXTitle("#it{p}_{T} (GeV/#it{c})");
d55bb5e1 3010 outputContainer->Add(fhMCChHad1EOverPR02);
f5036bcb 3011
78451bcd 3012 fhMCNeutral1EOverPR02 = new TH2F("hMCNeutral1EOverPR02","TRACK matches #it{E}/#it{p}, MC neutrals",nptbins,ptmin,ptmax, nPoverEbins,eOverPmin,eOverPmax);
3013 fhMCNeutral1EOverPR02->SetYTitle("#it{E}/#it{p}");
3014 fhMCNeutral1EOverPR02->SetXTitle("#it{p}_{T} (GeV/#it{c})");
d55bb5e1 3015 outputContainer->Add(fhMCNeutral1EOverPR02);
7206b21b 3016
78451bcd 3017 fhMCEle1EleEOverP = new TH2F("hMCEle1EleEOverP","Electron candidates #it{E}/#it{p} (60<dEdx<100), MC electrons",nptbins,ptmin,ptmax, nPoverEbins,eOverPmin,eOverPmax);
3018 fhMCEle1EleEOverP->SetYTitle("#it{E}/#it{p}");
3019 fhMCEle1EleEOverP->SetXTitle("#it{p}_{T} (GeV/#it{c})");
7206b21b 3020 outputContainer->Add(fhMCEle1EleEOverP);
3021
78451bcd 3022 fhMCChHad1EleEOverP = new TH2F("hMCEle1EleEOverP","Electron candidates #it{E}/#it{p} (60<dEdx<100), MC charged hadrons",nptbins,ptmin,ptmax, nPoverEbins,eOverPmin,eOverPmax);
3023 fhMCChHad1EleEOverP->SetYTitle("#it{E}/#it{p}");
3024 fhMCChHad1EleEOverP->SetXTitle("#it{p}_{T} (GeV/#it{c})");
7206b21b 3025 outputContainer->Add(fhMCChHad1EleEOverP);
3026
78451bcd 3027 fhMCNeutral1EleEOverP = new TH2F("hMCNeutral1EleEOverP","Electron candidates #it{E}/#it{p} (60<dEdx<100), MC neutrals",nptbins,ptmin,ptmax, nPoverEbins,eOverPmin,eOverPmax);
3028 fhMCNeutral1EleEOverP->SetYTitle("#it{E}/#it{p}");
3029 fhMCNeutral1EleEOverP->SetXTitle("#it{p}_{T} (GeV/#it{c})");
7206b21b 3030 outputContainer->Add(fhMCNeutral1EleEOverP);
3031
521636d2 3032 }
f5036bcb 3033
649b825d 3034 // for(Int_t i = 0; i < outputContainer->GetEntries() ; i++)
3035 // printf("i=%d, name= %s\n",i,outputContainer->At(i)->GetName());
521636d2 3036
649b825d 3037 return outputContainer;
902aa95c 3038}
3039
b94e038e 3040//______________________________________________________________________________________
3041Float_t AliAnaCalorimeterQA::GetECross(Int_t absID, AliVCaloCells* cells, Float_t dtcut)
1a72f6c5 3042{
3043 // Get energy in cross axis around maximum cell, for EMCAL only
3044
06f1b12a 3045 Int_t icol =-1, irow=-1,iRCU = -1;
0cea6003 3046 Int_t imod = GetModuleNumberCellIndexes(absID, GetCalorimeter(), icol, irow, iRCU);
e3300002 3047
f3da05f3 3048 if(GetCalorimeter()==kEMCAL)
57d8227a 3049 {
06f1b12a 3050 //Get close cells index, energy and time, not in corners
e3300002 3051
3052 Int_t absID1 = -1;
3053 Int_t absID2 = -1;
3054
3055 if( irow < AliEMCALGeoParams::fgkEMCALRows-1) absID1 = GetCaloUtils()->GetEMCALGeometry()->GetAbsCellIdFromCellIndexes(imod, irow+1, icol);
3056 if( irow > 0 ) absID2 = GetCaloUtils()->GetEMCALGeometry()->GetAbsCellIdFromCellIndexes(imod, irow-1, icol);
2747966a 3057
3058 // In case of cell in eta = 0 border, depending on SM shift the cross cell index
3059 Int_t absID3 = -1;
3060 Int_t absID4 = -1;
3061
3062 if ( icol == AliEMCALGeoParams::fgkEMCALCols - 1 && !(imod%2) )
3063 {
e3300002 3064 absID3 = GetCaloUtils()->GetEMCALGeometry()-> GetAbsCellIdFromCellIndexes(imod+1, irow, 0);
3065 absID4 = GetCaloUtils()->GetEMCALGeometry()-> GetAbsCellIdFromCellIndexes(imod , irow, icol-1);
2747966a 3066 }
3067 else if( icol == 0 && imod%2 )
3068 {
e3300002 3069 absID3 = GetCaloUtils()->GetEMCALGeometry()-> GetAbsCellIdFromCellIndexes(imod , irow, icol+1);
3070 absID4 = GetCaloUtils()->GetEMCALGeometry()-> GetAbsCellIdFromCellIndexes(imod-1, irow, AliEMCALGeoParams::fgkEMCALCols-1);
2747966a 3071 }
3072 else
3073 {
e3300002 3074 if( icol < AliEMCALGeoParams::fgkEMCALCols-1 )
3075 absID3 = GetCaloUtils()->GetEMCALGeometry()-> GetAbsCellIdFromCellIndexes(imod, irow, icol+1);
3076 if( icol > 0 )
3077 absID4 = GetCaloUtils()->GetEMCALGeometry()-> GetAbsCellIdFromCellIndexes(imod, irow, icol-1);
2747966a 3078 }
06f1b12a 3079
3080 //Recalibrate cell energy if needed
3081 //Float_t ecell = cells->GetCellAmplitude(absID);
0cea6003 3082 //GetCaloUtils()->RecalibrateCellAmplitude(ecell,GetCalorimeter(), absID);
06f1b12a 3083 Double_t tcell = cells->GetCellTime(absID);
0cea6003 3084 GetCaloUtils()->RecalibrateCellTime(tcell, GetCalorimeter(), absID,GetReader()->GetInputEvent()->GetBunchCrossNumber());
06f1b12a 3085
3086 Float_t ecell1 = 0, ecell2 = 0, ecell3 = 0, ecell4 = 0;
3087 Double_t tcell1 = 0, tcell2 = 0, tcell3 = 0, tcell4 = 0;
3088
45769d5b 3089 if(absID1 > 0 )
2747966a 3090 {
06f1b12a 3091 ecell1 = cells->GetCellAmplitude(absID1);
0cea6003 3092 GetCaloUtils()->RecalibrateCellAmplitude(ecell1, GetCalorimeter(), absID1);
06f1b12a 3093 tcell1 = cells->GetCellTime(absID1);
0cea6003 3094 GetCaloUtils()->RecalibrateCellTime (tcell1, GetCalorimeter(), absID1,GetReader()->GetInputEvent()->GetBunchCrossNumber());
06f1b12a 3095 }
45769d5b 3096 if(absID2 > 0 )
2747966a 3097 {
06f1b12a 3098 ecell2 = cells->GetCellAmplitude(absID2);
0cea6003 3099 GetCaloUtils()->RecalibrateCellAmplitude(ecell2, GetCalorimeter(), absID2);
06f1b12a 3100 tcell2 = cells->GetCellTime(absID2);
0cea6003 3101 GetCaloUtils()->RecalibrateCellTime (tcell2, GetCalorimeter(), absID2, GetReader()->GetInputEvent()->GetBunchCrossNumber());
06f1b12a 3102 }
45769d5b 3103 if(absID3 > 0 )
2747966a 3104 {
06f1b12a 3105 ecell3 = cells->GetCellAmplitude(absID3);
0cea6003 3106 GetCaloUtils()->RecalibrateCellAmplitude(ecell3, GetCalorimeter(), absID3);
06f1b12a 3107 tcell3 = cells->GetCellTime(absID3);
0cea6003 3108 GetCaloUtils()->RecalibrateCellTime (tcell3, GetCalorimeter(), absID3, GetReader()->GetInputEvent()->GetBunchCrossNumber());
06f1b12a 3109 }
45769d5b 3110 if(absID4 > 0 )
2747966a 3111 {
06f1b12a 3112 ecell4 = cells->GetCellAmplitude(absID4);
0cea6003 3113 GetCaloUtils()->RecalibrateCellAmplitude(ecell4, GetCalorimeter(), absID4);
06f1b12a 3114 tcell4 = cells->GetCellTime(absID4);
0cea6003 3115 GetCaloUtils()->RecalibrateCellTime (tcell4, GetCalorimeter(), absID4, GetReader()->GetInputEvent()->GetBunchCrossNumber());
06f1b12a 3116 }
f1538a5f 3117
3118 if(TMath::Abs(tcell-tcell1)*1.e9 > dtcut) ecell1 = 0 ;
3119 if(TMath::Abs(tcell-tcell2)*1.e9 > dtcut) ecell2 = 0 ;
3120 if(TMath::Abs(tcell-tcell3)*1.e9 > dtcut) ecell3 = 0 ;
3121 if(TMath::Abs(tcell-tcell4)*1.e9 > dtcut) ecell4 = 0 ;
06f1b12a 3122
3123 return ecell1+ecell2+ecell3+ecell4;
1a72f6c5 3124 }
57d8227a 3125 else //PHOS
3126 {
06f1b12a 3127
3128 Int_t absId1 = -1, absId2 = -1, absId3 = -1, absId4 = -1;
3129
3130 Int_t relId1[] = { imod+1, 0, irow+1, icol };
3131 Int_t relId2[] = { imod+1, 0, irow-1, icol };
3132 Int_t relId3[] = { imod+1, 0, irow , icol+1 };
3133 Int_t relId4[] = { imod+1, 0, irow , icol-1 };
3134
3135 GetCaloUtils()->GetPHOSGeometry()->RelToAbsNumbering(relId1, absId1);
3136 GetCaloUtils()->GetPHOSGeometry()->RelToAbsNumbering(relId2, absId2);
3137 GetCaloUtils()->GetPHOSGeometry()->RelToAbsNumbering(relId3, absId3);
3138 GetCaloUtils()->GetPHOSGeometry()->RelToAbsNumbering(relId4, absId4);
3139
3140 Float_t ecell1 = 0, ecell2 = 0, ecell3 = 0, ecell4 = 0;
3141
3142 if(absId1 > 0 ) ecell1 = cells->GetCellAmplitude(absId1);
3143 if(absId2 > 0 ) ecell2 = cells->GetCellAmplitude(absId2);
3144 if(absId3 > 0 ) ecell3 = cells->GetCellAmplitude(absId3);
3145 if(absId4 > 0 ) ecell4 = cells->GetCellAmplitude(absId4);
3146
3147 return ecell1+ecell2+ecell3+ecell4;
3148
1a72f6c5 3149 }
3150
1a72f6c5 3151}
3152
f4282b09 3153//___________________________________________________________________________________________________________
3154void AliAnaCalorimeterQA::InvariantMassHistograms(Int_t iclus, Int_t nModule, const TObjArray* caloClusters,
a82b4462 3155 AliVCaloCells * cells)
649b825d 3156{
3157 // Fill Invariant mass histograms
c8fe2783 3158
2db10729 3159 AliDebug(1,"Start");
3748ffb5 3160
649b825d 3161 //Get vertex for photon momentum calculation and event selection
3162 Double_t v[3] = {0,0,0}; //vertex ;
1a83b960 3163 //GetReader()->GetVertex(v);
a6f26052 3164
649b825d 3165 Int_t nModule2 = -1;
649b825d 3166 Int_t nCaloClusters = caloClusters->GetEntriesFast();
a6f26052 3167
d07278cf 3168 for(Int_t jclus = iclus + 1 ; jclus < nCaloClusters ; jclus++)
3169 {
649b825d 3170 AliVCluster* clus2 = (AliVCluster*)caloClusters->At(jclus);
a82b4462 3171
3172 Float_t maxCellFraction = 0.;
45769d5b 3173 Int_t absIdMax = GetCaloUtils()->GetMaxEnergyCell(cells, clus2, maxCellFraction);
55c05f8c 3174
d07278cf 3175 // Try to rediuce background with a mild shower shape cut and no more than 1 maxima
3176 // in cluster and remove low energy clusters
3177 if( clus2->GetNCells() <= 1 || !IsGoodCluster(absIdMax,cells) ||
3178 GetCaloUtils()->GetNumberOfLocalMaxima(clus2,cells) > 1 ||
07e4c878 3179 clus2->GetM02() > 0.5 || clus2->E() < fMinInvMassECut ) continue;
c8fe2783 3180
649b825d 3181 //Get cluster kinematics
f4282b09 3182 clus2->GetMomentum(fClusterMomentum2,v);
c8fe2783 3183
649b825d 3184 //Check only certain regions
3185 Bool_t in2 = kTRUE;
f4282b09 3186 if(IsFiducialCutOn()) in2 = GetFiducialCut()->IsInFiducialCut(fClusterMomentum2.Eta(),fClusterMomentum2.Phi(),GetCalorimeter()) ;
649b825d 3187 if(!in2) continue;
2302a644 3188
649b825d 3189 //Get module of cluster
3190 nModule2 = GetModuleNumber(clus2);
c8fe2783 3191
649b825d 3192 //Fill histograms
c8fe2783 3193
649b825d 3194 //All modules
f4282b09 3195 fhIM ->Fill((fClusterMomentum+fClusterMomentum2).Pt(),(fClusterMomentum+fClusterMomentum2).M());
49214ef9 3196
649b825d 3197 //Single module
45769d5b 3198 if(nModule == nModule2 && nModule >= 0 && nModule < fNModules)
f4282b09 3199 fhIMMod[nModule]->Fill((fClusterMomentum+fClusterMomentum2).Pt(),(fClusterMomentum+fClusterMomentum2).M());
d07278cf 3200
c8fe2783 3201
649b825d 3202 //Asymetry histograms
f4282b09 3203 fhAsym->Fill((fClusterMomentum+fClusterMomentum2).Pt(),
3204 TMath::Abs((fClusterMomentum.E()-fClusterMomentum2.E())/(fClusterMomentum.E()+fClusterMomentum2.E())));
649b825d 3205
3206 }// 2nd cluster loop
3207
3208}
3209
3210//______________________________
3211void AliAnaCalorimeterQA::Init()
3212{
3213 //Check if the data or settings are ok
c8fe2783 3214
f3da05f3 3215 if(GetCalorimeter() != kPHOS && GetCalorimeter() !=kEMCAL)
3216 AliFatal(Form("Wrong calorimeter name <%s>", GetCalorimeterString().Data()));
521636d2 3217
b6002a83 3218 //if(GetReader()->GetDataType()== AliCaloTrackReader::kMC)
3219 // AliFatal("Analysis of reconstructed data, MC reader not aplicable");
649b825d 3220
3221}
3bfc4732 3222
649b825d 3223//________________________________________
3224void AliAnaCalorimeterQA::InitParameters()
3225{
3226 //Initialize the parameters of the analysis.
3227 AddToHistogramsName("AnaCaloQA_");
3228
9f48b3f0 3229 fNModules = 22; // set maximum to maximum number of EMCAL modules
649b825d 3230 fNRCU = 2; // set maximum number of RCU in EMCAL per SM
07e4c878 3231
e6fec6f5 3232 fTimeCutMin = -9999999;
3233 fTimeCutMax = 9999999;
07e4c878 3234
3235 fEMCALCellAmpMin = 0.2; // 200 MeV
3236 fPHOSCellAmpMin = 0.2; // 200 MeV
3237 fCellAmpMin = 0.2; // 200 MeV
3238 fMinInvMassECut = 0.5; // 500 MeV
c8fe2783 3239
f1538a5f 3240 // Exotic studies
3241 fExoNECrossCuts = 10 ;
3242 fExoNDTimeCuts = 4 ;
3243
3244 fExoDTimeCuts [0] = 1.e4 ; fExoDTimeCuts [1] = 50.0 ; fExoDTimeCuts [2] = 25.0 ; fExoDTimeCuts [3] = 10.0 ;
3245 fExoECrossCuts[0] = 0.80 ; fExoECrossCuts[1] = 0.85 ; fExoECrossCuts[2] = 0.90 ; fExoECrossCuts[3] = 0.92 ; fExoECrossCuts[4] = 0.94 ;
3246 fExoECrossCuts[5] = 0.95 ; fExoECrossCuts[6] = 0.96 ; fExoECrossCuts[7] = 0.97 ; fExoECrossCuts[8] = 0.98 ; fExoECrossCuts[9] = 0.99 ;
3247
649b825d 3248}
c8fe2783 3249
b94e038e 3250//_____________________________________________________________________________
3251Bool_t AliAnaCalorimeterQA::IsGoodCluster(Int_t absIdMax, AliVCaloCells* cells)
649b825d 3252{
3253 //Identify cluster as exotic or not
3254
06f1b12a 3255 if(!fStudyBadClusters) return kTRUE;
a82b4462 3256
f3da05f3 3257 if(GetCalorimeter()==kEMCAL)
2747966a 3258 {
06f1b12a 3259 if(!GetCaloUtils()->GetEMCALRecoUtils()->IsRejectExoticCluster())
2747966a 3260 {
06f1b12a 3261 return !( GetCaloUtils()->GetEMCALRecoUtils()->IsExoticCell(absIdMax,cells,(GetReader()->GetInputEvent())->GetBunchCrossNumber()) );
2747966a 3262 }
3263 else
3264 {
06f1b12a 3265 return kTRUE;
2747966a 3266 }
649b825d 3267 }
06f1b12a 3268 else // PHOS
3269 {
1a83b960 3270 Float_t ampMax = cells->GetCellAmplitude(absIdMax);
0cea6003 3271 GetCaloUtils()->RecalibrateCellAmplitude(ampMax, GetCalorimeter(), absIdMax);
2747966a 3272
3273 if(ampMax < 0.01) return kFALSE;
3274
1a83b960 3275 if(1-GetECross(absIdMax,cells)/ampMax > 0.95) return kFALSE;
3276 else return kTRUE;
06f1b12a 3277 }
3278
649b825d 3279}
17708df9 3280
a82b4462 3281//_________________________________________________________
649b825d 3282void AliAnaCalorimeterQA::Print(const Option_t * opt) const
3283{
3284 //Print some relevant parameters set for the analysis
3285 if(! opt)
3286 return;
521636d2 3287
649b825d 3288 printf("**** Print %s %s ****\n", GetName(), GetTitle() ) ;
745913ae 3289 AliAnaCaloTrackCorrBaseClass::Print(" ");
2302a644 3290
f3da05f3 3291 printf("Select Calorimeter %s \n",GetCalorimeterString().Data());
649b825d 3292 printf("Time Cut: %3.1f < TOF < %3.1f\n", fTimeCutMin, fTimeCutMax);
3293 printf("EMCAL Min Amplitude : %2.1f GeV/c\n", fEMCALCellAmpMin) ;
3294 printf("PHOS Min Amplitude : %2.1f GeV/c\n", fPHOSCellAmpMin) ;
07e4c878 3295 printf("Inv. Mass min. E clus : %2.1f GeV/c\n", fMinInvMassECut) ;
c8fe2783 3296
649b825d 3297}
3298
649b825d 3299//_____________________________________________________
3300void AliAnaCalorimeterQA::MakeAnalysisFillHistograms()
3301{
3302 //Fill Calorimeter QA histograms
c8fe2783 3303
649b825d 3304 //Play with the MC stack if available
701cbf54 3305 if(IsDataMC()) MCHistograms();
798a9b04 3306
95aee5e1 3307 // Correlate Calorimeters and V0 and track Multiplicity
3308 if(fCorrelate) Correlate();
3309
701cbf54 3310 //Get List with CaloClusters , calo Cells, init min amplitude
3311 TObjArray * caloClusters = NULL;
3312 AliVCaloCells * cells = 0x0;
f3da05f3 3313 if (GetCalorimeter() == kPHOS)
701cbf54 3314 {
3315 fCellAmpMin = fPHOSCellAmpMin;
3316 caloClusters = GetPHOSClusters();
3317 cells = GetPHOSCells();
3318 }
f3da05f3 3319 else if (GetCalorimeter() == kEMCAL)
701cbf54 3320 {
3321 fCellAmpMin = fEMCALCellAmpMin;
3322 caloClusters = GetEMCALClusters();
3323 cells = GetEMCALCells();
3324 }
3325 else
2db10729 3326 AliFatal(Form("AliAnaCalorimeterQA::MakeAnalysisFillHistograms() - Wrong calorimeter name <%s>, END", GetCalorimeterString().Data()));
798a9b04 3327
07e4c878 3328 if( !caloClusters || !cells )
45769d5b 3329 {
2db10729 3330 AliFatal(Form("AliAnaCalorimeterQA::MakeAnalysisFillHistograms() - No CaloClusters or CaloCells available"));
649b825d 3331 return; // trick coverity
c8fe2783 3332 }
649b825d 3333
07e4c878 3334 if(caloClusters->GetEntriesFast() == 0) return ;
3335
1a72f6c5 3336 //printf("QA: N cells %d, N clusters %d \n",cells->GetNumberOfCells(),caloClusters->GetEntriesFast());
3337
95aee5e1 3338 // Clusters
649b825d 3339 ClusterLoopHistograms(caloClusters,cells);
3340
3341 // Cells
3342 CellHistograms(cells);
3343
2db10729 3344 AliDebug(1,"End");
649b825d 3345
a0bb4dc0 3346}
3347
649b825d 3348//______________________________________
3349void AliAnaCalorimeterQA::MCHistograms()
3350{
3351 //Get the MC arrays and do some checks before filling MC histograms
9e9f04cb 3352
95aee5e1 3353 Int_t pdg = 0 ;
3354 Int_t status = 0 ;
3355 Int_t nprim = 0 ;
3356
3357 TParticle * primStack = 0;
3358 AliAODMCParticle * primAOD = 0;
649b825d 3359
95aee5e1 3360 // Get the ESD MC particles container
3361 AliStack * stack = 0;
3362 if( GetReader()->ReadStack() )
45769d5b 3363 {
95aee5e1 3364 stack = GetMCStack();
3365 if(!stack ) return;
3366 nprim = stack->GetNtrack();
3367 }
2302a644 3368
95aee5e1 3369 // Get the AOD MC particles container
3370 TClonesArray * mcparticles = 0;
3371 if( GetReader()->ReadAODMCParticles() )
3372 {
3373 mcparticles = GetReader()->GetAODMCParticles();
3374 if( !mcparticles ) return;
3375 nprim = mcparticles->GetEntriesFast();
3376 }
35c71d5c 3377
95aee5e1 3378 //printf("N primaries %d\n",nprim);
3379 for(Int_t i=0 ; i < nprim; i++)
2747966a 3380 {
95aee5e1 3381 if(GetReader()->AcceptOnlyHIJINGLabels() && !GetReader()->IsHIJINGLabel(i)) continue ;
35c71d5c 3382
95aee5e1 3383 // Get the generated particles, check that it is primary (not parton, resonance)
3384 // and get its momentum. Different way to recover from ESD or AOD
3385 if(GetReader()->ReadStack())
2747966a 3386 {
95aee5e1 3387 primStack = stack->Particle(i) ;
3388 if(!primStack)
3389 {
2db10729 3390 AliWarning("ESD primaries pointer not available!!");
95aee5e1 3391 continue;
3392 }
3393
3394 pdg = primStack->GetPdgCode();
3395 status = primStack->GetStatusCode();
35c71d5c 3396
95aee5e1 3397 //printf("Input: i %d, %s, pdg %d, status %d \n",i, primStack->GetName(), pdg, status);
35c71d5c 3398
95aee5e1 3399 if ( status > 11 ) continue; //Working for PYTHIA and simple generators, check for HERWIG, HIJING?
3400
3401 if ( primStack->Energy() == TMath::Abs(primStack->Pz()) ) continue ; //Protection against floating point exception
3402
3403 //printf("Take : i %d, %s, pdg %d, status %d \n",i, primStack->GetName(), pdg, status);
3404
3405 //Photon kinematics
f4282b09 3406 primStack->Momentum(fPrimaryMomentum);
35c71d5c 3407 }
2747966a 3408 else
3409 {
95aee5e1 3410 primAOD = (AliAODMCParticle *) mcparticles->At(i);
3411 if(!primAOD)
3412 {
2db10729 3413 AliWarning("AOD primaries pointer not available!!");
95aee5e1 3414 continue;
3415 }
3416
3417 pdg = primAOD->GetPdgCode();
3418 status = primAOD->GetStatus();
3419
3420 //printf("Input: i %d, %s, pdg %d, status %d \n",i, primAOD->GetName(), pdg, status);
3421
3422 if (!primAOD->IsPrimary()) continue; //accept all which is not MC transport generated. Don't know how to avoid partons
3423
3424 if ( status > 11 ) continue; //Working for PYTHIA and simple generators, check for HERWIG
3425
3426 if ( primAOD->E() == TMath::Abs(primAOD->Pz()) ) continue ; //Protection against floating point exception
3427
3428 //printf("Take : i %d, %s, pdg %d, status %d \n",i, primAOD->GetName(), pdg, status);
3429
3430 //kinematics
f4282b09 3431 fPrimaryMomentum.SetPxPyPzE(primAOD->Px(),primAOD->Py(),primAOD->Pz(),primAOD->E());
2302a644 3432 }
95aee5e1 3433
f4282b09 3434 Float_t eMC = fPrimaryMomentum.E();
95aee5e1 3435 if(eMC < 0.2) continue;
f4282b09 3436 Float_t ptMC = fPrimaryMomentum.E();
95aee5e1 3437
f4282b09 3438 Float_t etaMC = fPrimaryMomentum.Eta();
95aee5e1 3439 // Only particles in |eta| < 1
3440 if (TMath::Abs(etaMC) > 1) continue;
3441
f4282b09 3442 Float_t phiMC = fPrimaryMomentum.Phi();
95aee5e1 3443 if(phiMC < 0)
3444 phiMC += TMath::TwoPi();
3445
3446 Int_t mcIndex = -1;
3447 if (pdg==22) mcIndex = kmcPhoton;
3448 else if (pdg==111) mcIndex = kmcPi0;
3449 else if (pdg==221) mcIndex = kmcEta;
3450 else if (TMath::Abs(pdg)==11) mcIndex = kmcElectron;
3451
3452 if( mcIndex >=0 )
2747966a 3453 {
95aee5e1 3454 fhGenMCE [mcIndex]->Fill( eMC);
3455 fhGenMCPt[mcIndex]->Fill(ptMC);
3456 if(eMC > 0.5) fhGenMCEtaPhi[mcIndex]->Fill(etaMC,phiMC);
3457
3458 Bool_t inacceptance = kTRUE;
3459 // Check same fidutial borders as in data analysis on top of real acceptance if real was requested.
f4282b09 3460 if( IsFiducialCutOn() && !GetFiducialCut()->IsInFiducialCut(fPrimaryMomentum.Eta(),fPrimaryMomentum.Phi(),GetCalorimeter()) )
3461 inacceptance = kFALSE ;
95aee5e1 3462
3463 if(IsRealCaloAcceptanceOn()) // defined on base class
3464 {
3465 if(GetReader()->ReadStack() &&
0cea6003 3466 !GetCaloUtils()->IsMCParticleInCalorimeterAcceptance(GetCalorimeter(), primStack)) inacceptance = kFALSE ;
95aee5e1 3467 if(GetReader()->ReadAODMCParticles() &&
0cea6003 3468 !GetCaloUtils()->IsMCParticleInCalorimeterAcceptance(GetCalorimeter(), primAOD )) inacceptance = kFALSE ;
95aee5e1 3469 }
3470
3471 if(!inacceptance) continue;
3472
3473 fhGenMCAccE [mcIndex]->Fill( eMC);
3474 fhGenMCAccPt[mcIndex]->Fill(ptMC);
3475 if(eMC > 0.5) fhGenMCAccEtaPhi[mcIndex]->Fill(etaMC,phiMC);
3476
2302a644 3477 }
95aee5e1 3478
2302a644 3479 }
902aa95c 3480}
c8fe2783 3481
649b825d 3482//_________________________________________________________________________________
3483void AliAnaCalorimeterQA::WeightHistograms(AliVCluster *clus, AliVCaloCells* cells)
3484{
3485 // Calculate weights
3486
3487 // First recalculate energy in case non linearity was applied
3488 Float_t energy = 0;
701cbf54 3489 Float_t ampMax = 0;
3490 Float_t energyOrg = clus->E();
3491
3492 // Do study when there are enough cells in cluster
3493 if(clus->GetNCells() < 3) return ;
3494
2747966a 3495 for (Int_t ipos = 0; ipos < clus->GetNCells(); ipos++)
3496 {
649b825d 3497 Int_t id = clus->GetCellsAbsId()[ipos];
3498
3499 //Recalibrate cell energy if needed
3500 Float_t amp = cells->GetCellAmplitude(id);
0cea6003 3501 GetCaloUtils()->RecalibrateCellAmplitude(amp, GetCalorimeter(), id);
649b825d 3502
3503 energy += amp;
3504
2db10729 3505 if ( amp > ampMax )
649b825d 3506 ampMax = amp;
3507
3508 } // energy loop
3509
2db10729 3510 if ( energy <=0 )
2747966a 3511 {
2db10729 3512 AliWarning(Form("Wrong calculated energy %f",energy));
649b825d 3513 return;
3514 }
3515
701cbf54 3516 //Remove non lin correction
3517 clus->SetE(energy);
3518
649b825d 3519 fhEMaxCellClusterRatio ->Fill(energy,ampMax/energy);
3520 fhEMaxCellClusterLogRatio->Fill(energy,TMath::Log(ampMax/energy));
3521
3522 //Get the ratio and log ratio to all cells in cluster
2747966a 3523 for (Int_t ipos = 0; ipos < clus->GetNCells(); ipos++)
3524 {
649b825d 3525 Int_t id = clus->GetCellsAbsId()[ipos];
3526
3527 //Recalibrate cell energy if needed
3528 Float_t amp = cells->GetCellAmplitude(id);
0cea6003 3529 GetCaloUtils()->RecalibrateCellAmplitude(amp, GetCalorimeter(), id);
649b825d 3530
3531 fhECellClusterRatio ->Fill(energy,amp/energy);
3532 fhECellClusterLogRatio->Fill(energy,TMath::Log(amp/energy));
3533 }
3534
3535 //Recalculate shower shape for different W0
f3da05f3 3536 if(GetCalorimeter()==kEMCAL)
2747966a 3537 {
649b825d 3538 Float_t l0org = clus->GetM02();
3539 Float_t l1org = clus->GetM20();
3540 Float_t dorg = clus->GetDispersion();
701cbf54 3541
3542 Int_t tagMC = -1;
3543 if(IsDataMC() && clus->GetNLabels() > 0)
45769d5b 3544 {
0cea6003 3545 Int_t tag = GetMCAnalysisUtils()->CheckOrigin(clus->GetLabels(),clus->GetNLabels(), GetReader(),GetCalorimeter());
701cbf54 3546
3547 if( GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCPhoton) &&
3548 !GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCPi0) &&
3549 !GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCEta) &&
3550 !GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCConversion) ){
3551 tagMC = 0;
3552 } // Pure Photon
3553 else if( GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCElectron) &&
3554 !GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCConversion) ){
3555 tagMC = 1;
3556 } // Electron
3557 else if( GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCConversion) ){
3558 tagMC = 2;
3559 } // Conversion
3560 else if( GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCPi0) ){
3561 tagMC = 3;
3562 }// Pi0
3563 else if(!GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCEta) &&
3564 !GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCPhoton) ){
3565 tagMC = 4;
3566 }// Hadron
3567 }// Is MC
3568
3569 for(Int_t iw = 0; iw < 12; iw++)
3570 {
3571 GetCaloUtils()->GetEMCALRecoUtils()->SetW0(3+iw*0.25);
649b825d 3572 GetCaloUtils()->GetEMCALRecoUtils()->RecalculateClusterShowerShapeParameters(GetEMCALGeometry(), cells, clus);
3573
3574 fhLambda0ForW0[iw]->Fill(energy,clus->GetM02());
1a72f6c5 3575 //fhLambda1ForW0[iw]->Fill(energy,clus->GetM20());
649b825d 3576
701cbf54 3577 if(IsDataMC() && tagMC >= 0)
45769d5b 3578 {
701cbf54 3579 fhLambda0ForW0MC[iw][tagMC]->Fill(energy,clus->GetM02());
3580 //fhLambda1ForW0MC[iw][tagMC]->Fill(energy,clus->GetM20());
3581 }
649b825d 3582 } // w0 loop
3583
3584 // Set the original values back
3585 clus->SetM02(l0org);
3586 clus->SetM20(l1org);
3587 clus->SetDispersion(dorg);
3588
3589 }// EMCAL
701cbf54 3590
3591 clus->SetE(energyOrg);
3592
649b825d 3593}
3594
3595
3596