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