You need cstdlib to declare exit()
[u/mrichter/AliRoot.git] / PWGGA / CaloTrackCorrelations / AliAnaPhoton.cxx
CommitLineData
85c4406e 1/**************************************************************************
1c5acb87 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 hereby granted *
cadbb0f3 9 * without fee, provided that the above copyright notice appears in all *
1c5acb87 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 **************************************************************************/
1c5acb87 15
16//_________________________________________________________________________
17//
18// Class for the photon identification.
19// Clusters from calorimeters are identified as photons
20// and kept in the AOD. Few histograms produced.
85c4406e 21// Produces input for other analysis classes like AliAnaPi0,
22// AliAnaParticleHadronCorrelation ...
1c5acb87 23//
85c4406e 24// -- Author: Gustavo Conesa (LNF-INFN)
1c5acb87 25//////////////////////////////////////////////////////////////////////////////
85c4406e 26
27
28// --- ROOT system ---
1c5acb87 29#include <TH2F.h>
477d6cee 30#include <TClonesArray.h>
0c1383b5 31#include <TObjString.h>
123fc3bd 32#include "TParticle.h"
6175da48 33#include "TDatabasePDG.h"
1c5acb87 34
85c4406e 35// --- Analysis system ---
36#include "AliAnaPhoton.h"
1c5acb87 37#include "AliCaloTrackReader.h"
123fc3bd 38#include "AliStack.h"
1c5acb87 39#include "AliCaloPID.h"
6639984f 40#include "AliMCAnalysisUtils.h"
ff45398a 41#include "AliFiducialCut.h"
0ae57829 42#include "AliVCluster.h"
591cc579 43#include "AliAODMCParticle.h"
c8fe2783 44#include "AliMixedEvent.h"
fc195fd0 45#include "AliAODEvent.h"
2ad19c3d 46#include "AliESDEvent.h"
c8fe2783 47
85c4406e 48// --- Detectors ---
c5693f62 49#include "AliPHOSGeoUtils.h"
50#include "AliEMCALGeometry.h"
1c5acb87 51
52ClassImp(AliAnaPhoton)
521636d2 53
85c4406e 54//____________________________
55AliAnaPhoton::AliAnaPhoton() :
56AliAnaCaloTrackCorrBaseClass(), fCalorimeter(""),
57fMinDist(0.), fMinDist2(0.), fMinDist3(0.),
58fRejectTrackMatch(0), fFillTMHisto(kFALSE),
59fTimeCutMin(-10000), fTimeCutMax(10000),
60fNCellsCut(0),
61fNLMCutMin(-1), fNLMCutMax(10),
62fFillSSHistograms(kFALSE), fFillOnlySimpleSSHisto(1),
bc41680b 63fFillPileUpHistograms(0),
85c4406e 64fNOriginHistograms(8), fNPrimaryHistograms(4),
85c4406e 65// Histograms
80b5ae80 66
67// Control histograms
68fhNCellsE(0), fhCellsE(0),
69fhMaxCellDiffClusterE(0), fhTimePt(0), fhEtaPhi(0),
70
85c4406e 71fhEPhoton(0), fhPtPhoton(0),
80b5ae80 72fhPhiPhoton(0), fhEtaPhoton(0),
73fhEtaPhiPhoton(0), fhEtaPhi05Photon(0),
85c4406e 74fhPtCentralityPhoton(0), fhPtEventPlanePhoton(0),
75
76// Shower shape histograms
77fhNLocMax(0),
78fhDispE(0), fhLam0E(0), fhLam1E(0),
79fhDispETRD(0), fhLam0ETRD(0), fhLam1ETRD(0),
80fhDispETM(0), fhLam0ETM(0), fhLam1ETM(0),
81fhDispETMTRD(0), fhLam0ETMTRD(0), fhLam1ETMTRD(0),
82
83fhNCellsLam0LowE(0), fhNCellsLam1LowE(0), fhNCellsDispLowE(0),
84fhNCellsLam0HighE(0), fhNCellsLam1HighE(0), fhNCellsDispHighE(0),
85
86fhEtaLam0LowE(0), fhPhiLam0LowE(0),
87fhEtaLam0HighE(0), fhPhiLam0HighE(0),
88fhLam0DispLowE(0), fhLam0DispHighE(0),
89fhLam1Lam0LowE(0), fhLam1Lam0HighE(0),
90fhDispLam1LowE(0), fhDispLam1HighE(0),
91fhDispEtaE(0), fhDispPhiE(0),
92fhSumEtaE(0), fhSumPhiE(0), fhSumEtaPhiE(0),
93fhDispEtaPhiDiffE(0), fhSphericityE(0),
94fhDispSumEtaDiffE(0), fhDispSumPhiDiffE(0),
95
96// MC histograms
97fhMCPhotonELambda0NoOverlap(0), fhMCPhotonELambda0TwoOverlap(0), fhMCPhotonELambda0NOverlap(0),
98// Embedding
99fhEmbeddedSignalFractionEnergy(0),
100fhEmbedPhotonELambda0FullSignal(0), fhEmbedPhotonELambda0MostlySignal(0),
101fhEmbedPhotonELambda0MostlyBkg(0), fhEmbedPhotonELambda0FullBkg(0),
102fhEmbedPi0ELambda0FullSignal(0), fhEmbedPi0ELambda0MostlySignal(0),
103fhEmbedPi0ELambda0MostlyBkg(0), fhEmbedPi0ELambda0FullBkg(0),
bc41680b 104
b2e375c7 105fhTimePtPhotonNoCut(0), fhTimePtPhotonSPD(0),
85c4406e 106fhTimeNPileUpVertSPD(0), fhTimeNPileUpVertTrack(0),
85c4406e 107fhPtPhotonNPileUpSPDVtx(0), fhPtPhotonNPileUpTrkVtx(0),
108fhPtPhotonNPileUpSPDVtxTimeCut(0), fhPtPhotonNPileUpTrkVtxTimeCut(0),
6df33fcb 109fhPtPhotonNPileUpSPDVtxTimeCut2(0), fhPtPhotonNPileUpTrkVtxTimeCut2(0),
bc41680b 110
6df33fcb 111fhEClusterSM(0), fhEPhotonSM(0),
112fhPtClusterSM(0), fhPtPhotonSM(0)
85c4406e 113{
1c5acb87 114 //default ctor
115
4bfeae64 116 for(Int_t i = 0; i < 14; i++)
117 {
4c8f7c2e 118 fhMCPt [i] = 0;
119 fhMCE [i] = 0;
120 fhMCPhi [i] = 0;
121 fhMCEta [i] = 0;
85c4406e 122 fhMCDeltaE [i] = 0;
4c8f7c2e 123 fhMCDeltaPt[i] = 0;
85c4406e 124 fhMC2E [i] = 0;
4c8f7c2e 125 fhMC2Pt [i] = 0;
521636d2 126 }
127
4bfeae64 128 for(Int_t i = 0; i < 7; i++)
129 {
3d5d5078 130 fhPtPrimMC [i] = 0;
131 fhEPrimMC [i] = 0;
132 fhPhiPrimMC[i] = 0;
4cf13296 133 fhEtaPrimMC[i] = 0;
3d5d5078 134 fhYPrimMC [i] = 0;
135
136 fhPtPrimMCAcc [i] = 0;
137 fhEPrimMCAcc [i] = 0;
138 fhPhiPrimMCAcc[i] = 0;
4cf13296 139 fhEtaPrimMCAcc[i] = 0;
3d5d5078 140 fhYPrimMCAcc [i] = 0;
d2655d46 141
142 fhDispEtaDispPhi[i] = 0;
143 fhLambda0DispPhi[i] = 0;
144 fhLambda0DispEta[i] = 0;
bc41680b 145
146 fhPtPhotonPileUp[i] = 0;
fad96885 147 fhClusterTimeDiffPhotonPileUp [i] = 0;
bc41680b 148
d2655d46 149 for(Int_t j = 0; j < 6; j++)
150 {
151 fhMCDispEtaDispPhi[i][j] = 0;
152 fhMCLambda0DispEta[i][j] = 0;
153 fhMCLambda0DispPhi[i][j] = 0;
154 }
85c4406e 155 }
3d5d5078 156
4bfeae64 157 for(Int_t i = 0; i < 6; i++)
158 {
f66d95af 159 fhMCELambda0 [i] = 0;
160 fhMCELambda1 [i] = 0;
161 fhMCEDispersion [i] = 0;
85c4406e 162 fhMCNCellsE [i] = 0;
163 fhMCMaxCellDiffClusterE[i] = 0;
bfdcf7fb 164 fhLambda0DispEta[i] = 0;
165 fhLambda0DispPhi[i] = 0;
85c4406e 166
f66d95af 167 fhMCLambda0vsClusterMaxCellDiffE0[i] = 0;
168 fhMCLambda0vsClusterMaxCellDiffE2[i] = 0;
169 fhMCLambda0vsClusterMaxCellDiffE6[i] = 0;
170 fhMCNCellsvsClusterMaxCellDiffE0 [i] = 0;
171 fhMCNCellsvsClusterMaxCellDiffE2 [i] = 0;
172 fhMCNCellsvsClusterMaxCellDiffE6 [i] = 0;
34c16486 173
174 fhMCEDispEta [i] = 0;
175 fhMCEDispPhi [i] = 0;
176 fhMCESumEtaPhi [i] = 0;
177 fhMCEDispEtaPhiDiff[i] = 0;
178 fhMCESphericity [i] = 0;
521636d2 179 }
180
85c4406e 181 for(Int_t i = 0; i < 5; i++)
182 {
58ea8ce5 183 fhClusterCutsE [i] = 0;
184 fhClusterCutsPt[i] = 0;
85c4406e 185 }
186
187 // Track matching residuals
188 for(Int_t i = 0; i < 2; i++)
189 {
126b8c62 190 fhTrackMatchedDEta [i] = 0; fhTrackMatchedDPhi [i] = 0; fhTrackMatchedDEtaDPhi [i] = 0;
191 fhTrackMatchedDEtaNeg[i] = 0; fhTrackMatchedDPhiNeg[i] = 0; fhTrackMatchedDEtaDPhiNeg[i] = 0;
192 fhTrackMatchedDEtaPos[i] = 0; fhTrackMatchedDPhiPos[i] = 0; fhTrackMatchedDEtaDPhiPos[i] = 0;
85c4406e 193 fhTrackMatchedDEtaTRD[i] = 0; fhTrackMatchedDPhiTRD[i] = 0;
194 fhTrackMatchedDEtaMCOverlap[i] = 0; fhTrackMatchedDPhiMCOverlap[i] = 0;
195 fhTrackMatchedDEtaMCNoOverlap[i] = 0; fhTrackMatchedDPhiMCNoOverlap[i] = 0;
196 fhTrackMatchedDEtaMCConversion[i] = 0; fhTrackMatchedDPhiMCConversion[i] = 0;
197 fhTrackMatchedMCParticle[i] = 0; fhTrackMatchedMCParticle[i] = 0;
198 fhdEdx[i] = 0; fhEOverP[i] = 0;
199 fhEOverPTRD[i] = 0;
200 }
201
1c5acb87 202 //Initialize parameters
203 InitParameters();
85c4406e 204
1c5acb87 205}
206
b2e375c7 207//_________________________________________________________________________________________
85c4406e 208Bool_t AliAnaPhoton::ClusterSelected(AliVCluster* calo, TLorentzVector mom, Int_t nMaxima)
c4a7d28a 209{
210 //Select clusters if they pass different cuts
fad96885 211
c2a62a94 212 Float_t ptcluster = mom.Pt();
213 Float_t ecluster = mom.E();
c2a62a94 214 Float_t etacluster = mom.Eta();
215 Float_t phicluster = mom.Phi();
58ea8ce5 216
6df33fcb 217 if(phicluster < 0) phicluster+=TMath::TwoPi();
c2a62a94 218
b2e375c7 219 Bool_t matched = IsTrackMatched(calo,GetReader()->GetInputEvent());
fad96885 220
221 if(GetDebug() > 2)
afb3af8a 222 printf("AliAnaPhoton::ClusterSelected() - Current Event %d; Before selection : E %2.2f, pT %2.2f, phi %2.2f, eta %2.2f\n",
c4a7d28a 223 GetReader()->GetEventNumber(),
c2a62a94 224 ecluster,ptcluster, phicluster*TMath::RadToDeg(),etacluster);
fad96885 225
58ea8ce5 226 fhClusterCutsE [1]->Fill( ecluster);
227 fhClusterCutsPt[1]->Fill(ptcluster);
c4a7d28a 228
c2a62a94 229 if(ecluster > 0.5) fhEtaPhi->Fill(etacluster, phicluster);
230
6df33fcb 231 Int_t nSM = GetModuleNumber(calo);
232 if(nSM < GetCaloUtils()->GetNumberOfSuperModulesUsed() && nSM >=0)
233 {
234 fhEClusterSM ->Fill(ecluster ,nSM);
235 fhPtClusterSM->Fill(ptcluster,nSM);
236 }
237
c4a7d28a 238 //.......................................
239 //If too small or big energy, skip it
afb3af8a 240 if(ecluster < GetMinEnergy() || ecluster > GetMaxEnergy() ) return kFALSE ;
09273901 241
c4a7d28a 242 if(GetDebug() > 2) printf("\t Cluster %d Pass E Cut \n",calo->GetID());
09273901 243
58ea8ce5 244 fhClusterCutsE [2]->Fill( ecluster);
245 fhClusterCutsPt[2]->Fill(ptcluster);
85c4406e 246
c4a7d28a 247 //.......................................
248 // TOF cut, BE CAREFUL WITH THIS CUT
249 Double_t tof = calo->GetTOF()*1e9;
250 if(tof < fTimeCutMin || tof > fTimeCutMax) return kFALSE;
09273901 251
c4a7d28a 252 if(GetDebug() > 2) printf("\t Cluster %d Pass Time Cut \n",calo->GetID());
09273901 253
58ea8ce5 254 fhClusterCutsE [3]->Fill( ecluster);
255 fhClusterCutsPt[3]->Fill(ptcluster);
85c4406e 256
c4a7d28a 257 //.......................................
258 if(calo->GetNCells() <= fNCellsCut && GetReader()->GetDataType() != AliCaloTrackReader::kMC) return kFALSE;
09273901 259
c4a7d28a 260 if(GetDebug() > 2) printf("\t Cluster %d Pass NCell Cut \n",calo->GetID());
09273901 261
58ea8ce5 262 fhClusterCutsE [4]->Fill( ecluster);
263 fhClusterCutsPt[4]->Fill(ptcluster);
85c4406e 264
9e51e29a 265 if(nMaxima < fNLMCutMin || nMaxima > fNLMCutMax) return kFALSE ;
266 if(GetDebug() > 2) printf(" \t Cluster %d pass NLM %d of out of range \n",calo->GetID(), nMaxima);
85c4406e 267
58ea8ce5 268 fhClusterCutsE [5]->Fill( ecluster);
269 fhClusterCutsPt[5]->Fill(ptcluster);
85c4406e 270
c4a7d28a 271 //.......................................
272 //Check acceptance selection
34c16486 273 if(IsFiducialCutOn())
274 {
c4a7d28a 275 Bool_t in = GetFiducialCut()->IsInFiducialCut(mom,fCalorimeter) ;
276 if(! in ) return kFALSE ;
277 }
09273901 278
6df33fcb 279 if(GetDebug() > 2) printf("\t Fiducial cut passed \n");
09273901 280
58ea8ce5 281 fhClusterCutsE [6]->Fill( ecluster);
282 fhClusterCutsPt[6]->Fill(ptcluster);
85c4406e 283
c4a7d28a 284 //.......................................
285 //Skip matched clusters with tracks
09273901 286
4bfeae64 287 // Fill matching residual histograms before PID cuts
288 if(fFillTMHisto) FillTrackMatchingResidualHistograms(calo,0);
09273901 289
34c16486 290 if(fRejectTrackMatch)
291 {
b2e375c7 292 if(matched)
34c16486 293 {
c4a7d28a 294 if(GetDebug() > 2) printf("\t Reject track-matched clusters\n");
295 return kFALSE ;
296 }
85c4406e 297 else
c4a7d28a 298 if(GetDebug() > 2) printf(" Track-matching cut passed \n");
299 }// reject matched clusters
09273901 300
58ea8ce5 301 fhClusterCutsE [7]->Fill( ecluster);
302 fhClusterCutsPt[7]->Fill(ptcluster);
85c4406e 303
c4a7d28a 304 //.......................................
305 //Check Distance to Bad channel, set bit.
306 Double_t distBad=calo->GetDistanceToBadChannel() ; //Distance to bad channel
307 if(distBad < 0.) distBad=9999. ; //workout strange convension dist = -1. ;
85c4406e 308 if(distBad < fMinDist)
34c16486 309 {//In bad channel (PHOS cristal size 2.2x2.2 cm), EMCAL ( cell units )
c4a7d28a 310 return kFALSE ;
311 }
312 else if(GetDebug() > 2) printf("\t Bad channel cut passed %4.2f > %2.2f \n",distBad, fMinDist);
fc195fd0 313
58ea8ce5 314 fhClusterCutsE [8]->Fill( ecluster);
315 fhClusterCutsPt[8]->Fill(ptcluster);
09273901 316
85c4406e 317 if(GetDebug() > 0)
fad96885 318 printf("AliAnaPhoton::ClusterSelected() Current Event %d; After selection : E %2.2f, pT %2.2f, phi %2.2f, eta %2.2f\n",
85c4406e 319 GetReader()->GetEventNumber(),
fad96885 320 ecluster, ptcluster,mom.Phi()*TMath::RadToDeg(),mom.Eta());
c4a7d28a 321
322 //All checks passed, cluster selected
323 return kTRUE;
85c4406e 324
c4a7d28a 325}
326
34c16486 327//___________________________________________
328void AliAnaPhoton::FillAcceptanceHistograms()
329{
3d5d5078 330 //Fill acceptance histograms if MC data is available
331
34c16486 332 Double_t photonY = -100 ;
333 Double_t photonE = -1 ;
334 Double_t photonPt = -1 ;
335 Double_t photonPhi = 100 ;
336 Double_t photonEta = -1 ;
85c4406e 337
34c16486 338 Int_t pdg = 0 ;
339 Int_t tag = 0 ;
f1c9c78f 340 Int_t status = 0 ;
34c16486 341 Int_t mcIndex = 0 ;
f1c9c78f 342 Int_t nprim = 0 ;
343 Bool_t inacceptance = kFALSE ;
85c4406e 344
f1c9c78f 345 TParticle * primStack = 0;
346 AliAODMCParticle * primAOD = 0;
347 TLorentzVector lv;
348
349 // Get the ESD MC particles container
350 AliStack * stack = 0;
351 if( GetReader()->ReadStack() )
85c4406e 352 {
f1c9c78f 353 stack = GetMCStack();
354 if(!stack ) return;
355 nprim = stack->GetNtrack();
356 }
357
358 // Get the AOD MC particles container
359 TClonesArray * mcparticles = 0;
360 if( GetReader()->ReadAODMCParticles() )
34c16486 361 {
f1c9c78f 362 mcparticles = GetReader()->GetAODMCParticles();
363 if( !mcparticles ) return;
364 nprim = mcparticles->GetEntriesFast();
365 }
366
367 for(Int_t i=0 ; i < nprim; i++)
368 {
369 if(GetReader()->AcceptOnlyHIJINGLabels() && !GetReader()->IsHIJINGLabel(i)) continue ;
370
371 if(GetReader()->ReadStack())
34c16486 372 {
f1c9c78f 373 primStack = stack->Particle(i) ;
28a99be3 374 if(!primStack)
375 {
376 printf("AliAnaPhoton::FillAcceptanceHistograms() - ESD primaries pointer not available!!\n");
377 continue;
378 }
379
f1c9c78f 380 pdg = primStack->GetPdgCode();
381 status = primStack->GetStatusCode();
3d5d5078 382
f1c9c78f 383 if(primStack->Energy() == TMath::Abs(primStack->Pz())) continue ; //Protection against floating point exception
384
385 //printf("i %d, %s %d %s %d \n",i, stack->Particle(i)->GetName(), stack->Particle(i)->GetPdgCode(),
386 // prim->GetName(), prim->GetPdgCode());
387
388 //Photon kinematics
389 primStack->Momentum(lv);
390
391 photonY = 0.5*TMath::Log((primStack->Energy()-primStack->Pz())/(primStack->Energy()+primStack->Pz())) ;
392 }
393 else
394 {
395 primAOD = (AliAODMCParticle *) mcparticles->At(i);
28a99be3 396 if(!primAOD)
397 {
398 printf("AliAnaPhoton::FillAcceptanceHistograms() - AOD primaries pointer not available!!\n");
399 continue;
400 }
401
f1c9c78f 402 pdg = primAOD->GetPdgCode();
403 status = primAOD->GetStatus();
404
405 if(primAOD->E() == TMath::Abs(primAOD->Pz())) continue ; //Protection against floating point exception
406
407 //Photon kinematics
408 lv.SetPxPyPzE(primAOD->Px(),primAOD->Py(),primAOD->Pz(),primAOD->E());
409
410 photonY = 0.5*TMath::Log((primAOD->E()-primAOD->Pz())/(primAOD->E()+primAOD->Pz())) ;
411 }
412
413 // Select only photons in the final state
414 if(pdg != 22 ) continue ;
415
416 // If too small or too large pt, skip, same cut as for data analysis
417 photonPt = lv.Pt () ;
418
419 if(photonPt < GetMinPt() || photonPt > GetMaxPt() ) continue ;
420
421 photonE = lv.E () ;
422 photonEta = lv.Eta() ;
423 photonPhi = lv.Phi() ;
424
425 if(photonPhi < 0) photonPhi+=TMath::TwoPi();
426
427 // Check if photons hit desired acceptance
428 inacceptance = kFALSE;
429
430 // Check same fidutial borders as in data analysis on top of real acceptance if real was requested.
431 if( GetFiducialCut()->IsInFiducialCut(lv,fCalorimeter)) inacceptance = kTRUE ;
432
433 // Check if photons hit the Calorimeter acceptance
434 if(IsRealCaloAcceptanceOn() && inacceptance) // defined on base class
435 {
436 if(GetReader()->ReadStack() &&
437 !GetCaloUtils()->IsMCParticleInCalorimeterAcceptance(fCalorimeter, primStack)) inacceptance = kFALSE ;
438 if(GetReader()->ReadAODMCParticles() &&
439 !GetCaloUtils()->IsMCParticleInCalorimeterAcceptance(fCalorimeter, primAOD )) inacceptance = kFALSE ;
440 }
441
442 // Get tag of this particle photon from fragmentation, decay, prompt ...
443 // Set the origin of the photon.
444 tag = GetMCAnalysisUtils()->CheckOrigin(i,GetReader());
445
446 if(!GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPhoton))
447 {
448 // A conversion photon from a hadron, skip this kind of photon
449 // printf("AliAnaPhoton::FillAcceptanceHistograms() - not a photon, weird!\n ");
450 // GetMCAnalysisUtils()->PrintMCTag(tag);
451
452 continue;
453 }
454
455 // Consider only final state particles, but this depends on generator,
456 // status 1 is the usual one, in case of not being ok, leave the possibility
457 // to not consider this.
458 if(status > 1) continue ; // Avoid "partonic" photons
459
460 Bool_t takeIt = kFALSE ;
461 if(status == 1 && GetMCAnalysisUtils()->GetMCGenerator()!="" ) takeIt = kTRUE ;
462
760b98f5 463 if (GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCConversion)) continue;
f1c9c78f 464
465 //Origin of photon
760b98f5 466 if (GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPrompt))
f1c9c78f 467 {
468 mcIndex = kmcPPrompt;
469 }
470 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCFragmentation))
471 {
472 mcIndex = kmcPFragmentation ;
473 }
474 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCISR))
475 {
476 mcIndex = kmcPISR;
477 }
478 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0Decay))
479 {
480 mcIndex = kmcPPi0Decay;
481 }
482 else if( (GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCEtaDecay) ||
483 GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCOtherDecay)))
484 {
485 mcIndex = kmcPOtherDecay;
486 }
487 else
488 {
489 // Other decay but from non final state particle
490 mcIndex = kmcPOtherDecay;
491 }//Other origin
492
493 if(!takeIt && (mcIndex == kmcPPi0Decay || mcIndex == kmcPOtherDecay)) takeIt = kTRUE ;
494
495 if(!takeIt) continue ;
3d5d5078 496
f1c9c78f 497 //Fill histograms
498 fhYPrimMC[kmcPPhoton]->Fill(photonPt, photonY) ;
499 if(TMath::Abs(photonY) < 1.0)
500 {
501 fhEPrimMC [kmcPPhoton]->Fill(photonE ) ;
502 fhPtPrimMC [kmcPPhoton]->Fill(photonPt) ;
503 fhPhiPrimMC[kmcPPhoton]->Fill(photonE , photonPhi) ;
504 fhEtaPrimMC[kmcPPhoton]->Fill(photonE , photonEta) ;
505 }
506 if(inacceptance)
507 {
508 fhEPrimMCAcc [kmcPPhoton]->Fill(photonE ) ;
509 fhPtPrimMCAcc [kmcPPhoton]->Fill(photonPt) ;
510 fhPhiPrimMCAcc[kmcPPhoton]->Fill(photonE , photonPhi) ;
511 fhEtaPrimMCAcc[kmcPPhoton]->Fill(photonE , photonEta) ;
512 fhYPrimMCAcc [kmcPPhoton]->Fill(photonE , photonY) ;
513 }//Accepted
514
760b98f5 515 if(mcIndex < fNPrimaryHistograms)
f1c9c78f 516 {
760b98f5 517 fhYPrimMC[mcIndex]->Fill(photonPt, photonY) ;
518 if(TMath::Abs(photonY) < 1.0)
519 {
520 fhEPrimMC [mcIndex]->Fill(photonE ) ;
521 fhPtPrimMC [mcIndex]->Fill(photonPt) ;
522 fhPhiPrimMC[mcIndex]->Fill(photonE , photonPhi) ;
523 fhEtaPrimMC[mcIndex]->Fill(photonE , photonEta) ;
524 }
525
526 if(inacceptance)
527 {
528 fhEPrimMCAcc [mcIndex]->Fill(photonE ) ;
529 fhPtPrimMCAcc [mcIndex]->Fill(photonPt) ;
530 fhPhiPrimMCAcc[mcIndex]->Fill(photonE , photonPhi) ;
531 fhEtaPrimMCAcc[mcIndex]->Fill(photonE , photonEta) ;
532 fhYPrimMCAcc [mcIndex]->Fill(photonE , photonY) ;
533 }//Accepted
f1c9c78f 534 }
f1c9c78f 535 }//loop on primaries
536
3d5d5078 537}
521636d2 538
bc41680b 539//________________________________________________________________________________
540void AliAnaPhoton::FillPileUpHistograms(AliVCluster* cluster, AliVCaloCells *cells)
b2e375c7 541{
bc41680b 542 // Fill some histograms to understand pile-up
b2e375c7 543
bc41680b 544 TLorentzVector mom;
545 cluster->GetMomentum(mom,GetVertex(0));
546 Float_t pt = mom.Pt();
547 Float_t time = cluster->GetTOF()*1.e9;
b2e375c7 548
bc41680b 549 AliVEvent * event = GetReader()->GetInputEvent();
b2e375c7 550
bc41680b 551 if(GetReader()->IsPileUpFromSPD()) fhPtPhotonPileUp[0]->Fill(pt);
552 if(GetReader()->IsPileUpFromEMCal()) fhPtPhotonPileUp[1]->Fill(pt);
553 if(GetReader()->IsPileUpFromSPDOrEMCal()) fhPtPhotonPileUp[2]->Fill(pt);
554 if(GetReader()->IsPileUpFromSPDAndEMCal()) fhPtPhotonPileUp[3]->Fill(pt);
555 if(GetReader()->IsPileUpFromSPDAndNotEMCal()) fhPtPhotonPileUp[4]->Fill(pt);
556 if(GetReader()->IsPileUpFromEMCalAndNotSPD()) fhPtPhotonPileUp[5]->Fill(pt);
557 if(GetReader()->IsPileUpFromNotSPDAndNotEMCal()) fhPtPhotonPileUp[6]->Fill(pt);
b2e375c7 558
bc41680b 559 fhTimePtPhotonNoCut->Fill(pt,time);
560 if(GetReader()->IsPileUpFromSPD()) fhTimePtPhotonSPD->Fill(pt,time);
b2e375c7 561
bc41680b 562 // cells inside the cluster
b2e375c7 563 Float_t maxCellFraction = 0.;
bc41680b 564 Int_t absIdMax = GetCaloUtils()->GetMaxEnergyCell( cells, cluster, maxCellFraction);
b2e375c7 565
126b8c62 566 //Loop on cells inside cluster, max cell must be over 100 MeV and time in BC=0
bc41680b 567 if(cells->GetCellAmplitude(absIdMax) > 0.1 && TMath::Abs(time) < 30)
b2e375c7 568 {
bc41680b 569 for (Int_t ipos = 0; ipos < cluster->GetNCells(); ipos++)
b2e375c7 570 {
bc41680b 571 Int_t absId = cluster->GetCellsAbsId()[ipos];
126b8c62 572
573 if( absId == absIdMax ) continue ;
574
bc41680b 575 Double_t tcell = cells->GetCellTime(absId);
b2e375c7 576 Float_t amp = cells->GetCellAmplitude(absId);
577 Int_t bc = GetReader()->GetInputEvent()->GetBunchCrossNumber();
b2e375c7 578
bc41680b 579 GetCaloUtils()->GetEMCALRecoUtils()->AcceptCalibrateCell(absId,bc,amp,tcell,cells);
580 tcell*=1e9;
b2e375c7 581
bc41680b 582 Float_t diff = (time-tcell);
b2e375c7 583
36769d30 584 if( cells->GetCellAmplitude(absIdMax) < 0.1 ) continue ;
126b8c62 585
bc41680b 586 if(GetReader()->IsPileUpFromSPD()) fhClusterTimeDiffPhotonPileUp[0]->Fill(pt, diff);
587 if(GetReader()->IsPileUpFromEMCal()) fhClusterTimeDiffPhotonPileUp[1]->Fill(pt, diff);
588 if(GetReader()->IsPileUpFromSPDOrEMCal()) fhClusterTimeDiffPhotonPileUp[2]->Fill(pt, diff);
589 if(GetReader()->IsPileUpFromSPDAndEMCal()) fhClusterTimeDiffPhotonPileUp[3]->Fill(pt, diff);
590 if(GetReader()->IsPileUpFromSPDAndNotEMCal()) fhClusterTimeDiffPhotonPileUp[4]->Fill(pt, diff);
591 if(GetReader()->IsPileUpFromEMCalAndNotSPD()) fhClusterTimeDiffPhotonPileUp[5]->Fill(pt, diff);
592 if(GetReader()->IsPileUpFromNotSPDAndNotEMCal()) fhClusterTimeDiffPhotonPileUp[6]->Fill(pt, diff);
b2e375c7 593
bc41680b 594 }//loop
85c4406e 595 }
acd56ca4 596
2ad19c3d 597 AliESDEvent* esdEv = dynamic_cast<AliESDEvent*> (event);
598 AliAODEvent* aodEv = dynamic_cast<AliAODEvent*> (event);
599
600 // N pile up vertices
0f7e7205 601 Int_t nVtxSPD = -1;
602 Int_t nVtxTrk = -1;
2ad19c3d 603
604 if (esdEv)
605 {
0f7e7205 606 nVtxSPD = esdEv->GetNumberOfPileupVerticesSPD();
607 nVtxTrk = esdEv->GetNumberOfPileupVerticesTracks();
85c4406e 608
2ad19c3d 609 }//ESD
610 else if (aodEv)
611 {
0f7e7205 612 nVtxSPD = aodEv->GetNumberOfPileupVerticesSPD();
613 nVtxTrk = aodEv->GetNumberOfPileupVerticesTracks();
2ad19c3d 614 }//AOD
615
bc41680b 616 if(pt < 8)
617 {
618 fhTimeNPileUpVertSPD ->Fill(time,nVtxSPD);
619 fhTimeNPileUpVertTrack->Fill(time,nVtxTrk);
620 }
2ad19c3d 621
85c4406e 622 fhPtPhotonNPileUpSPDVtx->Fill(pt,nVtxSPD);
0f7e7205 623 fhPtPhotonNPileUpTrkVtx->Fill(pt,nVtxTrk);
624
625 if(TMath::Abs(time) < 25)
85c4406e 626 {
627 fhPtPhotonNPileUpSPDVtxTimeCut->Fill(pt,nVtxSPD);
628 fhPtPhotonNPileUpTrkVtxTimeCut->Fill(pt,nVtxTrk);
0f7e7205 629 }
630
85c4406e 631 if(time < 75 && time > -25)
632 {
633 fhPtPhotonNPileUpSPDVtxTimeCut2->Fill(pt,nVtxSPD);
634 fhPtPhotonNPileUpTrkVtxTimeCut2->Fill(pt,nVtxTrk);
635 }
636
2ad19c3d 637}
638
34c16486 639//____________________________________________________________________________________
22ad7981 640void AliAnaPhoton::FillShowerShapeHistograms(AliVCluster* cluster, Int_t mcTag)
34c16486 641{
85c4406e 642 //Fill cluster Shower Shape histograms
521636d2 643
644 if(!fFillSSHistograms || GetMixedEvent()) return;
85c4406e 645
521636d2 646 Float_t energy = cluster->E();
647 Int_t ncells = cluster->GetNCells();
521636d2 648 Float_t lambda0 = cluster->GetM02();
649 Float_t lambda1 = cluster->GetM20();
650 Float_t disp = cluster->GetDispersion()*cluster->GetDispersion();
651
652 TLorentzVector mom;
34c16486 653 if(GetReader()->GetDataType() != AliCaloTrackReader::kMC)
654 {
655 cluster->GetMomentum(mom,GetVertex(0)) ;
656 }//Assume that come from vertex in straight line
657 else
658 {
521636d2 659 Double_t vertex[]={0,0,0};
660 cluster->GetMomentum(mom,vertex) ;
661 }
662
663 Float_t eta = mom.Eta();
664 Float_t phi = mom.Phi();
665 if(phi < 0) phi+=TMath::TwoPi();
666
667 fhLam0E ->Fill(energy,lambda0);
668 fhLam1E ->Fill(energy,lambda1);
669 fhDispE ->Fill(energy,disp);
85c4406e 670
4d1d8f00 671 if(fCalorimeter == "EMCAL" && GetFirstSMCoveredByTRD() >= 0 &&
672 GetModuleNumber(cluster) >= GetFirstSMCoveredByTRD() )
34c16486 673 {
521636d2 674 fhLam0ETRD->Fill(energy,lambda0);
675 fhLam1ETRD->Fill(energy,lambda1);
676 fhDispETRD->Fill(energy,disp);
521636d2 677 }
678
34c16486 679 Float_t l0 = 0., l1 = 0.;
85c4406e 680 Float_t dispp= 0., dEta = 0., dPhi = 0.;
681 Float_t sEta = 0., sPhi = 0., sEtaPhi = 0.;
764ab1f4 682 if(fCalorimeter == "EMCAL" && !fFillOnlySimpleSSHisto)
34c16486 683 {
684 GetCaloUtils()->GetEMCALRecoUtils()->RecalculateClusterShowerShapeParameters(GetEMCALGeometry(), GetReader()->GetInputEvent()->GetEMCALCells(), cluster,
685 l0, l1, dispp, dEta, dPhi, sEta, sPhi, sEtaPhi);
686 //printf("AliAnaPhoton::FillShowerShapeHistogram - l0 %2.6f, l1 %2.6f, disp %2.6f, dEta %2.6f, dPhi %2.6f, sEta %2.6f, sPhi %2.6f, sEtaPhi %2.6f \n",
687 // l0, l1, dispp, dEta, dPhi, sEta, sPhi, sEtaPhi );
688 //printf("AliAnaPhoton::FillShowerShapeHistogram - dispersion %f, dispersion eta+phi %f \n",
689 // disp, dPhi+dEta );
690 fhDispEtaE -> Fill(energy,dEta);
691 fhDispPhiE -> Fill(energy,dPhi);
692 fhSumEtaE -> Fill(energy,sEta);
693 fhSumPhiE -> Fill(energy,sPhi);
694 fhSumEtaPhiE -> Fill(energy,sEtaPhi);
695 fhDispEtaPhiDiffE -> Fill(energy,dPhi-dEta);
696 if(dEta+dPhi>0)fhSphericityE -> Fill(energy,(dPhi-dEta)/(dEta+dPhi));
697 if(dEta+sEta>0)fhDispSumEtaDiffE -> Fill(energy,(dEta-sEta)/((dEta+sEta)/2.));
85c4406e 698 if(dPhi+sPhi>0)fhDispSumPhiDiffE -> Fill(energy,(dPhi-sPhi)/((dPhi+sPhi)/2.));
34c16486 699
bfdcf7fb 700 Int_t ebin = -1;
701 if (energy < 2 ) ebin = 0;
702 else if (energy < 4 ) ebin = 1;
703 else if (energy < 6 ) ebin = 2;
704 else if (energy < 10) ebin = 3;
85c4406e 705 else if (energy < 15) ebin = 4;
706 else if (energy < 20) ebin = 5;
707 else ebin = 6;
bfdcf7fb 708
709 fhDispEtaDispPhi[ebin]->Fill(dEta ,dPhi);
710 fhLambda0DispEta[ebin]->Fill(lambda0,dEta);
711 fhLambda0DispPhi[ebin]->Fill(lambda0,dPhi);
712
34c16486 713 }
714
85c4406e 715 // if track-matching was of, check effect of track-matching residual cut
b5dbb99b 716
717 if(!fRejectTrackMatch)
718 {
719 Float_t dZ = cluster->GetTrackDz();
720 Float_t dR = cluster->GetTrackDx();
34c16486 721 if(cluster->IsEMCAL() && GetCaloUtils()->IsRecalculationOfClusterTrackMatchingOn())
722 {
b5dbb99b 723 dR = 2000., dZ = 2000.;
724 GetCaloUtils()->GetEMCALRecoUtils()->GetMatchedResiduals(cluster->GetID(),dZ,dR);
85c4406e 725 }
b5dbb99b 726
727 if(TMath::Abs(dZ) < 0.05 && TMath::Abs(dR) < 0.05)
728 {
729 fhLam0ETM ->Fill(energy,lambda0);
730 fhLam1ETM ->Fill(energy,lambda1);
731 fhDispETM ->Fill(energy,disp);
732
4d1d8f00 733 if(fCalorimeter == "EMCAL" && GetFirstSMCoveredByTRD() >= 0 &&
734 GetModuleNumber(cluster) >= GetFirstSMCoveredByTRD() )
34c16486 735 {
b5dbb99b 736 fhLam0ETMTRD->Fill(energy,lambda0);
737 fhLam1ETMTRD->Fill(energy,lambda1);
738 fhDispETMTRD->Fill(energy,disp);
739 }
740 }
85c4406e 741 }// if track-matching was of, check effect of matching residual cut
b5dbb99b 742
764ab1f4 743
4301207e 744 if(!fFillOnlySimpleSSHisto)
745 {
764ab1f4 746 if(energy < 2)
747 {
748 fhNCellsLam0LowE ->Fill(ncells,lambda0);
749 fhNCellsLam1LowE ->Fill(ncells,lambda1);
750 fhNCellsDispLowE ->Fill(ncells,disp);
751
752 fhLam1Lam0LowE ->Fill(lambda1,lambda0);
753 fhLam0DispLowE ->Fill(lambda0,disp);
754 fhDispLam1LowE ->Fill(disp,lambda1);
755 fhEtaLam0LowE ->Fill(eta,lambda0);
85c4406e 756 fhPhiLam0LowE ->Fill(phi,lambda0);
764ab1f4 757 }
85c4406e 758 else
764ab1f4 759 {
760 fhNCellsLam0HighE ->Fill(ncells,lambda0);
761 fhNCellsLam1HighE ->Fill(ncells,lambda1);
762 fhNCellsDispHighE ->Fill(ncells,disp);
763
764 fhLam1Lam0HighE ->Fill(lambda1,lambda0);
765 fhLam0DispHighE ->Fill(lambda0,disp);
766 fhDispLam1HighE ->Fill(disp,lambda1);
767 fhEtaLam0HighE ->Fill(eta, lambda0);
768 fhPhiLam0HighE ->Fill(phi, lambda0);
769 }
521636d2 770 }
4301207e 771
34c16486 772 if(IsDataMC())
773 {
f66d95af 774 AliVCaloCells* cells = 0;
775 if(fCalorimeter == "EMCAL") cells = GetEMCALCells();
776 else cells = GetPHOSCells();
3d5d5078 777
778 //Fill histograms to check shape of embedded clusters
779 Float_t fraction = 0;
85c4406e 780 // printf("check embedding %i\n",GetReader()->IsEmbeddedClusterSelectionOn());
781
34c16486 782 if(GetReader()->IsEmbeddedClusterSelectionOn())
783 {//Only working for EMCAL
85c4406e 784 // printf("embedded\n");
3d5d5078 785 Float_t clusterE = 0; // recalculate in case corrections applied.
786 Float_t cellE = 0;
34c16486 787 for(Int_t icell = 0; icell < cluster->GetNCells(); icell++)
788 {
3d5d5078 789 cellE = cells->GetCellAmplitude(cluster->GetCellAbsId(icell));
85c4406e 790 clusterE+=cellE;
3d5d5078 791 fraction+=cellE*cluster->GetCellAmplitudeFraction(icell);
792 }
793
794 //Fraction of total energy due to the embedded signal
795 fraction/=clusterE;
796
85c4406e 797 if(GetDebug() > 1 )
8d6b7f60 798 printf("AliAnaPhoton::FillShowerShapeHistogram() - Energy fraction of embedded signal %2.3f, Energy %2.3f\n",fraction, clusterE);
3d5d5078 799
800 fhEmbeddedSignalFractionEnergy->Fill(clusterE,fraction);
801
85c4406e 802 } // embedded fraction
3d5d5078 803
f66d95af 804 // Get the fraction of the cluster energy that carries the cell with highest energy
f66d95af 805 Float_t maxCellFraction = 0.;
4301207e 806 Int_t absID = GetCaloUtils()->GetMaxEnergyCell(cells, cluster,maxCellFraction);
f66d95af 807
4301207e 808 if( absID < 0 ) AliFatal("Wrong absID");
809
f66d95af 810 // Check the origin and fill histograms
34c16486 811
812 Int_t mcIndex = -1;
813
85c4406e 814 if( GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPhoton) &&
3d5d5078 815 !GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCConversion) &&
34c16486 816 !GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPi0) &&
817 !GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCEta))
818 {
819 mcIndex = kmcssPhoton ;
85c4406e 820
34c16486 821 if(!GetReader()->IsEmbeddedClusterSelectionOn())
822 {
3d5d5078 823 //Check particle overlaps in cluster
824
85c4406e 825 // Compare the primary depositing more energy with the rest,
8d6b7f60 826 // if no photon/electron as comon ancestor (conversions), count as other particle
4914e781 827 const UInt_t nlabels = cluster->GetNLabels();
828 Int_t overpdg[nlabels];
829 Int_t noverlaps = GetMCAnalysisUtils()->GetNOverlaps(cluster->GetLabels(), nlabels,mcTag,-1,GetReader(),overpdg);
830
8d6b7f60 831 //printf("N overlaps %d \n",noverlaps);
3d5d5078 832
f27fe026 833 if(noverlaps == 0)
34c16486 834 {
3d5d5078 835 fhMCPhotonELambda0NoOverlap ->Fill(energy, lambda0);
3d5d5078 836 }
f27fe026 837 else if(noverlaps == 1)
85c4406e 838 {
3d5d5078 839 fhMCPhotonELambda0TwoOverlap ->Fill(energy, lambda0);
3d5d5078 840 }
f27fe026 841 else if(noverlaps > 1)
85c4406e 842 {
3d5d5078 843 fhMCPhotonELambda0NOverlap ->Fill(energy, lambda0);
3d5d5078 844 }
85c4406e 845 else
34c16486 846 {
f27fe026 847 printf("AliAnaPhoton::FillShowerShapeHistogram() - n overlaps = %d!!\n", noverlaps);
3d5d5078 848 }
849 }//No embedding
850
851 //Fill histograms to check shape of embedded clusters
34c16486 852 if(GetReader()->IsEmbeddedClusterSelectionOn())
853 {
85c4406e 854 if (fraction > 0.9)
3d5d5078 855 {
856 fhEmbedPhotonELambda0FullSignal ->Fill(energy, lambda0);
3d5d5078 857 }
858 else if(fraction > 0.5)
859 {
860 fhEmbedPhotonELambda0MostlySignal ->Fill(energy, lambda0);
3d5d5078 861 }
862 else if(fraction > 0.1)
85c4406e 863 {
3d5d5078 864 fhEmbedPhotonELambda0MostlyBkg ->Fill(energy, lambda0);
3d5d5078 865 }
866 else
867 {
868 fhEmbedPhotonELambda0FullBkg ->Fill(energy, lambda0);
3d5d5078 869 }
870 } // embedded
871
521636d2 872 }//photon no conversion
4431f13a 873 else if ( GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPhoton) &&
874 GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCConversion) &&
875 !GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPi0) &&
876 !GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCEta))
877 {
878 mcIndex = kmcssConversion ;
879 }//conversion photon
880
34c16486 881 else if ( GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCElectron))
882 {
883 mcIndex = kmcssElectron ;
521636d2 884 }//electron
34c16486 885 else if ( GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPi0) )
886 {
887 mcIndex = kmcssPi0 ;
3d5d5078 888
889 //Fill histograms to check shape of embedded clusters
34c16486 890 if(GetReader()->IsEmbeddedClusterSelectionOn())
891 {
85c4406e 892 if (fraction > 0.9)
3d5d5078 893 {
894 fhEmbedPi0ELambda0FullSignal ->Fill(energy, lambda0);
3d5d5078 895 }
896 else if(fraction > 0.5)
897 {
898 fhEmbedPi0ELambda0MostlySignal ->Fill(energy, lambda0);
3d5d5078 899 }
900 else if(fraction > 0.1)
85c4406e 901 {
3d5d5078 902 fhEmbedPi0ELambda0MostlyBkg ->Fill(energy, lambda0);
3d5d5078 903 }
904 else
905 {
906 fhEmbedPi0ELambda0FullBkg ->Fill(energy, lambda0);
3d5d5078 907 }
85c4406e 908 } // embedded
3d5d5078 909
521636d2 910 }//pi0
34c16486 911 else if ( GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCEta) )
912 {
913 mcIndex = kmcssEta ;
85c4406e 914 }//eta
915 else
34c16486 916 {
85c4406e 917 mcIndex = kmcssOther ;
918 }//other particles
521636d2 919
34c16486 920 fhMCELambda0 [mcIndex]->Fill(energy, lambda0);
921 fhMCELambda1 [mcIndex]->Fill(energy, lambda1);
922 fhMCEDispersion [mcIndex]->Fill(energy, disp);
923 fhMCNCellsE [mcIndex]->Fill(energy, ncells);
924 fhMCMaxCellDiffClusterE[mcIndex]->Fill(energy, maxCellFraction);
925
764ab1f4 926 if(!fFillOnlySimpleSSHisto)
34c16486 927 {
764ab1f4 928 if (energy < 2.)
929 {
930 fhMCLambda0vsClusterMaxCellDiffE0[mcIndex]->Fill(lambda0, maxCellFraction);
931 fhMCNCellsvsClusterMaxCellDiffE0 [mcIndex]->Fill(ncells, maxCellFraction);
932 }
933 else if(energy < 6.)
934 {
935 fhMCLambda0vsClusterMaxCellDiffE2[mcIndex]->Fill(lambda0, maxCellFraction);
936 fhMCNCellsvsClusterMaxCellDiffE2 [mcIndex]->Fill(ncells, maxCellFraction);
937 }
938 else
939 {
940 fhMCLambda0vsClusterMaxCellDiffE6[mcIndex]->Fill(lambda0, maxCellFraction);
941 fhMCNCellsvsClusterMaxCellDiffE6 [mcIndex]->Fill(ncells, maxCellFraction);
942 }
943
944 if(fCalorimeter == "EMCAL")
945 {
946 fhMCEDispEta [mcIndex]-> Fill(energy,dEta);
947 fhMCEDispPhi [mcIndex]-> Fill(energy,dPhi);
948 fhMCESumEtaPhi [mcIndex]-> Fill(energy,sEtaPhi);
949 fhMCEDispEtaPhiDiff [mcIndex]-> Fill(energy,dPhi-dEta);
85c4406e 950 if(dEta+dPhi>0)fhMCESphericity[mcIndex]-> Fill(energy,(dPhi-dEta)/(dEta+dPhi));
764ab1f4 951
952 Int_t ebin = -1;
953 if (energy < 2 ) ebin = 0;
954 else if (energy < 4 ) ebin = 1;
955 else if (energy < 6 ) ebin = 2;
956 else if (energy < 10) ebin = 3;
85c4406e 957 else if (energy < 15) ebin = 4;
958 else if (energy < 20) ebin = 5;
959 else ebin = 6;
764ab1f4 960
961 fhMCDispEtaDispPhi[ebin][mcIndex]->Fill(dEta ,dPhi);
962 fhMCLambda0DispEta[ebin][mcIndex]->Fill(lambda0,dEta);
85c4406e 963 fhMCLambda0DispPhi[ebin][mcIndex]->Fill(lambda0,dPhi);
764ab1f4 964 }
34c16486 965 }
521636d2 966 }//MC data
967
968}
969
4bfeae64 970//__________________________________________________________________________
85c4406e 971void AliAnaPhoton::FillTrackMatchingResidualHistograms(AliVCluster* cluster,
22ad7981 972 Int_t cut)
4bfeae64 973{
974 // If selected, fill histograms with residuals of matched clusters, help to define track matching cut
975 // Residual filled for different cuts 0 (No cut), after 1 PID cut
85c4406e 976
4bfeae64 977 Float_t dZ = cluster->GetTrackDz();
978 Float_t dR = cluster->GetTrackDx();
979
980 if(cluster->IsEMCAL() && GetCaloUtils()->IsRecalculationOfClusterTrackMatchingOn())
981 {
982 dR = 2000., dZ = 2000.;
983 GetCaloUtils()->GetEMCALRecoUtils()->GetMatchedResiduals(cluster->GetID(),dZ,dR);
85c4406e 984 }
985
b2e375c7 986 AliVTrack *track = GetCaloUtils()->GetMatchedTrack(cluster, GetReader()->GetInputEvent());
987
988 Bool_t positive = kFALSE;
989 if(track) positive = (track->Charge()>0);
990
4bfeae64 991 if(fhTrackMatchedDEta[cut] && TMath::Abs(dR) < 999)
992 {
993 fhTrackMatchedDEta[cut]->Fill(cluster->E(),dZ);
994 fhTrackMatchedDPhi[cut]->Fill(cluster->E(),dR);
4bfeae64 995 if(cluster->E() > 0.5) fhTrackMatchedDEtaDPhi[cut]->Fill(dZ,dR);
b2e375c7 996
997 if(track)
998 {
999 if(positive)
1000 {
1001 fhTrackMatchedDEtaPos[cut]->Fill(cluster->E(),dZ);
1002 fhTrackMatchedDPhiPos[cut]->Fill(cluster->E(),dR);
1003 if(cluster->E() > 0.5) fhTrackMatchedDEtaDPhiPos[cut]->Fill(dZ,dR);
1004 }
1005 else
1006 {
1007 fhTrackMatchedDEtaNeg[cut]->Fill(cluster->E(),dZ);
1008 fhTrackMatchedDPhiNeg[cut]->Fill(cluster->E(),dR);
1009 if(cluster->E() > 0.5) fhTrackMatchedDEtaDPhiNeg[cut]->Fill(dZ,dR);
1010 }
1011 }
4bfeae64 1012
1013 Int_t nSMod = GetModuleNumber(cluster);
1014
4d1d8f00 1015 if(fCalorimeter=="EMCAL" && GetFirstSMCoveredByTRD() >= 0 &&
1016 nSMod >= GetFirstSMCoveredByTRD() )
4bfeae64 1017 {
1018 fhTrackMatchedDEtaTRD[cut]->Fill(cluster->E(),dZ);
1019 fhTrackMatchedDPhiTRD[cut]->Fill(cluster->E(),dR);
1020 }
1021
1022 // Check dEdx and E/p of matched clusters
b2e375c7 1023
4bfeae64 1024 if(TMath::Abs(dZ) < 0.05 && TMath::Abs(dR) < 0.05)
85c4406e 1025 {
85c4406e 1026 if(track)
4bfeae64 1027 {
4bfeae64 1028 Float_t dEdx = track->GetTPCsignal();
1029 Float_t eOverp = cluster->E()/track->P();
1030
1031 fhdEdx[cut] ->Fill(cluster->E(), dEdx);
1032 fhEOverP[cut]->Fill(cluster->E(), eOverp);
1033
4d1d8f00 1034 if(fCalorimeter=="EMCAL" && GetFirstSMCoveredByTRD() >= 0 &&
1035 nSMod >= GetFirstSMCoveredByTRD() )
4bfeae64 1036 fhEOverPTRD[cut]->Fill(cluster->E(), eOverp);
1037
1038
1039 }
1040 else
85c4406e 1041 printf("AliAnaPhoton::FillTrackMatchingResidualHistograms() - Residual OK but (dR, dZ)= (%2.4f,%2.4f) no track associated WHAT? \n", dR,dZ);
4bfeae64 1042
1043
1044
1045 if(IsDataMC())
1046 {
1047
2644ead9 1048 Int_t tag = GetMCAnalysisUtils()->CheckOrigin(cluster->GetLabels(),cluster->GetNLabels(),GetReader());
4bfeae64 1049
1050 if ( !GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCConversion) )
1051 {
1052 if ( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0) ||
85c4406e 1053 GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCEta) ) fhTrackMatchedMCParticle[cut]->Fill(cluster->E(), 2.5 );
4bfeae64 1054 else if ( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPhoton) ) fhTrackMatchedMCParticle[cut]->Fill(cluster->E(), 0.5 );
1055 else if ( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCElectron) ) fhTrackMatchedMCParticle[cut]->Fill(cluster->E(), 1.5 );
1056 else fhTrackMatchedMCParticle[cut]->Fill(cluster->E(), 3.5 );
1057
1058 // Check if several particles contributed to cluster and discard overlapped mesons
85c4406e 1059 if(!GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0) ||
34c16486 1060 !GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCEta))
1061 {
4bfeae64 1062 if(cluster->GetNLabels()==1)
1063 {
1064 fhTrackMatchedDEtaMCNoOverlap[cut]->Fill(cluster->E(),dZ);
1065 fhTrackMatchedDPhiMCNoOverlap[cut]->Fill(cluster->E(),dR);
1066 }
85c4406e 1067 else
4bfeae64 1068 {
1069 fhTrackMatchedDEtaMCOverlap[cut]->Fill(cluster->E(),dZ);
1070 fhTrackMatchedDPhiMCOverlap[cut]->Fill(cluster->E(),dR);
1071 }
1072
1073 }// Check overlaps
1074
1075 }
1076 else
1077 {
1078 if ( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0) ||
85c4406e 1079 GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCEta) ) fhTrackMatchedMCParticle[cut]->Fill(cluster->E(), 6.5 );
4bfeae64 1080 else if ( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPhoton) ) fhTrackMatchedMCParticle[cut]->Fill(cluster->E(), 4.5 );
1081 else if ( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCElectron) ) fhTrackMatchedMCParticle[cut]->Fill(cluster->E(), 5.5 );
1082 else fhTrackMatchedMCParticle[cut]->Fill(cluster->E(), 7.5 );
1083
1084 // Check if several particles contributed to cluster
85c4406e 1085 if(!GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0) ||
34c16486 1086 !GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCEta))
1087 {
4bfeae64 1088 fhTrackMatchedDEtaMCConversion[cut]->Fill(cluster->E(),dZ);
1089 fhTrackMatchedDPhiMCConversion[cut]->Fill(cluster->E(),dR);
1090
85c4406e 1091 }// Check overlaps
4bfeae64 1092
1093 }
1094
85c4406e 1095 } // MC
4bfeae64 1096
1097 } // residuals window
1098
1099 } // Small residual
1100
1101}
1102
1103//___________________________________________
0c1383b5 1104TObjString * AliAnaPhoton::GetAnalysisCuts()
85c4406e 1105{
0c1383b5 1106 //Save parameters used for analysis
1107 TString parList ; //this will be list of parameters used for this analysis.
5ae09196 1108 const Int_t buffersize = 255;
1109 char onePar[buffersize] ;
0c1383b5 1110
5ae09196 1111 snprintf(onePar,buffersize,"--- AliAnaPhoton ---\n") ;
85c4406e 1112 parList+=onePar ;
5ae09196 1113 snprintf(onePar,buffersize,"Calorimeter: %s\n",fCalorimeter.Data()) ;
0c1383b5 1114 parList+=onePar ;
5ae09196 1115 snprintf(onePar,buffersize,"fMinDist =%2.2f (Minimal distance to bad channel to accept cluster) \n",fMinDist) ;
0c1383b5 1116 parList+=onePar ;
5ae09196 1117 snprintf(onePar,buffersize,"fMinDist2=%2.2f (Cuts on Minimal distance to study acceptance evaluation) \n",fMinDist2) ;
0c1383b5 1118 parList+=onePar ;
5ae09196 1119 snprintf(onePar,buffersize,"fMinDist3=%2.2f (One more cut on distance used for acceptance-efficiency study) \n",fMinDist3) ;
0c1383b5 1120 parList+=onePar ;
5ae09196 1121 snprintf(onePar,buffersize,"fRejectTrackMatch: %d\n",fRejectTrackMatch) ;
85c4406e 1122 parList+=onePar ;
0c1383b5 1123
1124 //Get parameters set in base class.
1125 parList += GetBaseParametersList() ;
1126
1127 //Get parameters set in PID class.
1128 parList += GetCaloPID()->GetPIDParametersList() ;
1129
1130 //Get parameters set in FiducialCut class (not available yet)
85c4406e 1131 //parlist += GetFidCut()->GetFidCutParametersList()
0c1383b5 1132
1133 return new TObjString(parList) ;
1134}
1135
1c5acb87 1136//________________________________________________________________________
1137TList * AliAnaPhoton::GetCreateOutputObjects()
c2a62a94 1138{
85c4406e 1139 // Create histograms to be saved in output file and
477d6cee 1140 // store them in outputContainer
85c4406e 1141 TList * outputContainer = new TList() ;
1142 outputContainer->SetName("PhotonHistos") ;
4a745797 1143
85c4406e 1144 Int_t nptbins = GetHistogramRanges()->GetHistoPtBins(); Float_t ptmax = GetHistogramRanges()->GetHistoPtMax(); Float_t ptmin = GetHistogramRanges()->GetHistoPtMin();
1145 Int_t nphibins = GetHistogramRanges()->GetHistoPhiBins(); Float_t phimax = GetHistogramRanges()->GetHistoPhiMax(); Float_t phimin = GetHistogramRanges()->GetHistoPhiMin();
1146 Int_t netabins = GetHistogramRanges()->GetHistoEtaBins(); Float_t etamax = GetHistogramRanges()->GetHistoEtaMax(); Float_t etamin = GetHistogramRanges()->GetHistoEtaMin();
745913ae 1147 Int_t ssbins = GetHistogramRanges()->GetHistoShowerShapeBins(); Float_t ssmax = GetHistogramRanges()->GetHistoShowerShapeMax(); Float_t ssmin = GetHistogramRanges()->GetHistoShowerShapeMin();
85c4406e 1148 Int_t nbins = GetHistogramRanges()->GetHistoNClusterCellBins(); Int_t nmax = GetHistogramRanges()->GetHistoNClusterCellMax(); Int_t nmin = GetHistogramRanges()->GetHistoNClusterCellMin();
1149 Int_t ntimebins= GetHistogramRanges()->GetHistoTimeBins(); Float_t timemax = GetHistogramRanges()->GetHistoTimeMax(); Float_t timemin = GetHistogramRanges()->GetHistoTimeMin();
1150
1151 Int_t nresetabins = GetHistogramRanges()->GetHistoTrackResidualEtaBins();
1152 Float_t resetamax = GetHistogramRanges()->GetHistoTrackResidualEtaMax();
09273901 1153 Float_t resetamin = GetHistogramRanges()->GetHistoTrackResidualEtaMin();
85c4406e 1154 Int_t nresphibins = GetHistogramRanges()->GetHistoTrackResidualPhiBins();
1155 Float_t resphimax = GetHistogramRanges()->GetHistoTrackResidualPhiMax();
09273901 1156 Float_t resphimin = GetHistogramRanges()->GetHistoTrackResidualPhiMin();
1157
85c4406e 1158 Int_t ndedxbins = GetHistogramRanges()->GetHistodEdxBins();
1159 Float_t dedxmax = GetHistogramRanges()->GetHistodEdxMax();
31ae6d59 1160 Float_t dedxmin = GetHistogramRanges()->GetHistodEdxMin();
85c4406e 1161 Int_t nPoverEbins = GetHistogramRanges()->GetHistoPOverEBins();
1162 Float_t pOverEmax = GetHistogramRanges()->GetHistoPOverEMax();
31ae6d59 1163 Float_t pOverEmin = GetHistogramRanges()->GetHistoPOverEMin();
09273901 1164
d2655d46 1165 Int_t bin[] = {0,2,4,6,10,15,20,100}; // energy bins for SS studies
1166
9e51e29a 1167 TString cut[] = {"Open","Reader","E","Time","NCells","NLM","Fidutial","Matching","Bad","PID"};
85c4406e 1168 for (Int_t i = 0; i < 10 ; i++)
fc195fd0 1169 {
58ea8ce5 1170 fhClusterCutsE[i] = new TH1F(Form("hE_Cut_%d_%s", i, cut[i].Data()),
fc195fd0 1171 Form("Number of clusters that pass cuts <= %d, %s", i, cut[i].Data()),
85c4406e 1172 nptbins,ptmin,ptmax);
58ea8ce5 1173 fhClusterCutsE[i]->SetYTitle("d#it{N}/d#it{E} ");
1174 fhClusterCutsE[i]->SetXTitle("#it{E} (GeV)");
1175 outputContainer->Add(fhClusterCutsE[i]) ;
1176
1177 fhClusterCutsPt[i] = new TH1F(Form("hPt_Cut_%d_%s", i, cut[i].Data()),
1178 Form("Number of clusters that pass cuts <= %d, %s", i, cut[i].Data()),
1179 nptbins,ptmin,ptmax);
1180 fhClusterCutsPt[i]->SetYTitle("d#it{N}/d#it{E} ");
1181 fhClusterCutsPt[i]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
1182 outputContainer->Add(fhClusterCutsPt[i]) ;
fc195fd0 1183 }
1184
6df33fcb 1185 fhEClusterSM = new TH2F("hEClusterSM","Raw clusters E and super-module number",
1186 nptbins,ptmin,ptmax,
1187 GetCaloUtils()->GetNumberOfSuperModulesUsed(),0,GetCaloUtils()->GetNumberOfSuperModulesUsed());
1188 fhEClusterSM->SetYTitle("SuperModule ");
184ca640 1189 fhEClusterSM->SetXTitle("#it{E} (GeV)");
6df33fcb 1190 outputContainer->Add(fhEClusterSM) ;
1191
11baad66 1192 fhPtClusterSM = new TH2F("hPtClusterSM","Raw clusters #it{p}_{T} and super-module number",
6df33fcb 1193 nptbins,ptmin,ptmax,
1194 GetCaloUtils()->GetNumberOfSuperModulesUsed(),0,GetCaloUtils()->GetNumberOfSuperModulesUsed());
1195 fhPtClusterSM->SetYTitle("SuperModule ");
184ca640 1196 fhPtClusterSM->SetXTitle("#it{E} (GeV)");
6df33fcb 1197 outputContainer->Add(fhPtClusterSM) ;
1198
1199 fhEPhotonSM = new TH2F("hEPhotonSM","Selected clusters E and super-module number",
1200 nptbins,ptmin,ptmax,
1201 GetCaloUtils()->GetNumberOfSuperModulesUsed(),0,GetCaloUtils()->GetNumberOfSuperModulesUsed());
1202 fhEPhotonSM->SetYTitle("SuperModule ");
184ca640 1203 fhEPhotonSM->SetXTitle("#it{E} (GeV)");
6df33fcb 1204 outputContainer->Add(fhEPhotonSM) ;
1205
11baad66 1206 fhPtPhotonSM = new TH2F("hPtPhotonSM","Selected clusters #it{p}_{T} and super-module number",
6df33fcb 1207 nptbins,ptmin,ptmax,
1208 GetCaloUtils()->GetNumberOfSuperModulesUsed(),0,GetCaloUtils()->GetNumberOfSuperModulesUsed());
1209 fhPtPhotonSM->SetYTitle("SuperModule ");
184ca640 1210 fhPtPhotonSM->SetXTitle("#it{E} (GeV)");
6df33fcb 1211 outputContainer->Add(fhPtPhotonSM) ;
1212
85c4406e 1213 fhNCellsE = new TH2F ("hNCellsE","# of cells in cluster vs E of clusters", nptbins,ptmin,ptmax, nbins,nmin,nmax);
184ca640 1214 fhNCellsE->SetXTitle("#it{E} (GeV)");
c4a7d28a 1215 fhNCellsE->SetYTitle("# of cells in cluster");
85c4406e 1216 outputContainer->Add(fhNCellsE);
f15c25da 1217
85c4406e 1218 fhCellsE = new TH2F ("hCellsE","energy of cells in cluster vs E of clusters", nptbins,ptmin,ptmax, nptbins*2,ptmin,ptmax);
184ca640 1219 fhCellsE->SetXTitle("#it{E}_{cluster} (GeV)");
1220 fhCellsE->SetYTitle("#it{E}_{cell} (GeV)");
85c4406e 1221 outputContainer->Add(fhCellsE);
5c46c992 1222
b2e375c7 1223 fhTimePt = new TH2F ("hTimePt","time of cluster vs pT of clusters", nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
11baad66 1224 fhTimePt->SetXTitle("#it{p}_{T} (GeV/#it{c})");
184ca640 1225 fhTimePt->SetYTitle("#it{time} (ns)");
b2e375c7 1226 outputContainer->Add(fhTimePt);
6175da48 1227
f66d95af 1228 fhMaxCellDiffClusterE = new TH2F ("hMaxCellDiffClusterE","energy vs difference of cluster energy - max cell energy / cluster energy, good clusters",
85c4406e 1229 nptbins,ptmin,ptmax, 500,0,1.);
184ca640 1230 fhMaxCellDiffClusterE->SetXTitle("#it{E}_{cluster} (GeV) ");
1231 fhMaxCellDiffClusterE->SetYTitle("(#it{E}_{cluster} - #it{E}_{cell max})/ #it{E}_{cluster}");
85c4406e 1232 outputContainer->Add(fhMaxCellDiffClusterE);
f66d95af 1233
85c4406e 1234 fhEPhoton = new TH1F("hEPhoton","Number of #gamma over calorimeter vs energy",nptbins,ptmin,ptmax);
184ca640 1235 fhEPhoton->SetYTitle("#it{counts}");
1236 fhEPhoton->SetXTitle("#it{E}_{#gamma}(GeV)");
85c4406e 1237 outputContainer->Add(fhEPhoton) ;
20218aea 1238
11baad66 1239 fhPtPhoton = new TH1F("hPtPhoton","Number of #gamma over calorimeter vs #it{p}_{T}",nptbins,ptmin,ptmax);
184ca640 1240 fhPtPhoton->SetYTitle("#it{counts}");
1241 fhPtPhoton->SetXTitle("p_{T #gamma}(GeV/#it{c})");
85c4406e 1242 outputContainer->Add(fhPtPhoton) ;
1243
11baad66 1244 fhPtCentralityPhoton = new TH2F("hPtCentralityPhoton","centrality vs #it{p}_{T}",nptbins,ptmin,ptmax, 100,0,100);
c8710850 1245 fhPtCentralityPhoton->SetYTitle("Centrality");
11baad66 1246 fhPtCentralityPhoton->SetXTitle("#it{p}_{T}(GeV/#it{c})");
c8710850 1247 outputContainer->Add(fhPtCentralityPhoton) ;
1248
11baad66 1249 fhPtEventPlanePhoton = new TH2F("hPtEventPlanePhoton","centrality vs #it{p}_{T}",nptbins,ptmin,ptmax, 100,0,TMath::Pi());
c8710850 1250 fhPtEventPlanePhoton->SetYTitle("Event plane angle (rad)");
11baad66 1251 fhPtEventPlanePhoton->SetXTitle("#it{p}_{T} (GeV/#it{c})");
c8710850 1252 outputContainer->Add(fhPtEventPlanePhoton) ;
85c4406e 1253
c2a62a94 1254 fhEtaPhi = new TH2F
184ca640 1255 ("hEtaPhi","cluster,#it{E} > 0.5 GeV, #eta vs #phi",netabins,etamin,etamax,nphibins,phimin,phimax);
c2a62a94 1256 fhEtaPhi->SetYTitle("#phi (rad)");
1257 fhEtaPhi->SetXTitle("#eta");
1258 outputContainer->Add(fhEtaPhi) ;
85c4406e 1259
477d6cee 1260 fhPhiPhoton = new TH2F
11baad66 1261 ("hPhiPhoton","#phi_{#gamma} vs #it{p}_{T}",nptbins,ptmin,ptmax,nphibins,phimin,phimax);
6175da48 1262 fhPhiPhoton->SetYTitle("#phi (rad)");
184ca640 1263 fhPhiPhoton->SetXTitle("p_{T #gamma} (GeV/#it{c})");
85c4406e 1264 outputContainer->Add(fhPhiPhoton) ;
477d6cee 1265
1266 fhEtaPhoton = new TH2F
11baad66 1267 ("hEtaPhoton","#eta_{#gamma} vs #it{p}_{T}",nptbins,ptmin,ptmax,netabins,etamin,etamax);
477d6cee 1268 fhEtaPhoton->SetYTitle("#eta");
184ca640 1269 fhEtaPhoton->SetXTitle("p_{T #gamma} (GeV/#it{c})");
477d6cee 1270 outputContainer->Add(fhEtaPhoton) ;
1271
6175da48 1272 fhEtaPhiPhoton = new TH2F
85c4406e 1273 ("hEtaPhiPhoton","#eta vs #phi",netabins,etamin,etamax,nphibins,phimin,phimax);
6175da48 1274 fhEtaPhiPhoton->SetYTitle("#phi (rad)");
1275 fhEtaPhiPhoton->SetXTitle("#eta");
1276 outputContainer->Add(fhEtaPhiPhoton) ;
34c16486 1277 if(GetMinPt() < 0.5)
1278 {
20218aea 1279 fhEtaPhi05Photon = new TH2F
74e3eb22 1280 ("hEtaPhi05Photon","#eta vs #phi, E < 0.5",netabins,etamin,etamax,nphibins,phimin,phimax);
20218aea 1281 fhEtaPhi05Photon->SetYTitle("#phi (rad)");
1282 fhEtaPhi05Photon->SetXTitle("#eta");
1283 outputContainer->Add(fhEtaPhi05Photon) ;
1284 }
85c4406e 1285
9e51e29a 1286 fhNLocMax = new TH2F("hNLocMax","Number of local maxima in cluster",
85c4406e 1287 nptbins,ptmin,ptmax,10,0,10);
9e51e29a 1288 fhNLocMax ->SetYTitle("N maxima");
184ca640 1289 fhNLocMax ->SetXTitle("#it{E} (GeV)");
85c4406e 1290 outputContainer->Add(fhNLocMax) ;
9e51e29a 1291
521636d2 1292 //Shower shape
34c16486 1293 if(fFillSSHistograms)
1294 {
85c4406e 1295 fhLam0E = new TH2F ("hLam0E","#lambda_{0}^{2} vs E", nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
521636d2 1296 fhLam0E->SetYTitle("#lambda_{0}^{2}");
184ca640 1297 fhLam0E->SetXTitle("#it{E} (GeV)");
85c4406e 1298 outputContainer->Add(fhLam0E);
521636d2 1299
85c4406e 1300 fhLam1E = new TH2F ("hLam1E","#lambda_{1}^{2} vs E", nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
521636d2 1301 fhLam1E->SetYTitle("#lambda_{1}^{2}");
184ca640 1302 fhLam1E->SetXTitle("#it{E} (GeV)");
85c4406e 1303 outputContainer->Add(fhLam1E);
521636d2 1304
85c4406e 1305 fhDispE = new TH2F ("hDispE"," dispersion^{2} vs E", nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
521636d2 1306 fhDispE->SetYTitle("D^{2}");
184ca640 1307 fhDispE->SetXTitle("#it{E} (GeV) ");
521636d2 1308 outputContainer->Add(fhDispE);
85c4406e 1309
b5dbb99b 1310 if(!fRejectTrackMatch)
1311 {
85c4406e 1312 fhLam0ETM = new TH2F ("hLam0ETM","#lambda_{0}^{2} vs E, cut on track-matching residual |#Delta #eta| < 0.05, |#Delta #phi| < 0.05", nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
b5dbb99b 1313 fhLam0ETM->SetYTitle("#lambda_{0}^{2}");
184ca640 1314 fhLam0ETM->SetXTitle("#it{E} (GeV)");
85c4406e 1315 outputContainer->Add(fhLam0ETM);
b5dbb99b 1316
85c4406e 1317 fhLam1ETM = new TH2F ("hLam1ETM","#lambda_{1}^{2} vs E, cut on track-matching residual |#Delta #eta| < 0.05, |#Delta #phi| < 0.05", nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
b5dbb99b 1318 fhLam1ETM->SetYTitle("#lambda_{1}^{2}");
184ca640 1319 fhLam1ETM->SetXTitle("#it{E} (GeV)");
85c4406e 1320 outputContainer->Add(fhLam1ETM);
b5dbb99b 1321
85c4406e 1322 fhDispETM = new TH2F ("hDispETM"," dispersion^{2} vs E, cut on track-matching residual |#Delta #eta| < 0.05, |#Delta #phi| < 0.05", nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
b5dbb99b 1323 fhDispETM->SetYTitle("D^{2}");
184ca640 1324 fhDispETM->SetXTitle("#it{E} (GeV) ");
b5dbb99b 1325 outputContainer->Add(fhDispETM);
1326 }
521636d2 1327
4d1d8f00 1328 if(fCalorimeter == "EMCAL" && GetFirstSMCoveredByTRD() >= 0)
b5dbb99b 1329 {
85c4406e 1330 fhLam0ETRD = new TH2F ("hLam0ETRD","#lambda_{0}^{2} vs E, EMCAL SM covered by TRD", nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
521636d2 1331 fhLam0ETRD->SetYTitle("#lambda_{0}^{2}");
184ca640 1332 fhLam0ETRD->SetXTitle("#it{E} (GeV)");
85c4406e 1333 outputContainer->Add(fhLam0ETRD);
521636d2 1334
85c4406e 1335 fhLam1ETRD = new TH2F ("hLam1ETRD","#lambda_{1}^{2} vs E, EMCAL SM covered by TRD", nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
521636d2 1336 fhLam1ETRD->SetYTitle("#lambda_{1}^{2}");
184ca640 1337 fhLam1ETRD->SetXTitle("#it{E} (GeV)");
85c4406e 1338 outputContainer->Add(fhLam1ETRD);
521636d2 1339
85c4406e 1340 fhDispETRD = new TH2F ("hDispETRD"," dispersion^{2} vs E, EMCAL SM covered by TRD", nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
521636d2 1341 fhDispETRD->SetYTitle("Dispersion^{2}");
184ca640 1342 fhDispETRD->SetXTitle("#it{E} (GeV) ");
b5dbb99b 1343 outputContainer->Add(fhDispETRD);
1344
4d1d8f00 1345 if(!fRejectTrackMatch && GetFirstSMCoveredByTRD() >=0 )
b5dbb99b 1346 {
85c4406e 1347 fhLam0ETMTRD = new TH2F ("hLam0ETMTRD","#lambda_{0}^{2} vs E, EMCAL SM covered by TRD, cut on track-matching residual |#Delta #eta| < 0.05, |#Delta #phi| < 0.05", nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
b5dbb99b 1348 fhLam0ETMTRD->SetYTitle("#lambda_{0}^{2}");
184ca640 1349 fhLam0ETMTRD->SetXTitle("#it{E} (GeV)");
85c4406e 1350 outputContainer->Add(fhLam0ETMTRD);
b5dbb99b 1351
85c4406e 1352 fhLam1ETMTRD = new TH2F ("hLam1ETMTRD","#lambda_{1}^{2} vs E, EMCAL SM covered by TRD, cut on track-matching residual |#Delta #eta| < 0.05, |#Delta #phi| < 0.05", nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
b5dbb99b 1353 fhLam1ETMTRD->SetYTitle("#lambda_{1}^{2}");
184ca640 1354 fhLam1ETMTRD->SetXTitle("#it{E} (GeV)");
85c4406e 1355 outputContainer->Add(fhLam1ETMTRD);
b5dbb99b 1356
85c4406e 1357 fhDispETMTRD = new TH2F ("hDispETMTRD"," dispersion^{2} vs E, EMCAL SM covered by TRD, cut on track-matching residual |#Delta #eta| < 0.05, |#Delta #phi| < 0.05", nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
b5dbb99b 1358 fhDispETMTRD->SetYTitle("Dispersion^{2}");
184ca640 1359 fhDispETMTRD->SetXTitle("#it{E} (GeV) ");
85c4406e 1360 outputContainer->Add(fhDispETMTRD);
1361 }
1362 }
521636d2 1363
764ab1f4 1364 if(!fFillOnlySimpleSSHisto)
34c16486 1365 {
85c4406e 1366 fhNCellsLam0LowE = new TH2F ("hNCellsLam0LowE","N_{cells} in cluster vs #lambda_{0}^{2}, E < 2 GeV", nbins,nmin, nmax, ssbins,ssmin,ssmax);
764ab1f4 1367 fhNCellsLam0LowE->SetXTitle("N_{cells}");
1368 fhNCellsLam0LowE->SetYTitle("#lambda_{0}^{2}");
85c4406e 1369 outputContainer->Add(fhNCellsLam0LowE);
764ab1f4 1370
184ca640 1371 fhNCellsLam0HighE = new TH2F ("hNCellsLam0HighE","N_{cells} in cluster vs #lambda_{0}^{2}, #it{E} > 2 GeV", nbins,nmin, nmax, ssbins,ssmin,ssmax);
764ab1f4 1372 fhNCellsLam0HighE->SetXTitle("N_{cells}");
1373 fhNCellsLam0HighE->SetYTitle("#lambda_{0}^{2}");
85c4406e 1374 outputContainer->Add(fhNCellsLam0HighE);
764ab1f4 1375
85c4406e 1376 fhNCellsLam1LowE = new TH2F ("hNCellsLam1LowE","N_{cells} in cluster vs #lambda_{1}^{2}, E < 2 GeV", nbins,nmin, nmax, ssbins,ssmin,ssmax);
764ab1f4 1377 fhNCellsLam1LowE->SetXTitle("N_{cells}");
1378 fhNCellsLam1LowE->SetYTitle("#lambda_{0}^{2}");
85c4406e 1379 outputContainer->Add(fhNCellsLam1LowE);
764ab1f4 1380
184ca640 1381 fhNCellsLam1HighE = new TH2F ("hNCellsLam1HighE","N_{cells} in cluster vs #lambda_{1}^{2}, #it{E} > 2 GeV", nbins,nmin, nmax, ssbins,ssmin,ssmax);
764ab1f4 1382 fhNCellsLam1HighE->SetXTitle("N_{cells}");
1383 fhNCellsLam1HighE->SetYTitle("#lambda_{0}^{2}");
85c4406e 1384 outputContainer->Add(fhNCellsLam1HighE);
764ab1f4 1385
85c4406e 1386 fhNCellsDispLowE = new TH2F ("hNCellsDispLowE","N_{cells} in cluster vs dispersion^{2}, E < 2 GeV", nbins,nmin, nmax, ssbins,ssmin,ssmax);
764ab1f4 1387 fhNCellsDispLowE->SetXTitle("N_{cells}");
1388 fhNCellsDispLowE->SetYTitle("D^{2}");
85c4406e 1389 outputContainer->Add(fhNCellsDispLowE);
764ab1f4 1390
85c4406e 1391 fhNCellsDispHighE = new TH2F ("hNCellsDispHighE","N_{cells} in cluster vs dispersion^{2}, E < 2 GeV", nbins,nmin, nmax, ssbins,ssmin,ssmax);
764ab1f4 1392 fhNCellsDispHighE->SetXTitle("N_{cells}");
1393 fhNCellsDispHighE->SetYTitle("D^{2}");
85c4406e 1394 outputContainer->Add(fhNCellsDispHighE);
764ab1f4 1395
85c4406e 1396 fhEtaLam0LowE = new TH2F ("hEtaLam0LowE","#eta vs #lambda_{0}^{2}, E < 2 GeV", netabins,etamin,etamax, ssbins,ssmin,ssmax);
764ab1f4 1397 fhEtaLam0LowE->SetYTitle("#lambda_{0}^{2}");
1398 fhEtaLam0LowE->SetXTitle("#eta");
85c4406e 1399 outputContainer->Add(fhEtaLam0LowE);
764ab1f4 1400
85c4406e 1401 fhPhiLam0LowE = new TH2F ("hPhiLam0LowE","#phi vs #lambda_{0}^{2}, E < 2 GeV", nphibins,phimin,phimax, ssbins,ssmin,ssmax);
764ab1f4 1402 fhPhiLam0LowE->SetYTitle("#lambda_{0}^{2}");
1403 fhPhiLam0LowE->SetXTitle("#phi");
85c4406e 1404 outputContainer->Add(fhPhiLam0LowE);
764ab1f4 1405
184ca640 1406 fhEtaLam0HighE = new TH2F ("hEtaLam0HighE","#eta vs #lambda_{0}^{2}, #it{E} > 2 GeV", netabins,etamin,etamax, ssbins,ssmin,ssmax);
764ab1f4 1407 fhEtaLam0HighE->SetYTitle("#lambda_{0}^{2}");
1408 fhEtaLam0HighE->SetXTitle("#eta");
85c4406e 1409 outputContainer->Add(fhEtaLam0HighE);
764ab1f4 1410
184ca640 1411 fhPhiLam0HighE = new TH2F ("hPhiLam0HighE","#phi vs #lambda_{0}^{2}, #it{E} > 2 GeV", nphibins,phimin,phimax, ssbins,ssmin,ssmax);
764ab1f4 1412 fhPhiLam0HighE->SetYTitle("#lambda_{0}^{2}");
1413 fhPhiLam0HighE->SetXTitle("#phi");
85c4406e 1414 outputContainer->Add(fhPhiLam0HighE);
764ab1f4 1415
85c4406e 1416 fhLam1Lam0LowE = new TH2F ("hLam1Lam0LowE","#lambda_{0}^{2} vs #lambda_{1}^{2} in cluster of E < 2 GeV", ssbins,ssmin,ssmax, ssbins,ssmin,ssmax);
764ab1f4 1417 fhLam1Lam0LowE->SetYTitle("#lambda_{0}^{2}");
1418 fhLam1Lam0LowE->SetXTitle("#lambda_{1}^{2}");
85c4406e 1419 outputContainer->Add(fhLam1Lam0LowE);
764ab1f4 1420
184ca640 1421 fhLam1Lam0HighE = new TH2F ("hLam1Lam0HighE","#lambda_{0}^{2} vs #lambda_{1}^{2} in cluster of #it{E} > 2 GeV", ssbins,ssmin,ssmax, ssbins,ssmin,ssmax);
764ab1f4 1422 fhLam1Lam0HighE->SetYTitle("#lambda_{0}^{2}");
1423 fhLam1Lam0HighE->SetXTitle("#lambda_{1}^{2}");
85c4406e 1424 outputContainer->Add(fhLam1Lam0HighE);
764ab1f4 1425
85c4406e 1426 fhLam0DispLowE = new TH2F ("hLam0DispLowE","#lambda_{0}^{2} vs dispersion^{2} in cluster of E < 2 GeV", ssbins,ssmin,ssmax, ssbins,ssmin,ssmax);
764ab1f4 1427 fhLam0DispLowE->SetXTitle("#lambda_{0}^{2}");
1428 fhLam0DispLowE->SetYTitle("D^{2}");
85c4406e 1429 outputContainer->Add(fhLam0DispLowE);
764ab1f4 1430
184ca640 1431 fhLam0DispHighE = new TH2F ("hLam0DispHighE","#lambda_{0}^{2} vs dispersion^{2} in cluster of #it{E} > 2 GeV", ssbins,ssmin,ssmax, ssbins,ssmin,ssmax);
764ab1f4 1432 fhLam0DispHighE->SetXTitle("#lambda_{0}^{2}");
1433 fhLam0DispHighE->SetYTitle("D^{2}");
85c4406e 1434 outputContainer->Add(fhLam0DispHighE);
764ab1f4 1435
85c4406e 1436 fhDispLam1LowE = new TH2F ("hDispLam1LowE","Dispersion^{2} vs #lambda_{1}^{2} in cluster of E < 2 GeV", ssbins,ssmin,ssmax, ssbins,ssmin,ssmax);
764ab1f4 1437 fhDispLam1LowE->SetXTitle("D^{2}");
1438 fhDispLam1LowE->SetYTitle("#lambda_{1}^{2}");
85c4406e 1439 outputContainer->Add(fhDispLam1LowE);
764ab1f4 1440
184ca640 1441 fhDispLam1HighE = new TH2F ("hDispLam1HighE","Dispersion^{2} vs #lambda_{1^{2}} in cluster of #it{E} > 2 GeV", ssbins,ssmin,ssmax, ssbins,ssmin,ssmax);
764ab1f4 1442 fhDispLam1HighE->SetXTitle("D^{2}");
1443 fhDispLam1HighE->SetYTitle("#lambda_{1}^{2}");
85c4406e 1444 outputContainer->Add(fhDispLam1HighE);
764ab1f4 1445
1446 if(fCalorimeter == "EMCAL")
34c16486 1447 {
85c4406e 1448 fhDispEtaE = new TH2F ("hDispEtaE","#sigma^{2}_{#eta #eta} = #Sigma w_{i}(#eta_{i} - <#eta>)^{2}/ #Sigma w_{i} vs E", nptbins,ptmin,ptmax, ssbins,ssmin,ssmax);
184ca640 1449 fhDispEtaE->SetXTitle("#it{E} (GeV)");
764ab1f4 1450 fhDispEtaE->SetYTitle("#sigma^{2}_{#eta #eta}");
85c4406e 1451 outputContainer->Add(fhDispEtaE);
764ab1f4 1452
85c4406e 1453 fhDispPhiE = new TH2F ("hDispPhiE","#sigma^{2}_{#phi #phi} = #Sigma w_{i}(#phi_{i} - <#phi>)^{2} / #Sigma w_{i} vs E", nptbins,ptmin,ptmax, ssbins,ssmin,ssmax);
184ca640 1454 fhDispPhiE->SetXTitle("#it{E} (GeV)");
764ab1f4 1455 fhDispPhiE->SetYTitle("#sigma^{2}_{#phi #phi}");
85c4406e 1456 outputContainer->Add(fhDispPhiE);
764ab1f4 1457
85c4406e 1458 fhSumEtaE = new TH2F ("hSumEtaE","#delta^{2}_{#eta #eta} = #Sigma w_{i}(#eta_{i})^{2} / #Sigma w_{i} - <#eta>^{2} vs E", nptbins,ptmin,ptmax, ssbins,ssmin,ssmax);
184ca640 1459 fhSumEtaE->SetXTitle("#it{E} (GeV)");
764ab1f4 1460 fhSumEtaE->SetYTitle("#delta^{2}_{#eta #eta}");
85c4406e 1461 outputContainer->Add(fhSumEtaE);
764ab1f4 1462
85c4406e 1463 fhSumPhiE = new TH2F ("hSumPhiE","#delta^{2}_{#phi #phi} = #Sigma w_{i}(#phi_{i})^{2}/ #Sigma w_{i} - <#phi>^{2} vs E",
1464 nptbins,ptmin,ptmax, ssbins,ssmin,ssmax);
184ca640 1465 fhSumPhiE->SetXTitle("#it{E} (GeV)");
764ab1f4 1466 fhSumPhiE->SetYTitle("#delta^{2}_{#phi #phi}");
85c4406e 1467 outputContainer->Add(fhSumPhiE);
764ab1f4 1468
85c4406e 1469 fhSumEtaPhiE = new TH2F ("hSumEtaPhiE","#delta^{2}_{#eta #phi} = #Sigma w_{i}(#phi_{i} #eta_{i} ) / #Sigma w_{i} - <#phi><#eta> vs E",
1470 nptbins,ptmin,ptmax, 2*ssbins,-ssmax,ssmax);
184ca640 1471 fhSumEtaPhiE->SetXTitle("#it{E} (GeV)");
764ab1f4 1472 fhSumEtaPhiE->SetYTitle("#delta^{2}_{#eta #phi}");
1473 outputContainer->Add(fhSumEtaPhiE);
1474
85c4406e 1475 fhDispEtaPhiDiffE = new TH2F ("hDispEtaPhiDiffE","#sigma^{2}_{#phi #phi} - #sigma^{2}_{#eta #eta} vs E",
1476 nptbins,ptmin,ptmax,200, -10,10);
184ca640 1477 fhDispEtaPhiDiffE->SetXTitle("#it{E} (GeV)");
764ab1f4 1478 fhDispEtaPhiDiffE->SetYTitle("#sigma^{2}_{#phi #phi}-#sigma^{2}_{#eta #eta}");
85c4406e 1479 outputContainer->Add(fhDispEtaPhiDiffE);
bfdcf7fb 1480
85c4406e 1481 fhSphericityE = new TH2F ("hSphericityE","(#sigma^{2}_{#phi #phi} - #sigma^{2}_{#eta #eta}) / (#sigma^{2}_{#eta #eta} + #sigma^{2}_{#phi #phi}) vs E",
1482 nptbins,ptmin,ptmax, 200, -1,1);
184ca640 1483 fhSphericityE->SetXTitle("#it{E} (GeV)");
764ab1f4 1484 fhSphericityE->SetYTitle("s = (#sigma^{2}_{#phi #phi} - #sigma^{2}_{#eta #eta}) / (#sigma^{2}_{#eta #eta} + #sigma^{2}_{#phi #phi})");
1485 outputContainer->Add(fhSphericityE);
bfdcf7fb 1486
85c4406e 1487 fhDispSumEtaDiffE = new TH2F ("hDispSumEtaDiffE","#sigma^{2}_{#eta #eta} - #delta^{2}_{#eta #eta} / average vs E", nptbins,ptmin,ptmax, 200,-0.01,0.01);
184ca640 1488 fhDispSumEtaDiffE->SetXTitle("#it{E} (GeV)");
764ab1f4 1489 fhDispSumEtaDiffE->SetYTitle("#sigma^{2}_{#eta #eta} - #delta^{2}_{#eta #eta} / average");
85c4406e 1490 outputContainer->Add(fhDispSumEtaDiffE);
764ab1f4 1491
85c4406e 1492 fhDispSumPhiDiffE = new TH2F ("hDispSumPhiDiffE","#sigma^{2}_{#phi #phi} - #delta^{2}_{#phi #phi} / average vs E", nptbins,ptmin,ptmax, 200,-0.01,0.01);
184ca640 1493 fhDispSumPhiDiffE->SetXTitle("#it{E} (GeV)");
764ab1f4 1494 fhDispSumPhiDiffE->SetYTitle("#sigma^{2}_{#phi #phi} - #delta^{2}_{#phi #phi} / average");
85c4406e 1495 outputContainer->Add(fhDispSumPhiDiffE);
764ab1f4 1496
1497 for(Int_t i = 0; i < 7; i++)
1498 {
85c4406e 1499 fhDispEtaDispPhi[i] = new TH2F (Form("hDispEtaDispPhi_EBin%d",i),Form("#sigma^{2}_{#phi #phi} vs #sigma^{2}_{#eta #eta} for %d < E < %d GeV",bin[i],bin[i+1]),
1500 ssbins,ssmin,ssmax , ssbins,ssmin,ssmax);
764ab1f4 1501 fhDispEtaDispPhi[i]->SetXTitle("#sigma^{2}_{#eta #eta}");
1502 fhDispEtaDispPhi[i]->SetYTitle("#sigma^{2}_{#phi #phi}");
85c4406e 1503 outputContainer->Add(fhDispEtaDispPhi[i]);
764ab1f4 1504
85c4406e 1505 fhLambda0DispEta[i] = new TH2F (Form("hLambda0DispEta_EBin%d",i),Form("#lambda^{2}_{0} vs #sigma^{2}_{#eta #eta} for %d < E < %d GeV",bin[i],bin[i+1]),
1506 ssbins,ssmin,ssmax , ssbins,ssmin,ssmax);
764ab1f4 1507 fhLambda0DispEta[i]->SetXTitle("#lambda^{2}_{0}");
1508 fhLambda0DispEta[i]->SetYTitle("#sigma^{2}_{#eta #eta}");
85c4406e 1509 outputContainer->Add(fhLambda0DispEta[i]);
764ab1f4 1510
85c4406e 1511 fhLambda0DispPhi[i] = new TH2F (Form("hLambda0DispPhi_EBin%d",i),Form("#lambda^{2}_{0}} vs #sigma^{2}_{#phi #phi} for %d < E < %d GeV",bin[i],bin[i+1]),
1512 ssbins,ssmin,ssmax , ssbins,ssmin,ssmax);
764ab1f4 1513 fhLambda0DispPhi[i]->SetXTitle("#lambda^{2}_{0}");
1514 fhLambda0DispPhi[i]->SetYTitle("#sigma^{2}_{#phi #phi}");
85c4406e 1515 outputContainer->Add(fhLambda0DispPhi[i]);
764ab1f4 1516 }
34c16486 1517 }
1518 }
521636d2 1519 } // Shower shape
1520
09273901 1521 // Track Matching
1522
b5dbb99b 1523 if(fFillTMHisto)
1524 {
b2e375c7 1525 TString cutTM [] = {"NoCut",""};
b5dbb99b 1526
b2e375c7 1527 for(Int_t i = 0; i < 2; i++)
31ae6d59 1528 {
b2e375c7 1529 fhTrackMatchedDEta[i] = new TH2F
1530 (Form("hTrackMatchedDEta%s",cutTM[i].Data()),
1531 Form("d#eta of cluster-track vs cluster energy, %s",cutTM[i].Data()),
85c4406e 1532 nptbins,ptmin,ptmax,nresetabins,resetamin,resetamax);
b2e375c7 1533 fhTrackMatchedDEta[i]->SetYTitle("d#eta");
184ca640 1534 fhTrackMatchedDEta[i]->SetXTitle("#it{E}_{cluster} (GeV)");
4bfeae64 1535
b2e375c7 1536 fhTrackMatchedDPhi[i] = new TH2F
1537 (Form("hTrackMatchedDPhi%s",cutTM[i].Data()),
1538 Form("d#phi of cluster-track vs cluster energy, %s",cutTM[i].Data()),
85c4406e 1539 nptbins,ptmin,ptmax,nresphibins,resphimin,resphimax);
b2e375c7 1540 fhTrackMatchedDPhi[i]->SetYTitle("d#phi (rad)");
184ca640 1541 fhTrackMatchedDPhi[i]->SetXTitle("#it{E}_{cluster} (GeV)");
b2e375c7 1542
1543 fhTrackMatchedDEtaDPhi[i] = new TH2F
1544 (Form("hTrackMatchedDEtaDPhi%s",cutTM[i].Data()),
1545 Form("d#eta vs d#phi of cluster-track vs cluster energy, %s",cutTM[i].Data()),
1546 nresetabins,resetamin,resetamax,nresphibins,resphimin,resphimax);
1547 fhTrackMatchedDEtaDPhi[i]->SetYTitle("d#phi (rad)");
1548 fhTrackMatchedDEtaDPhi[i]->SetXTitle("d#eta");
1549
1550 fhTrackMatchedDEtaPos[i] = new TH2F
1551 (Form("hTrackMatchedDEtaPos%s",cutTM[i].Data()),
1552 Form("d#eta of cluster-track vs cluster energy, %s",cutTM[i].Data()),
85c4406e 1553 nptbins,ptmin,ptmax,nresetabins,resetamin,resetamax);
b2e375c7 1554 fhTrackMatchedDEtaPos[i]->SetYTitle("d#eta");
184ca640 1555 fhTrackMatchedDEtaPos[i]->SetXTitle("#it{E}_{cluster} (GeV)");
8d6b7f60 1556
b2e375c7 1557 fhTrackMatchedDPhiPos[i] = new TH2F
1558 (Form("hTrackMatchedDPhiPos%s",cutTM[i].Data()),
1559 Form("d#phi of cluster-track vs cluster energy, %s",cutTM[i].Data()),
85c4406e 1560 nptbins,ptmin,ptmax,nresphibins,resphimin,resphimax);
b2e375c7 1561 fhTrackMatchedDPhiPos[i]->SetYTitle("d#phi (rad)");
184ca640 1562 fhTrackMatchedDPhiPos[i]->SetXTitle("#it{E}_{cluster} (GeV)");
b2e375c7 1563
1564 fhTrackMatchedDEtaDPhiPos[i] = new TH2F
1565 (Form("hTrackMatchedDEtaDPhiPos%s",cutTM[i].Data()),
1566 Form("d#eta vs d#phi of cluster-track vs cluster energy, %s",cutTM[i].Data()),
1567 nresetabins,resetamin,resetamax,nresphibins,resphimin,resphimax);
1568 fhTrackMatchedDEtaDPhiPos[i]->SetYTitle("d#phi (rad)");
1569 fhTrackMatchedDEtaDPhiPos[i]->SetXTitle("d#eta");
1570
1571 fhTrackMatchedDEtaNeg[i] = new TH2F
1572 (Form("hTrackMatchedDEtaNeg%s",cutTM[i].Data()),
1573 Form("d#eta of cluster-track vs cluster energy, %s",cutTM[i].Data()),
85c4406e 1574 nptbins,ptmin,ptmax,nresetabins,resetamin,resetamax);
b2e375c7 1575 fhTrackMatchedDEtaNeg[i]->SetYTitle("d#eta");
184ca640 1576 fhTrackMatchedDEtaNeg[i]->SetXTitle("#it{E}_{cluster} (GeV)");
4bfeae64 1577
b2e375c7 1578 fhTrackMatchedDPhiNeg[i] = new TH2F
1579 (Form("hTrackMatchedDPhiNeg%s",cutTM[i].Data()),
1580 Form("d#phi of cluster-track vs cluster energy, %s",cutTM[i].Data()),
85c4406e 1581 nptbins,ptmin,ptmax,nresphibins,resphimin,resphimax);
b2e375c7 1582 fhTrackMatchedDPhiNeg[i]->SetYTitle("d#phi (rad)");
184ca640 1583 fhTrackMatchedDPhiNeg[i]->SetXTitle("#it{E}_{cluster} (GeV)");
b2e375c7 1584
1585 fhTrackMatchedDEtaDPhiNeg[i] = new TH2F
1586 (Form("hTrackMatchedDEtaDPhiNeg%s",cutTM[i].Data()),
1587 Form("d#eta vs d#phi of cluster-track vs cluster energy, %s",cutTM[i].Data()),
1588 nresetabins,resetamin,resetamax,nresphibins,resphimin,resphimax);
1589 fhTrackMatchedDEtaDPhiNeg[i]->SetYTitle("d#phi (rad)");
1590 fhTrackMatchedDEtaDPhiNeg[i]->SetXTitle("d#eta");
1591
1592 fhdEdx[i] = new TH2F (Form("hdEdx%s",cutTM[i].Data()),Form("matched track <dE/dx> vs cluster E, %s",cutTM[i].Data()),
1593 nptbins,ptmin,ptmax,ndedxbins, dedxmin, dedxmax);
184ca640 1594 fhdEdx[i]->SetXTitle("#it{E} (GeV)");
b2e375c7 1595 fhdEdx[i]->SetYTitle("<dE/dx>");
1596
1597 fhEOverP[i] = new TH2F (Form("hEOverP%s",cutTM[i].Data()),Form("matched track E/p vs cluster E, %s",cutTM[i].Data()),
1598 nptbins,ptmin,ptmax,nPoverEbins,pOverEmin,pOverEmax);
184ca640 1599 fhEOverP[i]->SetXTitle("#it{E} (GeV)");
b2e375c7 1600 fhEOverP[i]->SetYTitle("E/p");
1601
1602 outputContainer->Add(fhTrackMatchedDEta[i]) ;
1603 outputContainer->Add(fhTrackMatchedDPhi[i]) ;
1604 outputContainer->Add(fhTrackMatchedDEtaDPhi[i]) ;
1605 outputContainer->Add(fhTrackMatchedDEtaPos[i]) ;
1606 outputContainer->Add(fhTrackMatchedDPhiPos[i]) ;
1607 outputContainer->Add(fhTrackMatchedDEtaDPhiPos[i]) ;
1608 outputContainer->Add(fhTrackMatchedDEtaNeg[i]) ;
1609 outputContainer->Add(fhTrackMatchedDPhiNeg[i]) ;
1610 outputContainer->Add(fhTrackMatchedDEtaDPhiNeg[i]) ;
1611 outputContainer->Add(fhdEdx[i]);
1612 outputContainer->Add(fhEOverP[i]);
1613
4d1d8f00 1614 if(fCalorimeter=="EMCAL" && GetFirstSMCoveredByTRD() >=0 )
b2e375c7 1615 {
1616 fhTrackMatchedDEtaTRD[i] = new TH2F
1617 (Form("hTrackMatchedDEtaTRD%s",cutTM[i].Data()),
1618 Form("d#eta of cluster-track vs cluster energy, SM behind TRD, %s",cutTM[i].Data()),
1619 nptbins,ptmin,ptmax,nresetabins,resetamin,resetamax);
1620 fhTrackMatchedDEtaTRD[i]->SetYTitle("d#eta");
184ca640 1621 fhTrackMatchedDEtaTRD[i]->SetXTitle("#it{E}_{cluster} (GeV)");
b2e375c7 1622
1623 fhTrackMatchedDPhiTRD[i] = new TH2F
1624 (Form("hTrackMatchedDPhiTRD%s",cutTM[i].Data()),
1625 Form("d#phi of cluster-track vs cluster energy, SM behing TRD, %s",cutTM[i].Data()),
1626 nptbins,ptmin,ptmax,nresphibins,resphimin,resphimax);
1627 fhTrackMatchedDPhiTRD[i]->SetYTitle("d#phi (rad)");
184ca640 1628 fhTrackMatchedDPhiTRD[i]->SetXTitle("#it{E}_{cluster} (GeV)");
b2e375c7 1629
1630 fhEOverPTRD[i] = new TH2F
1631 (Form("hEOverPTRD%s",cutTM[i].Data()),
1632 Form("matched track E/p vs cluster E, behind TRD, %s",cutTM[i].Data()),
1633 nptbins,ptmin,ptmax,nPoverEbins,pOverEmin,pOverEmax);
184ca640 1634 fhEOverPTRD[i]->SetXTitle("#it{E} (GeV)");
b2e375c7 1635 fhEOverPTRD[i]->SetYTitle("E/p");
1636
1637 outputContainer->Add(fhTrackMatchedDEtaTRD[i]) ;
1638 outputContainer->Add(fhTrackMatchedDPhiTRD[i]) ;
1639 outputContainer->Add(fhEOverPTRD[i]);
1640 }
8d6b7f60 1641
b2e375c7 1642 if(IsDataMC())
1643 {
1644 fhTrackMatchedDEtaMCNoOverlap[i] = new TH2F
1645 (Form("hTrackMatchedDEtaMCNoOverlap%s",cutTM[i].Data()),
1646 Form("d#eta of cluster-track vs cluster energy, no other MC particles overlap %s",cutTM[i].Data()),
1647 nptbins,ptmin,ptmax,nresetabins,resetamin,resetamax);
1648 fhTrackMatchedDEtaMCNoOverlap[i]->SetYTitle("d#eta");
184ca640 1649 fhTrackMatchedDEtaMCNoOverlap[i]->SetXTitle("#it{E}_{cluster} (GeV)");
b2e375c7 1650
1651 fhTrackMatchedDPhiMCNoOverlap[i] = new TH2F
1652 (Form("hTrackMatchedDPhiMCNoOverlap%s",cutTM[i].Data()),
1653 Form("d#phi of cluster-track vs cluster energy, no other MC particles overlap %s",cutTM[i].Data()),
1654 nptbins,ptmin,ptmax,nresphibins,resphimin,resphimax);
1655 fhTrackMatchedDPhiMCNoOverlap[i]->SetYTitle("d#phi (rad)");
184ca640 1656 fhTrackMatchedDPhiMCNoOverlap[i]->SetXTitle("#it{E}_{cluster} (GeV)");
b2e375c7 1657
1658 outputContainer->Add(fhTrackMatchedDEtaMCNoOverlap[i]) ;
1659 outputContainer->Add(fhTrackMatchedDPhiMCNoOverlap[i]) ;
1660 fhTrackMatchedDEtaMCOverlap[i] = new TH2F
1661 (Form("hTrackMatchedDEtaMCOverlap%s",cutTM[i].Data()),
1662 Form("d#eta of cluster-track vs cluster energy, several MC particles overlap %s",cutTM[i].Data()),
1663 nptbins,ptmin,ptmax,nresetabins,resetamin,resetamax);
1664 fhTrackMatchedDEtaMCOverlap[i]->SetYTitle("d#eta");
184ca640 1665 fhTrackMatchedDEtaMCOverlap[i]->SetXTitle("#it{E}_{cluster} (GeV)");
b2e375c7 1666
1667 fhTrackMatchedDPhiMCOverlap[i] = new TH2F
1668 (Form("hTrackMatchedDPhiMCOverlap%s",cutTM[i].Data()),
1669 Form("d#phi of cluster-track vs cluster energy, several MC particles overlap %s",cutTM[i].Data()),
1670 nptbins,ptmin,ptmax,nresphibins,resphimin,resphimax);
1671 fhTrackMatchedDPhiMCOverlap[i]->SetYTitle("d#phi (rad)");
184ca640 1672 fhTrackMatchedDPhiMCOverlap[i]->SetXTitle("#it{E}_{cluster} (GeV)");
b2e375c7 1673
1674 outputContainer->Add(fhTrackMatchedDEtaMCOverlap[i]) ;
1675 outputContainer->Add(fhTrackMatchedDPhiMCOverlap[i]) ;
1676
1677 fhTrackMatchedDEtaMCConversion[i] = new TH2F
1678 (Form("hTrackMatchedDEtaMCConversion%s",cutTM[i].Data()),
1679 Form("d#eta of cluster-track vs cluster energy, no other MC particles overlap appart from conversions %s",cutTM[i].Data()),
1680 nptbins,ptmin,ptmax,nresetabins,resetamin,resetamax);
1681 fhTrackMatchedDEtaMCConversion[i]->SetYTitle("d#eta");
184ca640 1682 fhTrackMatchedDEtaMCConversion[i]->SetXTitle("#it{E}_{cluster} (GeV)");
b2e375c7 1683
1684 fhTrackMatchedDPhiMCConversion[i] = new TH2F
1685 (Form("hTrackMatchedDPhiMCConversion%s",cutTM[i].Data()),
1686 Form("d#phi of cluster-track vs cluster energy, no other MC particles overlap appart from conversions %s",cutTM[i].Data()),
1687 nptbins,ptmin,ptmax,nresphibins,resphimin,resphimax);
1688 fhTrackMatchedDPhiMCConversion[i]->SetYTitle("d#phi (rad)");
184ca640 1689 fhTrackMatchedDPhiMCConversion[i]->SetXTitle("#it{E}_{cluster} (GeV)");
b2e375c7 1690
1691 outputContainer->Add(fhTrackMatchedDEtaMCConversion[i]) ;
1692 outputContainer->Add(fhTrackMatchedDPhiMCConversion[i]) ;
1693
1694 fhTrackMatchedMCParticle[i] = new TH2F
1695 (Form("hTrackMatchedMCParticle%s",cutTM[i].Data()),
1696 Form("Origin of particle vs energy %s",cutTM[i].Data()),
1697 nptbins,ptmin,ptmax,8,0,8);
184ca640 1698 fhTrackMatchedMCParticle[i]->SetXTitle("#it{E} (GeV)");
b2e375c7 1699 //fhTrackMatchedMCParticle[i]->SetYTitle("Particle type");
1700
1701 fhTrackMatchedMCParticle[i]->GetYaxis()->SetBinLabel(1 ,"Photon");
1702 fhTrackMatchedMCParticle[i]->GetYaxis()->SetBinLabel(2 ,"Electron");
1703 fhTrackMatchedMCParticle[i]->GetYaxis()->SetBinLabel(3 ,"Meson Merged");
1704 fhTrackMatchedMCParticle[i]->GetYaxis()->SetBinLabel(4 ,"Rest");
1705 fhTrackMatchedMCParticle[i]->GetYaxis()->SetBinLabel(5 ,"Conv. Photon");
1706 fhTrackMatchedMCParticle[i]->GetYaxis()->SetBinLabel(6 ,"Conv. Electron");
1707 fhTrackMatchedMCParticle[i]->GetYaxis()->SetBinLabel(7 ,"Conv. Merged");
1708 fhTrackMatchedMCParticle[i]->GetYaxis()->SetBinLabel(8 ,"Conv. Rest");
1709
1710 outputContainer->Add(fhTrackMatchedMCParticle[i]);
1711 }
31ae6d59 1712 }
85c4406e 1713 }
09273901 1714
2ad19c3d 1715 if(fFillPileUpHistograms)
1716 {
5e5e056f 1717 TString pileUpName[] = {"SPD","EMCAL","SPDOrEMCAL","SPDAndEMCAL","SPDAndNotEMCAL","EMCALAndNotSPD","NotSPDAndNotEMCAL"} ;
1718
1719 for(Int_t i = 0 ; i < 7 ; i++)
1720 {
1721 fhPtPhotonPileUp[i] = new TH1F(Form("hPtPhotonPileUp%s",pileUpName[i].Data()),
11baad66 1722 Form("Selected photon #it{p}_{T} distribution, %s Pile-Up event",pileUpName[i].Data()), nptbins,ptmin,ptmax);
1723 fhPtPhotonPileUp[i]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
5e5e056f 1724 outputContainer->Add(fhPtPhotonPileUp[i]);
fad96885 1725
fad96885 1726 fhClusterTimeDiffPhotonPileUp[i] = new TH2F(Form("hClusterTimeDiffPhotonPileUp%s",pileUpName[i].Data()),
bc41680b 1727 Form("Photon cluster E vs #it{t}_{max}-#it{t}_{cell} in cluster, %s Pile-Up event",pileUpName[i].Data()),
1728 nptbins,ptmin,ptmax,400,-200,200);
184ca640 1729 fhClusterTimeDiffPhotonPileUp[i]->SetXTitle("#it{E} (GeV)");
bc41680b 1730 fhClusterTimeDiffPhotonPileUp[i]->SetYTitle("#it{t}_{max}-#it{t}_{cell} (ns)");
fad96885 1731 outputContainer->Add(fhClusterTimeDiffPhotonPileUp[i]);
5e5e056f 1732 }
1733
b2e375c7 1734 fhTimePtPhotonNoCut = new TH2F ("hTimePtPhoton_NoCut","time of photon cluster vs pT of clusters, no cut", nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
11baad66 1735 fhTimePtPhotonNoCut->SetXTitle("#it{p}_{T} (GeV/#it{c})");
184ca640 1736 fhTimePtPhotonNoCut->SetYTitle("#it{time} (ns)");
b2e375c7 1737 outputContainer->Add(fhTimePtPhotonNoCut);
1738
1739 fhTimePtPhotonSPD = new TH2F ("hTimePtPhoton_SPD","time of photon cluster vs pT of clusters, SPD cut", nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
11baad66 1740 fhTimePtPhotonSPD->SetXTitle("#it{p}_{T} (GeV/#it{c})");
184ca640 1741 fhTimePtPhotonSPD->SetYTitle("#it{time} (ns)");
b2e375c7 1742 outputContainer->Add(fhTimePtPhotonSPD);
bc41680b 1743
85c4406e 1744 fhTimeNPileUpVertSPD = new TH2F ("hTime_NPileUpVertSPD","time of cluster vs N pile-up SPD vertex", ntimebins,timemin,timemax,20,0,20);
2ad19c3d 1745 fhTimeNPileUpVertSPD->SetYTitle("# vertex ");
184ca640 1746 fhTimeNPileUpVertSPD->SetXTitle("#it{time} (ns)");
fad96885 1747 outputContainer->Add(fhTimeNPileUpVertSPD);
2ad19c3d 1748
85c4406e 1749 fhTimeNPileUpVertTrack = new TH2F ("hTime_NPileUpVertTracks","time of cluster vs N pile-up Tracks vertex", ntimebins,timemin,timemax, 20,0,20 );
2ad19c3d 1750 fhTimeNPileUpVertTrack->SetYTitle("# vertex ");
184ca640 1751 fhTimeNPileUpVertTrack->SetXTitle("#it{time} (ns)");
85c4406e 1752 outputContainer->Add(fhTimeNPileUpVertTrack);
2ad19c3d 1753
85c4406e 1754 fhPtPhotonNPileUpSPDVtx = new TH2F ("hPtPhoton_NPileUpVertSPD","pT of cluster vs N pile-up SPD vertex",
1755 nptbins,ptmin,ptmax,20,0,20);
1756 fhPtPhotonNPileUpSPDVtx->SetYTitle("# vertex ");
11baad66 1757 fhPtPhotonNPileUpSPDVtx->SetXTitle("#it{p}_{T} (GeV/#it{c})");
85c4406e 1758 outputContainer->Add(fhPtPhotonNPileUpSPDVtx);
0f7e7205 1759
85c4406e 1760 fhPtPhotonNPileUpTrkVtx = new TH2F ("hPtPhoton_NPileUpVertTracks","pT of cluster vs N pile-up Tracks vertex",
1761 nptbins,ptmin,ptmax, 20,0,20 );
1762 fhPtPhotonNPileUpTrkVtx->SetYTitle("# vertex ");
11baad66 1763 fhPtPhotonNPileUpTrkVtx->SetXTitle("#it{p}_{T} (GeV/#it{c})");
85c4406e 1764 outputContainer->Add(fhPtPhotonNPileUpTrkVtx);
0f7e7205 1765
85c4406e 1766 fhPtPhotonNPileUpSPDVtxTimeCut = new TH2F ("hPtPhoton_NPileUpVertSPD_TimeCut","pT of cluster vs N pile-up SPD vertex, |tof| < 25 ns",
1767 nptbins,ptmin,ptmax,20,0,20);
1768 fhPtPhotonNPileUpSPDVtxTimeCut->SetYTitle("# vertex ");
11baad66 1769 fhPtPhotonNPileUpSPDVtxTimeCut->SetXTitle("#it{p}_{T} (GeV/#it{c})");
85c4406e 1770 outputContainer->Add(fhPtPhotonNPileUpSPDVtxTimeCut);
0f7e7205 1771
85c4406e 1772 fhPtPhotonNPileUpTrkVtxTimeCut = new TH2F ("hPtPhoton_NPileUpVertTracks_TimeCut","pT of cluster vs N pile-up Tracks vertex, |tof| < 25 ns",
1773 nptbins,ptmin,ptmax, 20,0,20 );
1774 fhPtPhotonNPileUpTrkVtxTimeCut->SetYTitle("# vertex ");
11baad66 1775 fhPtPhotonNPileUpTrkVtxTimeCut->SetXTitle("#it{p}_{T} (GeV/#it{c})");
85c4406e 1776 outputContainer->Add(fhPtPhotonNPileUpTrkVtxTimeCut);
0f7e7205 1777
85c4406e 1778 fhPtPhotonNPileUpSPDVtxTimeCut2 = new TH2F ("hPtPhoton_NPileUpVertSPD_TimeCut2","pT of cluster vs N pile-up SPD vertex, -25 < tof < 75 ns",
1779 nptbins,ptmin,ptmax,20,0,20);
1780 fhPtPhotonNPileUpSPDVtxTimeCut2->SetYTitle("# vertex ");
11baad66 1781 fhPtPhotonNPileUpSPDVtxTimeCut2->SetXTitle("#it{p}_{T} (GeV/#it{c})");
85c4406e 1782 outputContainer->Add(fhPtPhotonNPileUpSPDVtxTimeCut2);
0f7e7205 1783
85c4406e 1784 fhPtPhotonNPileUpTrkVtxTimeCut2 = new TH2F ("hPtPhoton_NPileUpVertTracks_TimeCut2","pT of cluster vs N pile-up Tracks vertex, -25 < tof < 75 ns",
1785 nptbins,ptmin,ptmax, 20,0,20 );
1786 fhPtPhotonNPileUpTrkVtxTimeCut2->SetYTitle("# vertex ");
11baad66 1787 fhPtPhotonNPileUpTrkVtxTimeCut2->SetXTitle("#it{p}_{T} (GeV/#it{c})");
85c4406e 1788 outputContainer->Add(fhPtPhotonNPileUpTrkVtxTimeCut2);
1789
2ad19c3d 1790 }
bc41680b 1791
1792
2ad19c3d 1793
34c16486 1794 if(IsDataMC())
1795 {
f66d95af 1796 TString ptype[] = { "#gamma", "#gamma_{#pi decay}","#gamma_{other decay}", "#pi^{0}","#eta",
85c4406e 1797 "e^{#pm}","#gamma->e^{#pm}","hadron?","Anti-N","Anti-P",
1798 "#gamma_{prompt}","#gamma_{fragmentation}","#gamma_{ISR}","String" } ;
3d5d5078 1799
f66d95af 1800 TString pname[] = { "Photon","PhotonPi0Decay","PhotonOtherDecay","Pi0","Eta","Electron",
85c4406e 1801 "Conversion", "Hadron", "AntiNeutron","AntiProton",
1802 "PhotonPrompt","PhotonFragmentation","PhotonISR","String" } ;
521636d2 1803
34c16486 1804 for(Int_t i = 0; i < fNOriginHistograms; i++)
85c4406e 1805 {
3d5d5078 1806 fhMCE[i] = new TH1F(Form("hE_MC%s",pname[i].Data()),
85c4406e 1807 Form("cluster from %s : E ",ptype[i].Data()),
1808 nptbins,ptmin,ptmax);
184ca640 1809 fhMCE[i]->SetXTitle("#it{E} (GeV)");
85c4406e 1810 outputContainer->Add(fhMCE[i]) ;
521636d2 1811
4c8f7c2e 1812 fhMCPt[i] = new TH1F(Form("hPt_MC%s",pname[i].Data()),
11baad66 1813 Form("cluster from %s : #it{p}_{T} ",ptype[i].Data()),
85c4406e 1814 nptbins,ptmin,ptmax);
11baad66 1815 fhMCPt[i]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
4c8f7c2e 1816 outputContainer->Add(fhMCPt[i]) ;
521636d2 1817
4c8f7c2e 1818 fhMCEta[i] = new TH2F(Form("hEta_MC%s",pname[i].Data()),
85c4406e 1819 Form("cluster from %s : #eta ",ptype[i].Data()),
1820 nptbins,ptmin,ptmax,netabins,etamin,etamax);
4c8f7c2e 1821 fhMCEta[i]->SetYTitle("#eta");
184ca640 1822 fhMCEta[i]->SetXTitle("#it{E} (GeV)");
4c8f7c2e 1823 outputContainer->Add(fhMCEta[i]) ;
521636d2 1824
4c8f7c2e 1825 fhMCPhi[i] = new TH2F(Form("hPhi_MC%s",pname[i].Data()),
85c4406e 1826 Form("cluster from %s : #phi ",ptype[i].Data()),
1827 nptbins,ptmin,ptmax,nphibins,phimin,phimax);
4c8f7c2e 1828 fhMCPhi[i]->SetYTitle("#phi (rad)");
184ca640 1829 fhMCPhi[i]->SetXTitle("#it{E} (GeV)");
4c8f7c2e 1830 outputContainer->Add(fhMCPhi[i]) ;
1831
1832
d9105d92 1833 fhMCDeltaE[i] = new TH2F (Form("hDeltaE_MC%s",pname[i].Data()),
85c4406e 1834 Form("MC - Reco E from %s",pname[i].Data()),
1835 nptbins,ptmin,ptmax, 200,-50,50);
184ca640 1836 fhMCDeltaE[i]->SetYTitle("#Delta #it{E} (GeV)");
1837 fhMCDeltaE[i]->SetXTitle("#it{E} (GeV)");
4c8f7c2e 1838 outputContainer->Add(fhMCDeltaE[i]);
1839
d9105d92 1840 fhMCDeltaPt[i] = new TH2F (Form("hDeltaPt_MC%s",pname[i].Data()),
11baad66 1841 Form("MC - Reco #it{p}_{T} from %s",pname[i].Data()),
85c4406e 1842 nptbins,ptmin,ptmax, 200,-50,50);
184ca640 1843 fhMCDeltaPt[i]->SetXTitle("p_{T,rec} (GeV/#it{c})");
11baad66 1844 fhMCDeltaPt[i]->SetYTitle("#Delta #it{p}_{T} (GeV/#it{c})");
4c8f7c2e 1845 outputContainer->Add(fhMCDeltaPt[i]);
85c4406e 1846
4c8f7c2e 1847 fhMC2E[i] = new TH2F (Form("h2E_MC%s",pname[i].Data()),
85c4406e 1848 Form("E distribution, reconstructed vs generated from %s",pname[i].Data()),
1849 nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
184ca640 1850 fhMC2E[i]->SetXTitle("#it{E}_{rec} (GeV)");
1851 fhMC2E[i]->SetYTitle("#it{E}_{gen} (GeV)");
85c4406e 1852 outputContainer->Add(fhMC2E[i]);
4c8f7c2e 1853
1854 fhMC2Pt[i] = new TH2F (Form("h2Pt_MC%s",pname[i].Data()),
85c4406e 1855 Form("p_T distribution, reconstructed vs generated from %s",pname[i].Data()),
1856 nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
184ca640 1857 fhMC2Pt[i]->SetXTitle("p_{T,rec} (GeV/#it{c})");
1858 fhMC2Pt[i]->SetYTitle("p_{T,gen} (GeV/#it{c})");
4c8f7c2e 1859 outputContainer->Add(fhMC2Pt[i]);
1860
521636d2 1861
1862 }
3d5d5078 1863
f1c9c78f 1864 TString pptype[] = { "#gamma", "#gamma_{#pi decay}","#gamma_{other decay}",
85c4406e 1865 "#gamma_{prompt}","#gamma_{fragmentation}","#gamma_{ISR}"} ;
f66d95af 1866
f1c9c78f 1867 TString ppname[] = { "Photon","PhotonPi0Decay","PhotonOtherDecay",
85c4406e 1868 "PhotonPrompt","PhotonFragmentation","PhotonISR"} ;
f66d95af 1869
34c16486 1870 for(Int_t i = 0; i < fNPrimaryHistograms; i++)
85c4406e 1871 {
f66d95af 1872 fhEPrimMC[i] = new TH1F(Form("hEPrim_MC%s",ppname[i].Data()),
85c4406e 1873 Form("primary photon %s : E ",pptype[i].Data()),
1874 nptbins,ptmin,ptmax);
184ca640 1875 fhEPrimMC[i]->SetXTitle("#it{E} (GeV)");
85c4406e 1876 outputContainer->Add(fhEPrimMC[i]) ;
3d5d5078 1877
f66d95af 1878 fhPtPrimMC[i] = new TH1F(Form("hPtPrim_MC%s",ppname[i].Data()),
11baad66 1879 Form("primary photon %s : #it{p}_{T} ",pptype[i].Data()),
85c4406e 1880 nptbins,ptmin,ptmax);
11baad66 1881 fhPtPrimMC[i]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
3d5d5078 1882 outputContainer->Add(fhPtPrimMC[i]) ;
1883
f66d95af 1884 fhYPrimMC[i] = new TH2F(Form("hYPrim_MC%s",ppname[i].Data()),
85c4406e 1885 Form("primary photon %s : Rapidity ",pptype[i].Data()),
a421d7a3 1886 nptbins,ptmin,ptmax,200,-2,2);
3d5d5078 1887 fhYPrimMC[i]->SetYTitle("Rapidity");
184ca640 1888 fhYPrimMC[i]->SetXTitle("#it{E} (GeV)");
3d5d5078 1889 outputContainer->Add(fhYPrimMC[i]) ;
1890
4cf13296 1891 fhEtaPrimMC[i] = new TH2F(Form("hEtaPrim_MC%s",ppname[i].Data()),
1892 Form("primary photon %s : #eta",pptype[i].Data()),
a421d7a3 1893 nptbins,ptmin,ptmax,200,-2,2);
4cf13296 1894 fhEtaPrimMC[i]->SetYTitle("#eta");
184ca640 1895 fhEtaPrimMC[i]->SetXTitle("#it{E} (GeV)");
4cf13296 1896 outputContainer->Add(fhEtaPrimMC[i]) ;
1897
f66d95af 1898 fhPhiPrimMC[i] = new TH2F(Form("hPhiPrim_MC%s",ppname[i].Data()),
85c4406e 1899 Form("primary photon %s : #phi ",pptype[i].Data()),
a421d7a3 1900 nptbins,ptmin,ptmax,nphibins,0,TMath::TwoPi());
3d5d5078 1901 fhPhiPrimMC[i]->SetYTitle("#phi (rad)");
184ca640 1902 fhPhiPrimMC[i]->SetXTitle("#it{E} (GeV)");
3d5d5078 1903 outputContainer->Add(fhPhiPrimMC[i]) ;
85c4406e 1904
3d5d5078 1905
f66d95af 1906 fhEPrimMCAcc[i] = new TH1F(Form("hEPrimAcc_MC%s",ppname[i].Data()),
85c4406e 1907 Form("primary photon %s in acceptance: E ",pptype[i].Data()),
1908 nptbins,ptmin,ptmax);
184ca640 1909 fhEPrimMCAcc[i]->SetXTitle("#it{E} (GeV)");
85c4406e 1910 outputContainer->Add(fhEPrimMCAcc[i]) ;
3d5d5078 1911
f66d95af 1912 fhPtPrimMCAcc[i] = new TH1F(Form("hPtPrimAcc_MC%s",ppname[i].Data()),
11baad66 1913 Form("primary photon %s in acceptance: #it{p}_{T} ",pptype[i].Data()),
85c4406e 1914 nptbins,ptmin,ptmax);
11baad66 1915 fhPtPrimMCAcc[i]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
3d5d5078 1916 outputContainer->Add(fhPtPrimMCAcc[i]) ;
1917
f66d95af 1918 fhYPrimMCAcc[i] = new TH2F(Form("hYPrimAcc_MC%s",ppname[i].Data()),
85c4406e 1919 Form("primary photon %s in acceptance: Rapidity ",pptype[i].Data()),
1920 nptbins,ptmin,ptmax,100,-1,1);
3d5d5078 1921 fhYPrimMCAcc[i]->SetYTitle("Rapidity");
184ca640 1922 fhYPrimMCAcc[i]->SetXTitle("#it{E} (GeV)");
3d5d5078 1923 outputContainer->Add(fhYPrimMCAcc[i]) ;
4cf13296 1924
1925 fhEtaPrimMCAcc[i] = new TH2F(Form("hEtaPrimAcc_MC%s",ppname[i].Data()),
1926 Form("primary photon %s in acceptance: #eta ",pptype[i].Data()),
1927 nptbins,ptmin,ptmax,netabins,etamin,etamax);
1928 fhEtaPrimMCAcc[i]->SetYTitle("#eta");
184ca640 1929 fhEtaPrimMCAcc[i]->SetXTitle("#it{E} (GeV)");
667a3592 1930 outputContainer->Add(fhEtaPrimMCAcc[i]) ;
3d5d5078 1931
f66d95af 1932 fhPhiPrimMCAcc[i] = new TH2F(Form("hPhiPrimAcc_MC%s",ppname[i].Data()),
85c4406e 1933 Form("primary photon %s in acceptance: #phi ",pptype[i].Data()),
1934 nptbins,ptmin,ptmax,nphibins,phimin,phimax);
3d5d5078 1935 fhPhiPrimMCAcc[i]->SetYTitle("#phi (rad)");
184ca640 1936 fhPhiPrimMCAcc[i]->SetXTitle("#it{E} (GeV)");
3d5d5078 1937 outputContainer->Add(fhPhiPrimMCAcc[i]) ;
1938
1939 }
85c4406e 1940
34c16486 1941 if(fFillSSHistograms)
1942 {
85c4406e 1943 TString ptypess[] = { "#gamma","hadron?","#pi^{0}","#eta","#gamma->e^{#pm}","e^{#pm}"} ;
3d5d5078 1944
1945 TString pnamess[] = { "Photon","Hadron","Pi0","Eta","Conversion","Electron"} ;
1946
34c16486 1947 for(Int_t i = 0; i < 6; i++)
85c4406e 1948 {
3d5d5078 1949 fhMCELambda0[i] = new TH2F(Form("hELambda0_MC%s",pnamess[i].Data()),
1950 Form("cluster from %s : E vs #lambda_{0}^{2}",ptypess[i].Data()),
85c4406e 1951 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
3d5d5078 1952 fhMCELambda0[i]->SetYTitle("#lambda_{0}^{2}");
184ca640 1953 fhMCELambda0[i]->SetXTitle("#it{E} (GeV)");
85c4406e 1954 outputContainer->Add(fhMCELambda0[i]) ;
521636d2 1955
3d5d5078 1956 fhMCELambda1[i] = new TH2F(Form("hELambda1_MC%s",pnamess[i].Data()),
1957 Form("cluster from %s : E vs #lambda_{1}^{2}",ptypess[i].Data()),
85c4406e 1958 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
3d5d5078 1959 fhMCELambda1[i]->SetYTitle("#lambda_{1}^{2}");
184ca640 1960 fhMCELambda1[i]->SetXTitle("#it{E} (GeV)");
85c4406e 1961 outputContainer->Add(fhMCELambda1[i]) ;
34c16486 1962
3d5d5078 1963 fhMCEDispersion[i] = new TH2F(Form("hEDispersion_MC%s",pnamess[i].Data()),
1964 Form("cluster from %s : E vs dispersion^{2}",ptypess[i].Data()),
85c4406e 1965 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
3d5d5078 1966 fhMCEDispersion[i]->SetYTitle("D^{2}");
184ca640 1967 fhMCEDispersion[i]->SetXTitle("#it{E} (GeV)");
85c4406e 1968 outputContainer->Add(fhMCEDispersion[i]) ;
34c16486 1969
f66d95af 1970 fhMCNCellsE[i] = new TH2F (Form("hNCellsE_MC%s",pnamess[i].Data()),
85c4406e 1971 Form("# of cells in cluster from %s vs E of clusters",ptypess[i].Data()),
1972 nptbins,ptmin,ptmax, nbins,nmin,nmax);
184ca640 1973 fhMCNCellsE[i]->SetXTitle("#it{E} (GeV)");
f66d95af 1974 fhMCNCellsE[i]->SetYTitle("# of cells in cluster");
85c4406e 1975 outputContainer->Add(fhMCNCellsE[i]);
f66d95af 1976
1977 fhMCMaxCellDiffClusterE[i] = new TH2F (Form("hMaxCellDiffClusterE_MC%s",pnamess[i].Data()),
34c16486 1978 Form("energy vs difference of cluster energy from %s - max cell energy / cluster energy, good clusters",ptypess[i].Data()),
85c4406e 1979 nptbins,ptmin,ptmax, 500,0,1.);
184ca640 1980 fhMCMaxCellDiffClusterE[i]->SetXTitle("#it{E}_{cluster} (GeV) ");
1981 fhMCMaxCellDiffClusterE[i]->SetYTitle("(#it{E}_{cluster} - #it{E}_{cell max})/ #it{E}_{cluster}");
85c4406e 1982 outputContainer->Add(fhMCMaxCellDiffClusterE[i]);
f66d95af 1983
764ab1f4 1984 if(!fFillOnlySimpleSSHisto)
34c16486 1985 {
764ab1f4 1986 fhMCLambda0vsClusterMaxCellDiffE0[i] = new TH2F(Form("hLambda0vsClusterMaxCellDiffE0_MC%s",pnamess[i].Data()),
1987 Form("cluster from %s : #lambda^{2}_{0} vs fraction of energy carried by max cell, E < 2 GeV",ptypess[i].Data()),
85c4406e 1988 ssbins,ssmin,ssmax,500,0,1.);
764ab1f4 1989 fhMCLambda0vsClusterMaxCellDiffE0[i]->SetXTitle("#lambda_{0}^{2}");
184ca640 1990 fhMCLambda0vsClusterMaxCellDiffE0[i]->SetYTitle("(#it{E}_{cluster} - #it{E}_{cell max})/ #it{E}_{cluster}");
85c4406e 1991 outputContainer->Add(fhMCLambda0vsClusterMaxCellDiffE0[i]) ;
764ab1f4 1992
1993 fhMCLambda0vsClusterMaxCellDiffE2[i] = new TH2F(Form("hLambda0vsClusterMaxCellDiffE2_MC%s",pnamess[i].Data()),
1994 Form("cluster from %s : #lambda^{2}_{0} vs fraction of energy carried by max cell, 2< E < 6 GeV",ptypess[i].Data()),
85c4406e 1995 ssbins,ssmin,ssmax,500,0,1.);
764ab1f4 1996 fhMCLambda0vsClusterMaxCellDiffE2[i]->SetXTitle("#lambda_{0}^{2}");
184ca640 1997 fhMCLambda0vsClusterMaxCellDiffE2[i]->SetYTitle("(#it{E}_{cluster} - #it{E}_{cell max})/ #it{E}_{cluster}");
85c4406e 1998 outputContainer->Add(fhMCLambda0vsClusterMaxCellDiffE2[i]) ;
34c16486 1999
764ab1f4 2000 fhMCLambda0vsClusterMaxCellDiffE6[i] = new TH2F(Form("hLambda0vsClusterMaxCellDiffE6_MC%s",pnamess[i].Data()),
184ca640 2001 Form("cluster from %s : #lambda^{2}_{0} vs fraction of energy carried by max cell, #it{E} > 6 GeV",ptypess[i].Data()),
85c4406e 2002 ssbins,ssmin,ssmax,500,0,1.);
764ab1f4 2003 fhMCLambda0vsClusterMaxCellDiffE6[i]->SetXTitle("#lambda_{0}^{2}");
184ca640 2004 fhMCLambda0vsClusterMaxCellDiffE6[i]->SetYTitle("(#it{E}_{cluster} - #it{E}_{cell max})/ #it{E}_{cluster}");
85c4406e 2005 outputContainer->Add(fhMCLambda0vsClusterMaxCellDiffE6[i]) ;
34c16486 2006
764ab1f4 2007 fhMCNCellsvsClusterMaxCellDiffE0[i] = new TH2F(Form("hNCellsvsClusterMaxCellDiffE0_MC%s",pnamess[i].Data()),
2008 Form("cluster from %s : N cells in cluster vs fraction of energy carried by max cell, E < 2 GeV",ptypess[i].Data()),
85c4406e 2009 nbins/5,nmin,nmax/5,500,0,1.);
764ab1f4 2010 fhMCNCellsvsClusterMaxCellDiffE0[i]->SetXTitle("N cells in cluster");
184ca640 2011 fhMCNCellsvsClusterMaxCellDiffE0[i]->SetYTitle("(#it{E}_{cluster} - #it{E}_{cell max})/ #it{E}_{cluster}");
85c4406e 2012 outputContainer->Add(fhMCNCellsvsClusterMaxCellDiffE0[i]) ;
34c16486 2013
764ab1f4 2014 fhMCNCellsvsClusterMaxCellDiffE2[i] = new TH2F(Form("hNCellsvsClusterMaxCellDiffE2_MC%s",pnamess[i].Data()),
2015 Form("cluster from %s : N cells in cluster vs fraction of energy carried by max cell, 2< E < 6 GeV",ptypess[i].Data()),
85c4406e 2016 nbins/5,nmin,nmax/5,500,0,1.);
764ab1f4 2017 fhMCNCellsvsClusterMaxCellDiffE2[i]->SetXTitle("N cells in cluster");
184ca640 2018 fhMCNCellsvsClusterMaxCellDiffE2[i]->SetYTitle("(#it{E}_{cluster} - #it{E}_{cell max})/ #it{E}_{cluster}");
85c4406e 2019 outputContainer->Add(fhMCNCellsvsClusterMaxCellDiffE2[i]) ;
34c16486 2020
764ab1f4 2021 fhMCNCellsvsClusterMaxCellDiffE6[i] = new TH2F(Form("hNCellsvsClusterMaxCellDiffE6_MC%s",pnamess[i].Data()),
184ca640 2022 Form("cluster from %s : N cells in cluster vs fraction of energy carried by max cell, #it{E} > 6 GeV",ptypess[i].Data()),
85c4406e 2023 nbins/5,nmin,nmax/5,500,0,1.);
764ab1f4 2024 fhMCNCellsvsClusterMaxCellDiffE6[i]->SetXTitle("N cells in cluster");
184ca640 2025 fhMCNCellsvsClusterMaxCellDiffE6[i]->SetYTitle("#it{E} (GeV)");
85c4406e 2026 outputContainer->Add(fhMCNCellsvsClusterMaxCellDiffE6[i]) ;
34c16486 2027
764ab1f4 2028 if(fCalorimeter=="EMCAL")
34c16486 2029 {
764ab1f4 2030 fhMCEDispEta[i] = new TH2F (Form("hEDispEtaE_MC%s",pnamess[i].Data()),
2031 Form("cluster from %s : #sigma^{2}_{#eta #eta} = #Sigma w_{i}(#eta_{i} - <#eta>)^{2}/ #Sigma w_{i} vs E",ptypess[i].Data()),
85c4406e 2032 nptbins,ptmin,ptmax, ssbins,ssmin,ssmax);
184ca640 2033 fhMCEDispEta[i]->SetXTitle("#it{E} (GeV)");
764ab1f4 2034 fhMCEDispEta[i]->SetYTitle("#sigma^{2}_{#eta #eta}");
85c4406e 2035 outputContainer->Add(fhMCEDispEta[i]);
764ab1f4 2036
2037 fhMCEDispPhi[i] = new TH2F (Form("hEDispPhiE_MC%s",pnamess[i].Data()),
2038 Form("cluster from %s : #sigma^{2}_{#phi #phi} = #Sigma w_{i}(#phi_{i} - <#phi>)^{2} / #Sigma w_{i} vs E",ptypess[i].Data()),
85c4406e 2039 nptbins,ptmin,ptmax, ssbins,ssmin,ssmax);
184ca640 2040 fhMCEDispPhi[i]->SetXTitle("#it{E} (GeV)");
764ab1f4 2041 fhMCEDispPhi[i]->SetYTitle("#sigma^{2}_{#phi #phi}");
85c4406e 2042 outputContainer->Add(fhMCEDispPhi[i]);
764ab1f4 2043
2044 fhMCESumEtaPhi[i] = new TH2F (Form("hESumEtaPhiE_MC%s",pnamess[i].Data()),
85c4406e 2045 Form("cluster from %s : #delta^{2}_{#eta #phi} = #Sigma w_{i}(#phi_{i} #eta_{i} ) / #Sigma w_{i} - <#phi><#eta> vs E",ptypess[i].Data()),
2046 nptbins,ptmin,ptmax, 2*ssbins,-ssmax,ssmax);
184ca640 2047 fhMCESumEtaPhi[i]->SetXTitle("#it{E} (GeV)");
764ab1f4 2048 fhMCESumEtaPhi[i]->SetYTitle("#delta^{2}_{#eta #phi}");
2049 outputContainer->Add(fhMCESumEtaPhi[i]);
2050
2051 fhMCEDispEtaPhiDiff[i] = new TH2F (Form("hEDispEtaPhiDiffE_MC%s",pnamess[i].Data()),
85c4406e 2052 Form("cluster from %s : #sigma^{2}_{#phi #phi} - #sigma^{2}_{#eta #eta} vs E",ptypess[i].Data()),
2053 nptbins,ptmin,ptmax,200,-10,10);
184ca640 2054 fhMCEDispEtaPhiDiff[i]->SetXTitle("#it{E} (GeV)");
764ab1f4 2055 fhMCEDispEtaPhiDiff[i]->SetYTitle("#sigma^{2}_{#phi #phi}-#sigma^{2}_{#eta #eta}");
85c4406e 2056 outputContainer->Add(fhMCEDispEtaPhiDiff[i]);
764ab1f4 2057
2058 fhMCESphericity[i] = new TH2F (Form("hESphericity_MC%s",pnamess[i].Data()),
85c4406e 2059 Form("cluster from %s : (#sigma^{2}_{#phi #phi} - #sigma^{2}_{#eta #eta}) / (#sigma^{2}_{#eta #eta} + #sigma^{2}_{#phi #phi}) vs E",ptypess[i].Data()),
2060 nptbins,ptmin,ptmax, 200,-1,1);
184ca640 2061 fhMCESphericity[i]->SetXTitle("#it{E} (GeV)");
764ab1f4 2062 fhMCESphericity[i]->SetYTitle("s = (#sigma^{2}_{#phi #phi} - #sigma^{2}_{#eta #eta}) / (#sigma^{2}_{#eta #eta} + #sigma^{2}_{#phi #phi})");
2063 outputContainer->Add(fhMCESphericity[i]);
2064
2065 for(Int_t ie = 0; ie < 7; ie++)
2066 {
2067 fhMCDispEtaDispPhi[ie][i] = new TH2F (Form("hMCDispEtaDispPhi_EBin%d_MC%s",ie,pnamess[i].Data()),
85c4406e 2068 Form("cluster from %s : #sigma^{2}_{#phi #phi} vs #sigma^{2}_{#eta #eta} for %d < E < %d GeV",pnamess[i].Data(),bin[ie],bin[ie+1]),
2069 ssbins,ssmin,ssmax , ssbins,ssmin,ssmax);
764ab1f4 2070 fhMCDispEtaDispPhi[ie][i]->SetXTitle("#sigma^{2}_{#eta #eta}");
2071 fhMCDispEtaDispPhi[ie][i]->SetYTitle("#sigma^{2}_{#phi #phi}");
85c4406e 2072 outputContainer->Add(fhMCDispEtaDispPhi[ie][i]);
764ab1f4 2073
2074 fhMCLambda0DispEta[ie][i] = new TH2F (Form("hMCLambda0DispEta_EBin%d_MC%s",ie,pnamess[i].Data()),
85c4406e 2075 Form("cluster from %s : #lambda^{2}_{0} vs #sigma^{2}_{#eta #eta} for %d < E < %d GeV",pnamess[i].Data(),bin[ie],bin[ie+1]),
2076 ssbins,ssmin,ssmax , ssbins,ssmin,ssmax);
764ab1f4 2077 fhMCLambda0DispEta[ie][i]->SetXTitle("#lambda^{2}_{0}");
2078 fhMCLambda0DispEta[ie][i]->SetYTitle("#sigma^{2}_{#phi #phi}");
85c4406e 2079 outputContainer->Add(fhMCLambda0DispEta[ie][i]);
764ab1f4 2080
2081 fhMCLambda0DispPhi[ie][i] = new TH2F (Form("hMCLambda0DispPhi_EBin%d_MC%s",ie,pnamess[i].Data()),
85c4406e 2082 Form("cluster from %s :#lambda^{2}_{0} vs #sigma^{2}_{#phi #phi} for %d < E < %d GeV",pnamess[i].Data(),bin[ie],bin[ie+1]),
2083 ssbins,ssmin,ssmax , ssbins,ssmin,ssmax);
764ab1f4 2084 fhMCLambda0DispPhi[ie][i]->SetXTitle("#lambda^{2}_{0}");
2085 fhMCLambda0DispPhi[ie][i]->SetYTitle("#sigma^{2}_{#phi #phi}");
85c4406e 2086 outputContainer->Add(fhMCLambda0DispPhi[ie][i]);
764ab1f4 2087 }
34c16486 2088 }
34c16486 2089 }
85c4406e 2090 }// loop
3d5d5078 2091
2092 if(!GetReader()->IsEmbeddedClusterSelectionOn())
2093 {
2094 fhMCPhotonELambda0NoOverlap = new TH2F("hELambda0_MCPhoton_NoOverlap",
2095 "cluster from Photon : E vs #lambda_{0}^{2}",
85c4406e 2096 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
3d5d5078 2097 fhMCPhotonELambda0NoOverlap->SetYTitle("#lambda_{0}^{2}");
184ca640 2098 fhMCPhotonELambda0NoOverlap->SetXTitle("#it{E} (GeV)");
85c4406e 2099 outputContainer->Add(fhMCPhotonELambda0NoOverlap) ;
3d5d5078 2100
3d5d5078 2101 fhMCPhotonELambda0TwoOverlap = new TH2F("hELambda0_MCPhoton_TwoOverlap",
2102 "cluster from Photon : E vs #lambda_{0}^{2}",
85c4406e 2103 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
3d5d5078 2104 fhMCPhotonELambda0TwoOverlap->SetYTitle("#lambda_{0}^{2}");
184ca640 2105 fhMCPhotonELambda0TwoOverlap->SetXTitle("#it{E} (GeV)");
85c4406e 2106 outputContainer->Add(fhMCPhotonELambda0TwoOverlap) ;
3d5d5078 2107
3d5d5078 2108 fhMCPhotonELambda0NOverlap = new TH2F("hELambda0_MCPhoton_NOverlap",
2109 "cluster from Photon : E vs #lambda_{0}^{2}",
85c4406e 2110 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
3d5d5078 2111 fhMCPhotonELambda0NOverlap->SetYTitle("#lambda_{0}^{2}");
184ca640 2112 fhMCPhotonELambda0NOverlap->SetXTitle("#it{E} (GeV)");
85c4406e 2113 outputContainer->Add(fhMCPhotonELambda0NOverlap) ;
521636d2 2114
85c4406e 2115 } //No embedding
3d5d5078 2116
3d5d5078 2117 if(GetReader()->IsEmbeddedClusterSelectionOn())
2118 {
2119
2120 fhEmbeddedSignalFractionEnergy = new TH2F("hEmbeddedSignal_FractionEnergy",
34c16486 2121 "Energy Fraction of embedded signal versus cluster energy",
85c4406e 2122 nptbins,ptmin,ptmax,100,0.,1.);
3d5d5078 2123 fhEmbeddedSignalFractionEnergy->SetYTitle("Fraction");
184ca640 2124 fhEmbeddedSignalFractionEnergy->SetXTitle("#it{E} (GeV)");
85c4406e 2125 outputContainer->Add(fhEmbeddedSignalFractionEnergy) ;
3d5d5078 2126
2127 fhEmbedPhotonELambda0FullSignal = new TH2F("hELambda0_EmbedPhoton_FullSignal",
34c16486 2128 "cluster from Photon embedded with more than 90% energy in cluster : E vs #lambda_{0}^{2}",
85c4406e 2129 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
3d5d5078 2130 fhEmbedPhotonELambda0FullSignal->SetYTitle("#lambda_{0}^{2}");
184ca640 2131 fhEmbedPhotonELambda0FullSignal->SetXTitle("#it{E} (GeV)");
85c4406e 2132 outputContainer->Add(fhEmbedPhotonELambda0FullSignal) ;
34c16486 2133
3d5d5078 2134 fhEmbedPhotonELambda0MostlySignal = new TH2F("hELambda0_EmbedPhoton_MostlySignal",
34c16486 2135 "cluster from Photon embedded with 50% to 90% energy in cluster : E vs #lambda_{0}^{2}",
85c4406e 2136 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
3d5d5078 2137 fhEmbedPhotonELambda0MostlySignal->SetYTitle("#lambda_{0}^{2}");
184ca640 2138 fhEmbedPhotonELambda0MostlySignal->SetXTitle("#it{E} (GeV)");
85c4406e 2139 outputContainer->Add(fhEmbedPhotonELambda0MostlySignal) ;
3d5d5078 2140
3d5d5078 2141 fhEmbedPhotonELambda0MostlyBkg = new TH2F("hELambda0_EmbedPhoton_MostlyBkg",
34c16486 2142 "cluster from Photon embedded with 10% to 50% energy in cluster : E vs #lambda_{0}^{2}",
85c4406e 2143 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
3d5d5078 2144 fhEmbedPhotonELambda0MostlyBkg->SetYTitle("#lambda_{0}^{2}");
184ca640 2145 fhEmbedPhotonELambda0MostlyBkg->SetXTitle("#it{E} (GeV)");
85c4406e 2146 outputContainer->Add(fhEmbedPhotonELambda0MostlyBkg) ;
34c16486 2147
3d5d5078 2148 fhEmbedPhotonELambda0FullBkg = new TH2F("hELambda0_EmbedPhoton_FullBkg",
34c16486 2149 "cluster from Photonm embedded with 0% to 10% energy in cluster : E vs #lambda_{0}^{2}",
85c4406e 2150 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
3d5d5078 2151 fhEmbedPhotonELambda0FullBkg->SetYTitle("#lambda_{0}^{2}");
184ca640 2152 fhEmbedPhotonELambda0FullBkg->SetXTitle("#it{E} (GeV)");
85c4406e 2153 outputContainer->Add(fhEmbedPhotonELambda0FullBkg) ;
3d5d5078 2154
3d5d5078 2155 fhEmbedPi0ELambda0FullSignal = new TH2F("hELambda0_EmbedPi0_FullSignal",
34c16486 2156 "cluster from Pi0 embedded with more than 90% energy in cluster : E vs #lambda_{0}^{2}",
85c4406e 2157 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
3d5d5078 2158 fhEmbedPi0ELambda0FullSignal->SetYTitle("#lambda_{0}^{2}");
184ca640 2159 fhEmbedPi0ELambda0FullSignal->SetXTitle("#it{E} (GeV)");
85c4406e 2160 outputContainer->Add(fhEmbedPi0ELambda0FullSignal) ;
34c16486 2161
3d5d5078 2162 fhEmbedPi0ELambda0MostlySignal = new TH2F("hELambda0_EmbedPi0_MostlySignal",
34c16486 2163 "cluster from Pi0 embedded with 50% to 90% energy in cluster : E vs #lambda_{0}^{2}",
85c4406e 2164 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
3d5d5078 2165 fhEmbedPi0ELambda0MostlySignal->SetYTitle("#lambda_{0}^{2}");
184ca640 2166 fhEmbedPi0ELambda0MostlySignal->SetXTitle("#it{E} (GeV)");
85c4406e 2167 outputContainer->Add(fhEmbedPi0ELambda0MostlySignal) ;
3d5d5078 2168
3d5d5078 2169 fhEmbedPi0ELambda0MostlyBkg = new TH2F("hELambda0_EmbedPi0_MostlyBkg",
34c16486 2170 "cluster from Pi0 embedded with 10% to 50% energy in cluster : E vs #lambda_{0}^{2}",
85c4406e 2171 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
3d5d5078 2172 fhEmbedPi0ELambda0MostlyBkg->SetYTitle("#lambda_{0}^{2}");
184ca640 2173 fhEmbedPi0ELambda0MostlyBkg->SetXTitle("#it{E} (GeV)");
85c4406e 2174 outputContainer->Add(fhEmbedPi0ELambda0MostlyBkg) ;
3d5d5078 2175
3d5d5078 2176 fhEmbedPi0ELambda0FullBkg = new TH2F("hELambda0_EmbedPi0_FullBkg",
34c16486 2177 "cluster from Pi0 embedded with 0% to 10% energy in cluster : E vs #lambda_{0}^{2}",
85c4406e 2178 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
3d5d5078 2179 fhEmbedPi0ELambda0FullBkg->SetYTitle("#lambda_{0}^{2}");
184ca640 2180 fhEmbedPi0ELambda0FullBkg->SetXTitle("#it{E} (GeV)");
85c4406e 2181 outputContainer->Add(fhEmbedPi0ELambda0FullBkg) ;
34c16486 2182
3d5d5078 2183 }// embedded histograms
2184
521636d2 2185
2186 }// Fill SS MC histograms
2187
477d6cee 2188 }//Histos with MC
1035a8d9 2189
477d6cee 2190 return outputContainer ;
2191
1c5acb87 2192}
2193
34c16486 2194//_______________________
6639984f 2195void AliAnaPhoton::Init()
2196{
2197
2198 //Init
2199 //Do some checks
34c16486 2200 if(fCalorimeter == "PHOS" && !GetReader()->IsPHOSSwitchedOn() && NewOutputAOD())
2201 {
591cc579 2202 printf("AliAnaPhoton::Init() - !!STOP: You want to use PHOS in analysis but it is not read!! \n!!Check the configuration file!!\n");
6639984f 2203 abort();
2204 }
34c16486 2205 else if(fCalorimeter == "EMCAL" && !GetReader()->IsEMCALSwitchedOn() && NewOutputAOD())
2206 {
591cc579 2207 printf("AliAnaPhoton::Init() - !!STOP: You want to use EMCAL in analysis but it is not read!! \n!!Check the configuration file!!\n");
6639984f 2208 abort();
2209 }
2210
49b5c49b 2211 if(GetReader()->GetDataType() == AliCaloTrackReader::kMC) GetCaloPID()->SwitchOnBayesian();
2212
6639984f 2213}
2214
6639984f 2215//____________________________________________________________________________
1c5acb87 2216void AliAnaPhoton::InitParameters()
2217{
2218
2219 //Initialize the parameters of the analysis.
a3aebfff 2220 AddToHistogramsName("AnaPhoton_");
521636d2 2221
6175da48 2222 fCalorimeter = "EMCAL" ;
2223 fMinDist = 2.;
2224 fMinDist2 = 4.;
2225 fMinDist3 = 5.;
1e86c71e 2226
caa8a222 2227 fTimeCutMin =-1000000;
2228 fTimeCutMax = 1000000;
6175da48 2229 fNCellsCut = 0;
2ac125bf 2230
1e86c71e 2231 fRejectTrackMatch = kTRUE ;
1e86c71e 2232
1c5acb87 2233}
2234
2235//__________________________________________________________________
85c4406e 2236void AliAnaPhoton::MakeAnalysisFillAOD()
1c5acb87 2237{
f8006433 2238 //Do photon analysis and fill aods
f37fa8d2 2239
85c4406e 2240 //Get the vertex
5025c139 2241 Double_t v[3] = {0,0,0}; //vertex ;
2242 GetReader()->GetVertex(v);
f8006433 2243
f37fa8d2 2244 //Select the Calorimeter of the photon
85c4406e 2245 TObjArray * pl = 0x0;
2246 AliVCaloCells* cells = 0;
71e3889f 2247 if (fCalorimeter == "PHOS" )
2248 {
2249 pl = GetPHOSClusters();
2250 cells = GetPHOSCells();
2251 }
477d6cee 2252 else if (fCalorimeter == "EMCAL")
71e3889f 2253 {
2254 pl = GetEMCALClusters();
2255 cells = GetEMCALCells();
2256 }
5ae09196 2257
85c4406e 2258 if(!pl)
34c16486 2259 {
5ae09196 2260 Info("MakeAnalysisFillAOD","TObjArray with %s clusters is NULL!\n",fCalorimeter.Data());
2261 return;
2262 }
521636d2 2263
58ea8ce5 2264 TLorentzVector mom;
2265
fc195fd0 2266 // Loop on raw clusters before filtering in the reader and fill control histogram
34c16486 2267 if((GetReader()->GetEMCALClusterListName()=="" && fCalorimeter=="EMCAL") || fCalorimeter=="PHOS")
2268 {
2269 for(Int_t iclus = 0; iclus < GetReader()->GetInputEvent()->GetNumberOfCaloClusters(); iclus++ )
2270 {
fc195fd0 2271 AliVCluster * clus = GetReader()->GetInputEvent()->GetCaloCluster(iclus);
58ea8ce5 2272 if (fCalorimeter == "PHOS" && clus->IsPHOS() && clus->E() > GetReader()->GetPHOSPtMin() )
2273 {
2274 fhClusterCutsE [0]->Fill(clus->E());
2275
2276 clus->GetMomentum(mom,GetVertex(0)) ;
2277 fhClusterCutsPt[0]->Fill(mom.Pt());
2278 }
2279 else if(fCalorimeter == "EMCAL" && clus->IsEMCAL() && clus->E() > GetReader()->GetEMCALPtMin())
2280 {
2281 fhClusterCutsE [0]->Fill(clus->E());
2282
2283 clus->GetMomentum(mom,GetVertex(0)) ;
2284 fhClusterCutsPt[0]->Fill(mom.Pt());
2285 }
fc195fd0 2286 }
2287 }
85c4406e 2288 else
34c16486 2289 { // reclusterized
fc195fd0 2290 TClonesArray * clusterList = 0;
85c4406e 2291
7d650cb7 2292 if(GetReader()->GetInputEvent()->FindListObject(GetReader()->GetEMCALClusterListName()))
2293 clusterList = dynamic_cast<TClonesArray*> (GetReader()->GetInputEvent()->FindListObject(GetReader()->GetEMCALClusterListName()));
2294 else if(GetReader()->GetOutputEvent())
4a9e1073 2295 clusterList = dynamic_cast<TClonesArray*> (GetReader()->GetOutputEvent()->FindListObject(GetReader()->GetEMCALClusterListName()));
7d650cb7 2296
34c16486 2297 if(clusterList)
2298 {
fc195fd0 2299 Int_t nclusters = clusterList->GetEntriesFast();
85c4406e 2300 for (Int_t iclus = 0; iclus < nclusters; iclus++)
34c16486 2301 {
85c4406e 2302 AliVCluster * clus = dynamic_cast<AliVCluster*> (clusterList->At(iclus));
58ea8ce5 2303 if(clus)
2304 {
2305 fhClusterCutsE [0]->Fill(clus->E());
2306
2307 clus->GetMomentum(mom,GetVertex(0)) ;
2308 fhClusterCutsPt[0]->Fill(mom.Pt());
2309 }
6265ad55 2310 }
fc195fd0 2311 }
2312 }
fc195fd0 2313
6175da48 2314 //Init arrays, variables, get number of clusters
58ea8ce5 2315 TLorentzVector mom2 ;
1e86c71e 2316 Int_t nCaloClusters = pl->GetEntriesFast();
20218aea 2317
6175da48 2318 if(GetDebug() > 0) printf("AliAnaPhoton::MakeAnalysisFillAOD() - input %s cluster entries %d\n", fCalorimeter.Data(), nCaloClusters);
521636d2 2319
6175da48 2320 //----------------------------------------------------
2321 // Fill AOD with PHOS/EMCAL AliAODPWG4Particle objects
2322 //----------------------------------------------------
2323 // Loop on clusters
34c16486 2324 for(Int_t icalo = 0; icalo < nCaloClusters; icalo++)
85c4406e 2325 {
2326 AliVCluster * calo = (AliVCluster*) (pl->At(icalo));
0ae57829 2327 //printf("calo %d, %f\n",icalo,calo->E());
521636d2 2328
f8006433 2329 //Get the index where the cluster comes, to retrieve the corresponding vertex
85c4406e 2330 Int_t evtIndex = 0 ;
2331 if (GetMixedEvent())
34c16486 2332 {
85c4406e 2333 evtIndex=GetMixedEvent()->EventIndexForCaloCluster(calo->GetID()) ;
5025c139 2334 //Get the vertex and check it is not too large in z
96539743 2335 if(TMath::Abs(GetVertex(evtIndex)[2])> GetZvertexCut()) continue;
c8fe2783 2336 }
521636d2 2337
85c4406e 2338 //Cluster selection, not charged, with photon id and in fiducial cut
34c16486 2339 if(GetReader()->GetDataType() != AliCaloTrackReader::kMC)
2340 {
afb3af8a 2341 calo->GetMomentum(mom,GetVertex(evtIndex)) ;
2342 }//Assume that come from vertex in straight line
34c16486 2343 else
2344 {
f8006433 2345 Double_t vertex[]={0,0,0};
2346 calo->GetMomentum(mom,vertex) ;
2347 }
c2a62a94 2348
bc41680b 2349 //-----------------------------
6175da48 2350 // Cluster selection
bc41680b 2351 //-----------------------------
9e51e29a 2352 Int_t nMaxima = GetCaloUtils()->GetNumberOfLocalMaxima(calo, cells); // NLM
2353 if(!ClusterSelected(calo,mom,nMaxima)) continue;
85c4406e 2354
6175da48 2355 //----------------------------
2356 //Create AOD for analysis
2357 //----------------------------
2358 AliAODPWG4Particle aodph = AliAODPWG4Particle(mom);
2359
2360 //...............................................
85c4406e 2361 //Set the indeces of the original caloclusters (MC, ID), and calorimeter
6175da48 2362 Int_t label = calo->GetLabel();
2363 aodph.SetLabel(label);
6175da48 2364 aodph.SetCaloLabel(calo->GetID(),-1);
2365 aodph.SetDetector(fCalorimeter);
c4a7d28a 2366 //printf("Index %d, Id %d, iaod %d\n",icalo, calo->GetID(),GetOutputAODBranch()->GetEntriesFast());
521636d2 2367
6175da48 2368 //...............................................
2369 //Set bad channel distance bit
c4a7d28a 2370 Double_t distBad=calo->GetDistanceToBadChannel() ; //Distance to bad channel
f37fa8d2 2371 if (distBad > fMinDist3) aodph.SetDistToBad(2) ;
85c4406e 2372 else if(distBad > fMinDist2) aodph.SetDistToBad(1) ;
f37fa8d2 2373 else aodph.SetDistToBad(0) ;
af7b3903 2374 //printf("DistBad %f Bit %d\n",distBad, aodph.DistToBad());
c8fe2783 2375
bc41680b 2376 //-------------------------------------
8d6b7f60 2377 // Play with the MC stack if available
bc41680b 2378 //-------------------------------------
8d6b7f60 2379
2380 //Check origin of the candidates
2381 Int_t tag = -1;
2382
34c16486 2383 if(IsDataMC())
2384 {
2644ead9 2385 tag = GetMCAnalysisUtils()->CheckOrigin(calo->GetLabels(),calo->GetNLabels(),GetReader());
8d6b7f60 2386 aodph.SetTag(tag);
2387
2388 if(GetDebug() > 0)
2389 printf("AliAnaPhoton::MakeAnalysisFillAOD() - Origin of candidate, bit map %d\n",aodph.GetTag());
85c4406e 2390 }//Work with stack also
2391
bc41680b 2392 //--------------------------------------------------------
521636d2 2393 //Fill some shower shape histograms before PID is applied
bc41680b 2394 //--------------------------------------------------------
521636d2 2395
8d6b7f60 2396 FillShowerShapeHistograms(calo,tag);
6175da48 2397
2398 //-------------------------------------
f37fa8d2 2399 //PID selection or bit setting
6175da48 2400 //-------------------------------------
49b5c49b 2401
6175da48 2402 //...............................................
2403 // Data, PID check on
3c1d9afb 2404 if(IsCaloPIDOn())
2405 {
49b5c49b 2406 // Get most probable PID, 2 options check bayesian PID weights or redo PID
2407 // By default, redo PID
09273901 2408
3c1d9afb 2409 aodph.SetIdentifiedParticleType(GetCaloPID()->GetIdentifiedParticleType(calo));
477d6cee 2410
21a4b1c0 2411 if(GetDebug() > 1) printf("AliAnaPhoton::MakeAnalysisFillAOD() - PDG of identified particle %d\n",aodph.GetIdentifiedParticleType());
477d6cee 2412
f37fa8d2 2413 //If cluster does not pass pid, not photon, skip it.
85c4406e 2414 if(aodph.GetIdentifiedParticleType() != AliCaloPID::kPhoton) continue ;
477d6cee 2415
2416 }
85c4406e 2417
6175da48 2418 //...............................................
2419 // Data, PID check off
3c1d9afb 2420 else
2421 {
f37fa8d2 2422 //Set PID bits for later selection (AliAnaPi0 for example)
49b5c49b 2423 //GetIdentifiedParticleType already called in SetPIDBits.
2424
3c1d9afb 2425 GetCaloPID()->SetPIDBits(calo,&aodph, GetCaloUtils(),GetReader()->GetInputEvent());
49b5c49b 2426
85c4406e 2427 if(GetDebug() > 1) printf("AliAnaPhoton::MakeAnalysisFillAOD() - PID Bits set \n");
477d6cee 2428 }
2429
3c1d9afb 2430 if(GetDebug() > 1) printf("AliAnaPhoton::MakeAnalysisFillAOD() - Photon selection cuts passed: pT %3.2f, pdg %d\n",
2431 aodph.Pt(), aodph.GetIdentifiedParticleType());
09273901 2432
58ea8ce5 2433 fhClusterCutsE [9]->Fill(calo->E());
2434 fhClusterCutsPt[9]->Fill(mom.Pt());
85c4406e 2435
6df33fcb 2436 Int_t nSM = GetModuleNumber(calo);
2437 if(nSM < GetCaloUtils()->GetNumberOfSuperModulesUsed() && nSM >=0)
2438 {
2439 fhEPhotonSM ->Fill(mom.E (),nSM);
2440 fhPtPhotonSM->Fill(mom.Pt(),nSM);
2441 }
2442
9e51e29a 2443 fhNLocMax->Fill(calo->E(),nMaxima);
85c4406e 2444
09273901 2445 // Matching after cuts
bc41680b 2446 if( fFillTMHisto ) FillTrackMatchingResidualHistograms(calo,1);
09273901 2447
2ad19c3d 2448 // Fill histograms to undertand pile-up before other cuts applied
2449 // Remember to relax time cuts in the reader
bc41680b 2450 if( fFillPileUpHistograms ) FillPileUpHistograms(calo,cells);
2ad19c3d 2451
5c46c992 2452 // Add number of local maxima to AOD, method name in AOD to be FIXED
9e51e29a 2453 aodph.SetFiducialArea(nMaxima);
5c46c992 2454
f37fa8d2 2455 //Add AOD with photon object to aod branch
477d6cee 2456 AddAODParticle(aodph);
2457
2458 }//loop
85c4406e 2459
2460 if(GetDebug() > 1) printf("AliAnaPhoton::MakeAnalysisFillAOD() End fill AODs, with %d entries \n",GetOutputAODBranch()->GetEntriesFast());
477d6cee 2461
1c5acb87 2462}
2463
bc41680b 2464//______________________________________________
85c4406e 2465void AliAnaPhoton::MakeAnalysisFillHistograms()
1c5acb87 2466{
6175da48 2467 //Fill histograms
85c4406e 2468
f27fe026 2469 //In case of simulated data, fill acceptance histograms
2470 if(IsDataMC()) FillAcceptanceHistograms();
2471
6175da48 2472 // Get vertex
2244659d 2473 Double_t v[3] = {0,0,0}; //vertex ;
2474 GetReader()->GetVertex(v);
85c4406e 2475 //fhVertex->Fill(v[0],v[1],v[2]);
6175da48 2476 if(TMath::Abs(v[2]) > GetZvertexCut()) return ; // done elsewhere for Single Event analysis, but there for mixed event
2477
2478 //----------------------------------
577d9801 2479 //Loop on stored AOD photons
2480 Int_t naod = GetOutputAODBranch()->GetEntriesFast();
577d9801 2481 if(GetDebug() > 0) printf("AliAnaPhoton::MakeAnalysisFillHistograms() - aod branch entries %d\n", naod);
521636d2 2482
c8710850 2483 Float_t cen = GetEventCentrality();
85c4406e 2484 // printf("++++++++++ GetEventCentrality() %f\n",cen);
2485
c8710850 2486 Float_t ep = GetEventPlaneAngle();
2487
3c1d9afb 2488 for(Int_t iaod = 0; iaod < naod ; iaod++)
2489 {
577d9801 2490 AliAODPWG4Particle* ph = (AliAODPWG4Particle*) (GetOutputAODBranch()->At(iaod));
2491 Int_t pdg = ph->GetIdentifiedParticleType();
521636d2 2492
85c4406e 2493 if(GetDebug() > 3)
2494 printf("AliAnaPhoton::MakeAnalysisFillHistograms() - PDG %d, MC TAG %d, Calorimeter %s\n",
3c1d9afb 2495 ph->GetIdentifiedParticleType(),ph->GetTag(), (ph->GetDetector()).Data()) ;
521636d2 2496
577d9801 2497 //If PID used, fill histos with photons in Calorimeter fCalorimeter
85c4406e 2498 if(IsCaloPIDOn() && pdg != AliCaloPID::kPhoton) continue;
bc41680b 2499
577d9801 2500 if(ph->GetDetector() != fCalorimeter) continue;
521636d2