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