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