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