AliPhysics  cdeda5a (cdeda5a)
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
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 #include "AliMCEvent.h"
36 
37 // --- Detectors ---
38 #include "AliPHOSGeoUtils.h"
39 #include "AliEMCALGeometry.h"
40 
44 
45 //____________________________
48 //____________________________
51 fMinDist(0.), fMinDist2(0.), fMinDist3(0.),
52 fRejectTrackMatch(0), fFillTMHisto(kFALSE),
53 fTimeCutMin(-10000), fTimeCutMax(10000),
54 fNCellsCut(0),
55 fNLMCutMin(-1), fNLMCutMax(10),
56 fFillSSHistograms(0), fFillEMCALRegionSSHistograms(0),
57 fFillConversionVertexHisto(0),fFillOnlySimpleSSHisto(1),
58 fFillSSNLocMaxHisto(0),
59 fNOriginHistograms(9), fNPrimaryHistograms(5),
60 fMomentum(), fMomentum2(),
61 fPrimaryMom(), fProdVertex(),
62 fConstantTimeShift(0), fFillEBinAcceptanceHisto(0), fNEBinCuts(0),
63 fStudyActivityNearCluster(0),
64 // Histograms
65 
66 // Control histograms
67 fhNCellsE(0), fhCellsE(0),
68 fhMaxCellDiffClusterE(0), fhTimePt(0), fhEtaPhi(0),
69 
70 fhEPhoton(0), fhPtPhoton(0),
71 fhPhiPhoton(0), fhEtaPhoton(0),
72 fhEtaPhiPhoton(0), fhEtaPhi05Photon(0),
73 fhPtCentralityPhoton(0), fhPtEventPlanePhoton(0),
74 
75 // Shower shape histograms
76 fhNLocMax(0),
77 fhDispE(0), fhDispPt(0),
78 fhLam0E(0), fhLam0Pt(0),
79 fhLam1E(0), fhLam1Pt(0),
80 fhLam0PtNLM1(0), fhLam0PtNLM2(0),
81 fhLam1PtNLM1(0), fhLam1PtNLM2(0),
82 fhDispETRD(0), fhLam0ETRD(0), fhLam0PtTRD(0), fhLam1ETRD(0),
83 fhDispETM(0), fhLam0ETM(0), fhLam0PtTM(0), fhLam1ETM(0),
84 fhDispETMTRD(0), fhLam0ETMTRD(0), fhLam0PtTMTRD(0), fhLam1ETMTRD(0),
85 
86 fhNCellsLam0LowE(0), fhNCellsLam1LowE(0), fhNCellsDispLowE(0),
87 fhNCellsLam0HighE(0), fhNCellsLam1HighE(0), fhNCellsDispHighE(0),
88 
89 fhEtaLam0LowE(0), fhPhiLam0LowE(0),
90 fhEtaLam0HighE(0), fhPhiLam0HighE(0),
91 fhLam0DispLowE(0), fhLam0DispHighE(0),
92 fhLam1Lam0LowE(0), fhLam1Lam0HighE(0),
93 fhDispLam1LowE(0), fhDispLam1HighE(0),
94 fhDispEtaE(0), fhDispPhiE(0),
95 fhSumEtaE(0), fhSumPhiE(0), fhSumEtaPhiE(0),
96 fhDispEtaPhiDiffE(0), fhSphericityE(0),
97 fhDispSumEtaDiffE(0), fhDispSumPhiDiffE(0),
98 
99 // MC histograms
100 fhMCPhotonELambda0NoOverlap(0), fhMCPhotonELambda0TwoOverlap(0), fhMCPhotonELambda0NOverlap(0),
101 // Embedding
102 fhEmbeddedSignalFractionEnergy(0),
103 fhEmbedPhotonELambda0FullSignal(0), fhEmbedPhotonELambda0MostlySignal(0),
104 fhEmbedPhotonELambda0MostlyBkg(0), fhEmbedPhotonELambda0FullBkg(0),
105 fhEmbedPi0ELambda0FullSignal(0), fhEmbedPi0ELambda0MostlySignal(0),
106 fhEmbedPi0ELambda0MostlyBkg(0), fhEmbedPi0ELambda0FullBkg(0),
107 
108 fhTimePtPhotonNoCut(0), fhTimePtPhotonSPD(0),
109 fhTimeNPileUpVertSPD(0), fhTimeNPileUpVertTrack(0),
110 fhPtPhotonNPileUpSPDVtx(0), fhPtPhotonNPileUpTrkVtx(0),
111 fhPtPhotonNPileUpSPDVtxTimeCut(0), fhPtPhotonNPileUpTrkVtxTimeCut(0),
112 fhPtPhotonNPileUpSPDVtxTimeCut2(0), fhPtPhotonNPileUpTrkVtxTimeCut2(0),
113 
114 fhEClusterSM(0), fhEPhotonSM(0),
115 fhPtClusterSM(0), fhPtPhotonSM(0),
116 fhMCConversionVertex(0), fhMCConversionVertexTRD(0)
117 {
118  for(Int_t i = 0; i < fgkNmcTypes; i++)
119  {
120  fhMCPt [i] = 0;
121  fhMCE [i] = 0;
122  fhMCPhi [i] = 0;
123  fhMCEta [i] = 0;
124  fhMCDeltaE [i] = 0;
125  fhMCDeltaPt[i] = 0;
126  fhMC2E [i] = 0;
127  fhMC2Pt [i] = 0;
128  }
129 
130  for(Int_t i = 0; i < fgkNmcPrimTypes; i++)
131  {
132  fhPtPrimMC [i] = 0;
133  fhEPrimMC [i] = 0;
134  fhPhiPrimMC[i] = 0;
135  fhEtaPrimMC[i] = 0;
136  fhYPrimMC [i] = 0;
137 
138  fhPtPrimMCAcc [i] = 0;
139  fhEPrimMCAcc [i] = 0;
140  fhPhiPrimMCAcc[i] = 0;
141  fhEtaPrimMCAcc[i] = 0;
142  fhYPrimMCAcc [i] = 0;
143  }
144 
145  for(Int_t i = 0; i < 7; i++)
146  {
147  fhDispEtaDispPhi[i] = 0;
148  fhLambda0DispPhi[i] = 0;
149  fhLambda0DispEta[i] = 0;
150 
151  fhPtPhotonPileUp[i] = 0;
153 
154  for(Int_t j = 0; j < fgkNssTypes; j++)
155  {
156  fhMCDispEtaDispPhi[i][j] = 0;
157  fhMCLambda0DispEta[i][j] = 0;
158  fhMCLambda0DispPhi[i][j] = 0;
159  }
160  }
161 
162  for(Int_t i = 0; i < fgkNssTypes; i++)
163  {
164  fhMCELambda0 [i] = 0;
165  fhMCPtLambda0 [i] = 0;
166  fhMCELambda1 [i] = 0;
167  fhMCEDispersion [i] = 0;
168  fhMCNCellsE [i] = 0;
170  fhLambda0DispEta[i] = 0;
171  fhLambda0DispPhi[i] = 0;
172 
179 
180  fhMCEDispEta [i] = 0;
181  fhMCEDispPhi [i] = 0;
182  fhMCESumEtaPhi [i] = 0;
183  fhMCEDispEtaPhiDiff[i] = 0;
184  fhMCESphericity [i] = 0;
185  }
186 
187  for(Int_t i = 0; i < 5; i++)
188  {
189  fhClusterCutsE [i] = 0;
190  fhClusterCutsPt[i] = 0;
191  }
192 
193  // Track matching residuals
194  for(Int_t i = 0; i < 2; i++)
195  {
204  fhdEdx[i] = 0; fhEOverP[i] = 0;
205  fhEOverPTRD[i] = 0;
206  }
207 
208  for(Int_t i = 0; i < 6; i++)
209  {
210  fhMCConversionLambda0Rcut [i] = 0;
212  fhMCConversionLambda1Rcut [i] = 0;
214  }
215 
216  for(Int_t ieta = 0; ieta < 4; ieta++)
217  {
218  for(Int_t iphi = 0; iphi < 3; iphi++)
219  {
220 // fhLam0EMCALRegion [ieta][iphi] = 0;
221 // fhLam0EMCALRegionTRD[ieta][iphi] = 0;
222 //
223 // for(Int_t i = 0; i < 6; i++)
224 // {
225 // fhLam0EMCALRegionMCConvRcut [ieta][iphi][i] = 0;
226 // fhLam0EMCALRegionTRDMCConvRcut[ieta][iphi][i] = 0;
227 // }
228  for(Int_t ism =0; ism < 20; ism++)
229  {
230  fhLam0EMCALRegionPerSM[ieta][iphi][ism] = 0;
231  fhLam1EMCALRegionPerSM[ieta][iphi][ism] = 0;
232  }
233  }
234  }
235 
236  for(Int_t il0 = 0; il0 < 2; il0++)
237  {
238  for(Int_t i = 0; i < 7; i++)
239  {
240  fhEtaPhiLam0BinPtBin [il0][i] = 0 ;
241  fhEtaPhiLargeTimeInClusterCell[il0][i] = 0 ;
242 
243  fhColRowLam0BinPtBin [il0][i] = 0 ;
244  fhColRowLam0BinPtBinWeighted [il0][i] = 0 ;
245  fhColRowLam0BinPtBinLargeTime [il0][i] = 0 ;
246  fhCellClusterIndexEAndTime [il0][i] = 0 ;
247  fhCellClusterEAndTime [il0][i] = 0 ;
248  fhCellClusterEFracAndTime [il0][i] = 0 ;
249  }
250 
251  for(Int_t ism =0; ism < 20; ism++)
252  {
253  fhLam1Lam0BinPerSM [il0][ism] = 0;
254  fhTimeLam0BinPerSM [il0][ism] = 0;
255  fhTimeLam0BinPerSMWeighted [il0][ism] = 0;
256  fhDTimeLam0BinPerSM [il0][ism] = 0;
257  fhDTimeLam0BinPerSMWeighted [il0][ism] = 0;
258  fhCellClusterEFracLam0BinPerSM [il0][ism] = 0;
259 // fhCellClusterEFracLam0BinPerSMWeighted[il0][ism] = 0;
260  fhCellClusterELam0BinPerSM [il0][ism] = 0;
262 
265  fhCellClusterELam0BinPerSMLargeTime [il0][ism] = 0;
267  fhNCellsWithLargeTimeInCluster [il0][ism] = 0;
268  }
269  }
270 
271 // for(Int_t il0 = 0; il0 < 2; il0++)
272 // {
273 // for(Int_t i = 0; i < 7; i++)
274 // {
275 // fhEtaPhiLam0BinPtBinSMShared [il0][i] = 0 ;
276 // fhColRowLam0BinPtBinSMShared [il0][i] = 0 ;
277 // }
278 //
279 // for(Int_t ism =0; ism < 12; ism++)
280 // {
281 // fhTimeLam0BinPerSMShared[il0][ism] = 0;
282 // fhLam0PerSMShared [ism] = 0;
283 // fhLam1PerSMShared [ism] = 0;
284 // }
285 // }
286 
287  for(Int_t ism = 0; ism < 20; ism++)
288  {
289  fhLam0PerSM [ism] = 0;
290  fhLam1PerSM [ism] = 0;
293 // fhLam0PerSMSPDPileUp [ism] = 0;
294 // fhLam1PerSMSPDPileUp [ism] = 0;
295  }
296 
297  for(Int_t ilarge = 0; ilarge < 5; ilarge++)
298  {
301  }
302 
303  for(Int_t i = 0; i < 14; i++)
304  {
305  fhEBinClusterEtaPhi[i] = 0 ;
306  fhEBinClusterColRow[i] = 0 ;
307  fhEBinClusterEtaPhiPID[i] = 0 ;
308  fhEBinClusterColRowPID[i] = 0 ;
309  }
310 
311  for(Int_t igen = 0; igen < 10; igen++)
312  {
313  for(Int_t ip = 0; ip < fgkNGenTypes; ip++)
314  {
315  fhMergeGeneratorCluster [igen][ip] = 0;
316  fhMergeGeneratorClusterNotHijingBkg [igen][ip] = 0;
318  fhMergeGeneratorClusterHijingBkg [igen][ip] = 0;
319  fhCleanGeneratorCluster [igen][ip] = 0;
320 
326 
327  fhMergeGeneratorClusterEPrimRecoDiff [igen][ip] = 0;
331  fhCleanGeneratorClusterEPrimRecoDiff [igen][ip] = 0;
332  }
333  }
334 
335  for(Int_t icase = 0; icase < 6; icase++)
336  {
337  fhLocalRegionClusterEtaPhi [icase] = 0;
338  fhLocalRegionClusterEnergySum [icase] = 0;
350  }
351 
352  // Initialize parameters
353  InitParameters();
354 }
355 
356 //_____________________________________________________________________
365 //__________________________________________________________________________
367  Float_t phi, TObjArray *clusterList)
368 {
369  Float_t radius = 0.2; // Hardcoded unless real use.
370 
371  // Accept cluster on EMCal limits - radius
372  // CAREFUL if used for Run2 with DCal.
373  if(phi < 3.15 - radius && phi > 1.4 + radius && TMath::Abs(eta) < 0.7-radius)
374  {
375  TString genName, genNameBkg;
376  Int_t genBkgTag = -1;
378  {
379  AliVCluster * calo = (AliVCluster*) (clusterList->At(icalo));
380  genBkgTag = GetCocktailGeneratorBackgroundTag(calo, genName, genNameBkg);
381  }
382 
383  fhLocalRegionClusterEtaPhi[0] ->Fill(eta,phi, GetEventWeight());
384 
385  Float_t sumE = 0;
386  Int_t sumM = 0;
387  Float_t sumEHi = 0;
388  Int_t sumMHi = 0;
389  for(Int_t icalo2 = 0; icalo2 < clusterList->GetEntriesFast(); icalo2++)
390  {
391  if ( icalo2 == icalo ) continue;
392 
393  AliVCluster * calo2 = (AliVCluster*) (clusterList->At(icalo2));
394 
395  // Select clusters in a radius of R=0.2
396  calo2->GetMomentum(fMomentum2,GetVertex(0)) ;
397 
398  Float_t dEta = eta-fMomentum2.Eta();
399  Float_t dPhi = phi-GetPhi(fMomentum2.Phi());
400 
401  if(TMath::Abs(dPhi) >= TMath::Pi())
402  dPhi = TMath::TwoPi()-TMath::Abs(dPhi);
403 
404  if (TMath::Sqrt( dEta*dEta + dPhi*dPhi ) > 0.2) continue;
405 
406  sumM++;
407  sumE += calo2->E();
408 
409  if( IsDataMC() && IsStudyClusterOverlapsPerGeneratorOn() && calo2->GetNLabels() > 0)
410  {
411  TString genName2;
412  (GetReader()->GetMC())->GetCocktailGenerator(calo2->GetLabel(),genName2);
413 
414  if(genName2.Contains("ijing"))
415  {
416  sumMHi++;
417  sumEHi += calo2->E();
418  }
419  }
420  }
421 
422  Float_t sumMAdd = sumM-sumMHi;
423  Float_t sumEAdd = sumE-sumEHi;
424 
425  fhLocalRegionClusterEnergySum [0]->Fill(en,sumE,GetEventWeight());
427 
429  {
432  }
433 
435  {
436  fhLocalRegionClusterEnergySumHijing [0]->Fill(en,sumEHi ,GetEventWeight());
438 
439  fhLocalRegionClusterEnergySumAdded [0]->Fill(en,sumEAdd,GetEventWeight());
441 
443  {
446 
449  }
450 
451  if(genBkgTag < 0) return;
452 
453  fhLocalRegionClusterEtaPhi [genBkgTag+1] ->Fill(eta,phi, GetEventWeight());
454 
455  fhLocalRegionClusterEnergySum [genBkgTag+1]->Fill(en,sumE,GetEventWeight());
456  fhLocalRegionClusterMultiplicity [genBkgTag+1]->Fill(en,sumM,GetEventWeight());
457 
458  fhLocalRegionClusterEnergySumHijing [genBkgTag+1]->Fill(en,sumEHi ,GetEventWeight());
459  fhLocalRegionClusterMultiplicityHijing[genBkgTag+1]->Fill(en,sumMHi ,GetEventWeight());
460 
461  fhLocalRegionClusterEnergySumAdded [genBkgTag+1]->Fill(en,sumEAdd,GetEventWeight());
462  fhLocalRegionClusterMultiplicityAdded [genBkgTag+1]->Fill(en,sumMAdd,GetEventWeight());
463 
465  {
468 
471 
474  }
475 
476  if(genBkgTag > 0) // Not clean, all merged
477  {
478  fhLocalRegionClusterEtaPhi [5] ->Fill(eta,phi, GetEventWeight());
479 
480  fhLocalRegionClusterEnergySum [5]->Fill(en,sumE,GetEventWeight());
482 
483  fhLocalRegionClusterEnergySumHijing [5]->Fill(en,sumEHi ,GetEventWeight());
485 
486  fhLocalRegionClusterEnergySumAdded [5]->Fill(en,sumEAdd,GetEventWeight());
488 
490  {
493 
496 
499  }
500 
501  } // not clean, all merged
502 
503  } // check generator overlaps
504 
505  } // EMCal acceptance
506 }
507 
508 //_____________________________________________________________________
514 //__________________________________________________________________________
516 {
517  //
518  // Check the generators inside the cluster
519  TString genName = "", genNameBkg = "";
520  Int_t genBkgTag = GetCocktailGeneratorBackgroundTag(calo, genName, genNameBkg);
521  if (genBkgTag == -1) return;
522  else if(genBkgTag > 3) printf("Bkg generator tag larger than 3\n");
523 
524  //
525  // Get primary particle info of main particle contributing to the cluster
526  Float_t eprim = 0;
527  //Float_t ptprim = 0;
528  Bool_t ok = kFALSE;
529  Int_t pdg = 0, status = 0, momLabel = -1;
530 
531  //fPrimaryMom = GetMCAnalysisUtils()->GetMother(label,GetReader(),ok);
532  fPrimaryMom = GetMCAnalysisUtils()->GetMother(calo->GetLabel(),GetReader(), pdg, status, ok, momLabel);
533 
534  if(ok)
535  {
536  eprim = fPrimaryMom.Energy();
537  //ptprim = fPrimaryMom.Pt();
538  }
539 
540  //
541  Float_t en = calo->E();
542 
543  if ( eprim < 0.1 || en < 0.5 ) return;
544 
545  //
546  Int_t partType = -1;
547  if ( GetMCAnalysisUtils()->CheckTagBit(mctag,AliMCAnalysisUtils::kMCPi0) ) partType = kmcGenPi0Merged;
548  else if( GetMCAnalysisUtils()->CheckTagBit(mctag,AliMCAnalysisUtils::kMCPi0Decay) ) partType = kmcGenPi0Decay;
549  else if( GetMCAnalysisUtils()->CheckTagBit(mctag,AliMCAnalysisUtils::kMCEtaDecay) ) partType = kmcGenEtaDecay;
550  else if( GetMCAnalysisUtils()->CheckTagBit(mctag,AliMCAnalysisUtils::kMCPhoton) ) partType = kmcGenPhoton;
551  else if( GetMCAnalysisUtils()->CheckTagBit(mctag,AliMCAnalysisUtils::kMCElectron) ) partType = kmcGenElectron;
552  else partType = kmcGenOther;
553 
554  Int_t genType = GetNCocktailGenNamesToCheck()-1;
555  for(Int_t igen = 1; igen < GetNCocktailGenNamesToCheck(); igen++)
556  {
557  if ( genName.Contains(GetCocktailGenNameToCheck(igen)) )
558  {
559  genType = igen;
560  break;
561  }
562  }
563 
564  Float_t ratio = en / eprim;
565  Float_t diff = en - eprim;
566 
567  if ( genBkgTag > 0 )
568  {
569  fhMergeGeneratorCluster[0] [0] ->Fill(en,GetEventWeight());
570  fhMergeGeneratorCluster[0] [partType] ->Fill(en,GetEventWeight());
571  fhMergeGeneratorCluster[genType][0] ->Fill(en,GetEventWeight());
572  fhMergeGeneratorCluster[genType][partType] ->Fill(en,GetEventWeight());
573 
574  fhMergeGeneratorClusterEPrimRecoRatio[0] [0] ->Fill(en,ratio,GetEventWeight());
575  fhMergeGeneratorClusterEPrimRecoRatio[0] [partType] ->Fill(en,ratio,GetEventWeight());
576  fhMergeGeneratorClusterEPrimRecoRatio[genType][0] ->Fill(en,ratio,GetEventWeight());
577  fhMergeGeneratorClusterEPrimRecoRatio[genType][partType] ->Fill(en,ratio,GetEventWeight());
578 
579  fhMergeGeneratorClusterEPrimRecoDiff[0] [0] ->Fill(en,diff,GetEventWeight());
580  fhMergeGeneratorClusterEPrimRecoDiff[0] [partType] ->Fill(en,diff,GetEventWeight());
581  fhMergeGeneratorClusterEPrimRecoDiff[genType][0] ->Fill(en,diff,GetEventWeight());
582  fhMergeGeneratorClusterEPrimRecoDiff[genType][partType] ->Fill(en,diff,GetEventWeight());
583 
584  if(!genName.Contains("ijing"))
585  {
586  if(genBkgTag == 2)
587  {
589  fhMergeGeneratorClusterNotHijingBkg[0] [partType] ->Fill(en,GetEventWeight());
590  fhMergeGeneratorClusterNotHijingBkg[genType][0] ->Fill(en,GetEventWeight());
591  fhMergeGeneratorClusterNotHijingBkg[genType][partType] ->Fill(en,GetEventWeight());
592 
594  fhMergeGeneratorClusterNotHijingBkgEPrimRecoRatio[0] [partType] ->Fill(en,ratio,GetEventWeight());
596  fhMergeGeneratorClusterNotHijingBkgEPrimRecoRatio[genType][partType] ->Fill(en,ratio,GetEventWeight());
597 
601  fhMergeGeneratorClusterNotHijingBkgEPrimRecoDiff[genType][partType] ->Fill(en,diff,GetEventWeight());
602  }
603 
604  if(genBkgTag == 1)
605  {
607  fhMergeGeneratorClusterHijingBkg[0] [partType] ->Fill(en,GetEventWeight());
608  fhMergeGeneratorClusterHijingBkg[genType][0] ->Fill(en,GetEventWeight());
609  fhMergeGeneratorClusterHijingBkg[genType][partType] ->Fill(en,GetEventWeight());
610 
612  fhMergeGeneratorClusterHijingBkgEPrimRecoRatio[0] [partType] ->Fill(en,ratio,GetEventWeight());
613  fhMergeGeneratorClusterHijingBkgEPrimRecoRatio[genType][0] ->Fill(en,ratio,GetEventWeight());
614  fhMergeGeneratorClusterHijingBkgEPrimRecoRatio[genType][partType] ->Fill(en,ratio,GetEventWeight());
615 
617  fhMergeGeneratorClusterHijingBkgEPrimRecoDiff[0] [partType] ->Fill(en,diff,GetEventWeight());
619  fhMergeGeneratorClusterHijingBkgEPrimRecoDiff[genType][partType] ->Fill(en,diff,GetEventWeight());
620  }
621 
622  if ( genBkgTag == 3 )
623  {
625  fhMergeGeneratorClusterHijingAndOtherBkg[0] [partType] ->Fill(en,GetEventWeight());
627  fhMergeGeneratorClusterHijingAndOtherBkg[genType][partType] ->Fill(en,GetEventWeight());
628 
632  fhMergeGeneratorClusterHijingAndOtherBkgEPrimRecoRatio[genType][partType] ->Fill(en,ratio,GetEventWeight());
633 
637  fhMergeGeneratorClusterHijingAndOtherBkgEPrimRecoDiff[genType][partType] ->Fill(en,diff,GetEventWeight());
638  }
639  }
640  }
641  else
642  {
643  fhCleanGeneratorCluster[0] [0] ->Fill(en,GetEventWeight());
644  fhCleanGeneratorCluster[0] [partType] ->Fill(en,GetEventWeight());
645  fhCleanGeneratorCluster[genType][0] ->Fill(en,GetEventWeight());
646  fhCleanGeneratorCluster[genType][partType] ->Fill(en,GetEventWeight());
647 
648  fhCleanGeneratorClusterEPrimRecoRatio[0] [0] ->Fill(en,ratio,GetEventWeight());
649  fhCleanGeneratorClusterEPrimRecoRatio[0] [partType] ->Fill(en,ratio,GetEventWeight());
650  fhCleanGeneratorClusterEPrimRecoRatio[genType][0] ->Fill(en,ratio,GetEventWeight());
651  fhCleanGeneratorClusterEPrimRecoRatio[genType][partType] ->Fill(en,ratio,GetEventWeight());
652 
653  fhCleanGeneratorClusterEPrimRecoDiff[0] [0] ->Fill(en,diff,GetEventWeight());
654  fhCleanGeneratorClusterEPrimRecoDiff[0] [partType] ->Fill(en,diff,GetEventWeight());
655  fhCleanGeneratorClusterEPrimRecoDiff[genType][0] ->Fill(en,diff,GetEventWeight());
656  fhCleanGeneratorClusterEPrimRecoDiff[genType][partType] ->Fill(en,diff,GetEventWeight());
657  }
658 }
659 
660 
661 //_____________________________________________________________________
677 //_____________________________________________________________________
678 Bool_t AliAnaPhoton::ClusterSelected(AliVCluster* calo, Int_t nMaxima)
679 {
680  Float_t ptcluster = fMomentum.Pt();
681  Float_t ecluster = fMomentum.E();
682  Float_t etacluster = fMomentum.Eta();
683  Float_t phicluster = fMomentum.Phi();
684 
685  if(phicluster < 0) phicluster+=TMath::TwoPi();
686 
687  Bool_t matched = IsTrackMatched(calo,GetReader()->GetInputEvent());
688 
689  AliDebug(2,Form("Current Event %d; Before selection : E %2.2f, pT %2.2f, phi %2.2f, eta %2.2f",
691  ecluster,ptcluster, phicluster*TMath::RadToDeg(),etacluster));
692 
693  fhClusterCutsE [1]->Fill( ecluster, GetEventWeight());
694  fhClusterCutsPt[1]->Fill(ptcluster, GetEventWeight());
695 
696  if(ecluster > 0.5) fhEtaPhi->Fill(etacluster, phicluster, GetEventWeight());
697 
698  Int_t nSM = GetModuleNumber(calo);
699  if(nSM < GetCaloUtils()->GetNumberOfSuperModulesUsed() && nSM >=0)
700  {
701  fhEClusterSM ->Fill(ecluster , nSM, GetEventWeight());
702  fhPtClusterSM->Fill(ptcluster, nSM, GetEventWeight());
703  }
704 
705  //.......................................
706  //If too small or big energy, skip it
707  if(ecluster < GetMinEnergy() || ecluster > GetMaxEnergy() ) return kFALSE ;
708 
709  AliDebug(2,Form("\t Cluster %d Pass E Cut",calo->GetID()));
710 
711  fhClusterCutsE [2]->Fill( ecluster, GetEventWeight());
712  fhClusterCutsPt[2]->Fill(ptcluster, GetEventWeight());
713 
714  //.......................................
715  // TOF cut, BE CAREFUL WITH THIS CUT
716  Double_t tof = calo->GetTOF()*1e9;
717  if(tof > 400) tof-=fConstantTimeShift; // only for MC, rest shift = 0
718 
719  if(tof < fTimeCutMin || tof > fTimeCutMax) return kFALSE;
720 
721  AliDebug(2,Form("\t Cluster %d Pass Time Cut",calo->GetID()));
722 
723  fhClusterCutsE [3]->Fill( ecluster, GetEventWeight());
724  fhClusterCutsPt[3]->Fill(ptcluster, GetEventWeight());
725 
726  //.......................................
727  if(calo->GetNCells() <= fNCellsCut && GetReader()->GetDataType() != AliCaloTrackReader::kMC) return kFALSE;
728 
729  AliDebug(2,Form("\t Cluster %d Pass NCell Cut",calo->GetID()));
730 
731  fhClusterCutsE [4]->Fill( ecluster, GetEventWeight());
732  fhClusterCutsPt[4]->Fill(ptcluster, GetEventWeight());
733 
734  if(nMaxima < fNLMCutMin || nMaxima > fNLMCutMax) return kFALSE ;
735  AliDebug(2,Form("\t Cluster %d pass NLM %d of out of range",calo->GetID(), nMaxima));
736 
737  fhClusterCutsE [5]->Fill( ecluster, GetEventWeight());
738  fhClusterCutsPt[5]->Fill(ptcluster, GetEventWeight());
739 
740  //.......................................
741  //Check acceptance selection
742  if(IsFiducialCutOn())
743  {
745  if(! in ) return kFALSE ;
746  }
747 
748  AliDebug(2,Form("\t Fiducial cut passed"));
749 
750  fhClusterCutsE [6]->Fill( ecluster, GetEventWeight());
751  fhClusterCutsPt[6]->Fill(ptcluster, GetEventWeight());
752 
753  //.......................................
754  //Skip matched clusters with tracks
755 
756  // Fill matching residual histograms before PID cuts
758 
760  {
761  if(matched)
762  {
763  AliDebug(2,"\t Reject track-matched clusters");
764  return kFALSE ;
765  }
766  else
767  AliDebug(2,"\t Track-matching cut passed");
768  }// reject matched clusters
769 
770  fhClusterCutsE [7]->Fill( ecluster, GetEventWeight());
771  fhClusterCutsPt[7]->Fill(ptcluster, GetEventWeight());
772 
773  //.......................................
774  //Check Distance to Bad channel, set bit.
775  Double_t distBad=calo->GetDistanceToBadChannel() ; //Distance to bad channel
776  if(distBad < 0.) distBad=9999. ; //workout strange convension dist = -1. ;
777  if(distBad < fMinDist)
778  {//In bad channel (PHOS cristal size 2.2x2.2 cm), EMCAL ( cell units )
779  return kFALSE ;
780  }
781  else AliDebug(2,Form("\t Bad channel cut passed %4.2f > %2.2f",distBad, fMinDist));
782 
783  fhClusterCutsE [8]->Fill( ecluster, GetEventWeight());
784  fhClusterCutsPt[8]->Fill(ptcluster, GetEventWeight());
785 
786  AliDebug(1,Form("Current Event %d; After selection : E %2.2f, pT %2.2f, phi %2.2f, eta %2.2f",
788  ecluster, ptcluster,fMomentum.Phi()*TMath::RadToDeg(),fMomentum.Eta()));
789 
790  //All checks passed, cluster selected
791  return kTRUE;
792 
793 }
794 
795 //___________________________________________
800 //___________________________________________
802 {
803  Double_t photonY = -100 ;
804  Double_t photonE = -1 ;
805  Double_t photonPt = -1 ;
806  Double_t photonPhi = 100 ;
807  Double_t photonEta = -1 ;
808 
809  Int_t pdg = 0 ;
810  Int_t tag = 0 ;
811  Int_t status = 0 ;
812  Int_t mcIndex = 0 ;
813  Int_t nprim = 0 ;
814  Bool_t inacceptance = kFALSE ;
815 
816  TParticle * primStack = 0;
817  AliAODMCParticle * primAOD = 0;
818 
819  // Get the ESD MC particles container
820  AliStack * stack = 0;
821  if( GetReader()->ReadStack() )
822  {
823  stack = GetMCStack();
824  if( !stack )
825  {
826  AliFatal("Stack not available, is the MC handler called? STOP");
827  return;
828  }
829  nprim = stack->GetNtrack();
830  }
831 
832  // Get the AOD MC particles container
833  TClonesArray * mcparticles = 0;
834  if( GetReader()->ReadAODMCParticles() )
835  {
836  mcparticles = GetReader()->GetAODMCParticles();
837  if( !mcparticles )
838  {
839  AliFatal("Standard MCParticles not available!");
840  return;
841  }
842  nprim = mcparticles->GetEntriesFast();
843  }
844 
845  for(Int_t i=0 ; i < nprim; i++)
846  {
847  if ( !GetReader()->AcceptParticleMCLabel( i ) ) continue ;
848 
849  if(GetReader()->ReadStack())
850  {
851  primStack = stack->Particle(i) ;
852  if(!primStack)
853  {
854  AliWarning("ESD primaries pointer not available!!");
855  continue;
856  }
857 
858  pdg = primStack->GetPdgCode();
859  status = primStack->GetStatusCode();
860 
861  // Protection against floating point exception
862  if ( primStack->Energy() == TMath::Abs(primStack->Pz()) ||
863  (primStack->Energy() - primStack->Pz()) < 1e-3 ||
864  (primStack->Energy() + primStack->Pz()) < 0 ) continue ;
865 
866  //printf("i %d, %s %d %s %d \n",i, stack->Particle(i)->GetName(), stack->Particle(i)->GetPdgCode(),
867  // prim->GetName(), prim->GetPdgCode());
868 
869  // Photon kinematics
870  primStack->Momentum(fMomentum);
871 
872  photonY = 0.5*TMath::Log((primStack->Energy()+primStack->Pz())/(primStack->Energy()-primStack->Pz())) ;
873  }
874  else
875  {
876  primAOD = (AliAODMCParticle *) mcparticles->At(i);
877  if(!primAOD)
878  {
879  AliWarning("AOD primaries pointer not available!!");
880  continue;
881  }
882 
883  pdg = primAOD->GetPdgCode();
884  status = primAOD->GetStatus();
885 
886  // Protection against floating point exception
887  if ( primAOD->E() == TMath::Abs(primAOD->Pz()) ||
888  (primAOD->E() - primAOD->Pz()) < 1e-3 ||
889  (primAOD->E() + primAOD->Pz()) < 0 ) continue ;
890 
891  // Photon kinematics
892  fMomentum.SetPxPyPzE(primAOD->Px(),primAOD->Py(),primAOD->Pz(),primAOD->E());
893 
894  photonY = 0.5*TMath::Log((primAOD->E()+primAOD->Pz())/(primAOD->E()-primAOD->Pz())) ;
895  }
896 
897  // Select only photons in the final state
898  if(pdg != 22 ) continue ;
899 
900  // If too small or too large pt, skip, same cut as for data analysis
901  photonPt = fMomentum.Pt () ;
902 
903  if(photonPt < GetMinPt() || photonPt > GetMaxPt() ) continue ;
904 
905  photonE = fMomentum.E () ;
906  photonEta = fMomentum.Eta() ;
907  photonPhi = fMomentum.Phi() ;
908 
909  if(photonPhi < 0) photonPhi+=TMath::TwoPi();
910 
911  // Check if photons hit desired acceptance
912  inacceptance = kTRUE;
913 
914  // Check same fidutial borders as in data analysis on top of real acceptance if real was requested.
915  if( IsFiducialCutOn() && !GetFiducialCut()->IsInFiducialCut(fMomentum.Eta(),fMomentum.Phi(),GetCalorimeter())) inacceptance = kFALSE ;
916 
917  // Check if photons hit the Calorimeter acceptance
918  if(IsRealCaloAcceptanceOn()) // defined on base class
919  {
920  if(GetReader()->ReadStack() &&
921  !GetCaloUtils()->IsMCParticleInCalorimeterAcceptance(GetCalorimeter(), primStack)) inacceptance = kFALSE ;
922  if(GetReader()->ReadAODMCParticles() &&
923  !GetCaloUtils()->IsMCParticleInCalorimeterAcceptance(GetCalorimeter(), primAOD )) inacceptance = kFALSE ;
924  }
925 
926  // Get tag of this particle photon from fragmentation, decay, prompt ...
927  // Set the origin of the photon.
929 
930  if(!GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPhoton))
931  {
932  // A conversion photon from a hadron, skip this kind of photon
933  // printf("AliAnaPhoton::FillAcceptanceHistograms() - not a photon, weird!\n ");
934  // GetMCAnalysisUtils()->PrintMCTag(tag);
935 
936  continue;
937  }
938 
939  // Consider only final state particles, but this depends on generator,
940  // status 1 is the usual one, in case of not being ok, leave the possibility
941  // to not consider this.
942  if(status > 1) continue ; // Avoid "partonic" photons
943 
944  Bool_t takeIt = kFALSE ;
945  if(status == 1 && GetMCAnalysisUtils()->GetMCGenerator() != AliMCAnalysisUtils::kBoxLike ) takeIt = kTRUE ;
946 
947  if (GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCConversion)) continue;
948 
949  // Origin of photon
950  if (GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPrompt))
951  {
952  mcIndex = kmcPPrompt;
953  }
954  else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCFragmentation))
955  {
956  mcIndex = kmcPFragmentation ;
957  }
958  else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCISR))
959  {
960  mcIndex = kmcPISR;
961  }
962  else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0Decay))
963  {
964  mcIndex = kmcPPi0Decay;
965  }
966  else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCEtaDecay))
967  {
968  mcIndex = kmcPEtaDecay;
969  }
970  else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCOtherDecay))
971  {
972  mcIndex = kmcPOtherDecay;
973  }
974  else
975  {
976  // Other decay but from non final state particle
977  mcIndex = kmcPOtherDecay;
978  } // Other origin
979 
980  if(!takeIt && (mcIndex == kmcPPi0Decay || mcIndex == kmcPOtherDecay)) takeIt = kTRUE ;
981 
982  if(!takeIt) continue ;
983 
984  // Fill histograms for all photons
985  fhYPrimMC[kmcPPhoton]->Fill(photonPt, photonY, GetEventWeight()) ;
986  if(TMath::Abs(photonY) < 1.0)
987  {
988  fhEPrimMC [kmcPPhoton]->Fill(photonE , GetEventWeight()) ;
989  fhPtPrimMC [kmcPPhoton]->Fill(photonPt, GetEventWeight()) ;
990  fhPhiPrimMC[kmcPPhoton]->Fill(photonE , photonPhi, GetEventWeight()) ;
991  fhEtaPrimMC[kmcPPhoton]->Fill(photonE , photonEta, GetEventWeight()) ;
992  }
993 
994  if(inacceptance)
995  {
996  fhEPrimMCAcc [kmcPPhoton]->Fill(photonE , GetEventWeight()) ;
997  fhPtPrimMCAcc [kmcPPhoton]->Fill(photonPt, GetEventWeight()) ;
998  fhPhiPrimMCAcc[kmcPPhoton]->Fill(photonE , photonPhi, GetEventWeight()) ;
999  fhEtaPrimMCAcc[kmcPPhoton]->Fill(photonE , photonEta, GetEventWeight()) ;
1000  fhYPrimMCAcc [kmcPPhoton]->Fill(photonE , photonY , GetEventWeight()) ;
1001  } // Accepted
1002 
1003  // Fill histograms for photons origin
1004  if(mcIndex < fNPrimaryHistograms)
1005  {
1006  fhYPrimMC[mcIndex]->Fill(photonPt, photonY, GetEventWeight()) ;
1007  if(TMath::Abs(photonY) < 1.0)
1008  {
1009  fhEPrimMC [mcIndex]->Fill(photonE , GetEventWeight()) ;
1010  fhPtPrimMC [mcIndex]->Fill(photonPt, GetEventWeight()) ;
1011  fhPhiPrimMC[mcIndex]->Fill(photonE , photonPhi, GetEventWeight()) ;
1012  fhEtaPrimMC[mcIndex]->Fill(photonE , photonEta, GetEventWeight()) ;
1013  }
1014 
1015  if(inacceptance)
1016  {
1017  fhEPrimMCAcc [mcIndex]->Fill(photonE , GetEventWeight()) ;
1018  fhPtPrimMCAcc [mcIndex]->Fill(photonPt, GetEventWeight()) ;
1019  fhPhiPrimMCAcc[mcIndex]->Fill(photonE , photonPhi, GetEventWeight()) ;
1020  fhEtaPrimMCAcc[mcIndex]->Fill(photonE , photonEta, GetEventWeight()) ;
1021  fhYPrimMCAcc [mcIndex]->Fill(photonE , photonY , GetEventWeight()) ;
1022  } // Accepted
1023  }
1024  } // loop on primaries
1025 }
1026 
1027 //________________________________________________________________________________
1029 //________________________________________________________________________________
1030 void AliAnaPhoton::FillPileUpHistograms(AliVCluster* cluster, AliVCaloCells *cells,
1031  Int_t absIdMax)
1032 {
1033  Float_t pt = fMomentum.Pt();
1034  Float_t time = cluster->GetTOF()*1.e9;
1035  if ( time > 400 ) time-=fConstantTimeShift; // Only for MC, rest shift = 0
1036 
1037  AliVEvent * event = GetReader()->GetInputEvent();
1038 
1039  if(GetReader()->IsPileUpFromSPD()) fhPtPhotonPileUp[0]->Fill(pt, GetEventWeight()) ;
1040  if(GetReader()->IsPileUpFromEMCal()) fhPtPhotonPileUp[1]->Fill(pt, GetEventWeight()) ;
1041  if(GetReader()->IsPileUpFromSPDOrEMCal()) fhPtPhotonPileUp[2]->Fill(pt, GetEventWeight()) ;
1042  if(GetReader()->IsPileUpFromSPDAndEMCal()) fhPtPhotonPileUp[3]->Fill(pt, GetEventWeight()) ;
1043  if(GetReader()->IsPileUpFromSPDAndNotEMCal()) fhPtPhotonPileUp[4]->Fill(pt, GetEventWeight()) ;
1044  if(GetReader()->IsPileUpFromEMCalAndNotSPD()) fhPtPhotonPileUp[5]->Fill(pt, GetEventWeight()) ;
1045  if(GetReader()->IsPileUpFromNotSPDAndNotEMCal()) fhPtPhotonPileUp[6]->Fill(pt, GetEventWeight()) ;
1046 
1047  fhTimePtPhotonNoCut->Fill(pt, time, GetEventWeight());
1048  if(GetReader()->IsPileUpFromSPD()) fhTimePtPhotonSPD->Fill(pt, time, GetEventWeight());
1049 
1050  //
1051  // Cells inside the cluster
1052  //
1053 
1054  // Loop on cells inside cluster, max cell must be over 100 MeV and time in BC=0
1055  if(cells->GetCellAmplitude(absIdMax) > 0.1 && TMath::Abs(time) < 30)
1056  {
1057  for (Int_t ipos = 0; ipos < cluster->GetNCells(); ipos++)
1058  {
1059  Int_t absId = cluster->GetCellsAbsId()[ipos];
1060 
1061  if( absId == absIdMax ) continue ;
1062 
1063  Double_t tcell = cells->GetCellTime(absId);
1064  Float_t amp = cells->GetCellAmplitude(absId);
1065  Int_t bc = GetReader()->GetInputEvent()->GetBunchCrossNumber();
1066 
1067  GetCaloUtils()->GetEMCALRecoUtils()->AcceptCalibrateCell(absId,bc,amp,tcell,cells);
1068  tcell*=1e9;
1069  tcell-=fConstantTimeShift;
1070 
1071  Float_t diff = (time-tcell);
1072 
1073  if( cells->GetCellAmplitude(absIdMax) < 0.1 ) continue ;
1074 
1082 
1083  } // loop
1084  }
1085 
1086  AliESDEvent* esdEv = dynamic_cast<AliESDEvent*> (event);
1087  AliAODEvent* aodEv = dynamic_cast<AliAODEvent*> (event);
1088 
1089  //
1090  // N pile up vertices
1091  //
1092 
1093  Int_t nVtxSPD = -1;
1094  Int_t nVtxTrk = -1;
1095 
1096  if (esdEv)
1097  {
1098  nVtxSPD = esdEv->GetNumberOfPileupVerticesSPD();
1099  nVtxTrk = esdEv->GetNumberOfPileupVerticesTracks();
1100 
1101  } // ESD
1102  else if (aodEv)
1103  {
1104  nVtxSPD = aodEv->GetNumberOfPileupVerticesSPD();
1105  nVtxTrk = aodEv->GetNumberOfPileupVerticesTracks();
1106  } // AOD
1107 
1108  if(pt < 8)
1109  {
1110  fhTimeNPileUpVertSPD ->Fill(time, nVtxSPD, GetEventWeight());
1111  fhTimeNPileUpVertTrack->Fill(time, nVtxTrk, GetEventWeight());
1112  }
1113 
1114  fhPtPhotonNPileUpSPDVtx->Fill(pt, nVtxSPD, GetEventWeight());
1115  fhPtPhotonNPileUpTrkVtx->Fill(pt, nVtxTrk, GetEventWeight());
1116 
1117  if(TMath::Abs(time) < 25)
1118  {
1119  fhPtPhotonNPileUpSPDVtxTimeCut->Fill(pt, nVtxSPD, GetEventWeight());
1120  fhPtPhotonNPileUpTrkVtxTimeCut->Fill(pt, nVtxTrk, GetEventWeight());
1121  }
1122 
1123  if(time < 75 && time > -25)
1124  {
1125  fhPtPhotonNPileUpSPDVtxTimeCut2->Fill(pt, nVtxSPD, GetEventWeight());
1126  fhPtPhotonNPileUpTrkVtxTimeCut2->Fill(pt, nVtxTrk, GetEventWeight());
1127  }
1128 }
1129 
1130 //_________________________________________________________________________________
1132 //_________________________________________________________________________________
1133 void AliAnaPhoton::FillShowerShapeHistograms(AliVCluster* cluster, Int_t mcTag, Int_t nlm,
1134  Float_t maxCellFraction, Int_t & largeTime)
1135 {
1136  if(!fFillSSHistograms || GetMixedEvent()) return;
1137 
1138  Float_t energy = cluster->E();
1139  Int_t ncells = cluster->GetNCells();
1140  Float_t lambda0 = cluster->GetM02();
1141  Float_t lambda1 = cluster->GetM20();
1142  Float_t disp = cluster->GetDispersion()*cluster->GetDispersion();
1143 
1144  Float_t pt = fMomentum.Pt();
1145  Float_t eta = fMomentum.Eta();
1146  Float_t phi = fMomentum.Phi();
1147  if(phi < 0) phi+=TMath::TwoPi();
1148 
1149  fhLam0E ->Fill(energy, lambda0, GetEventWeight());
1150  fhLam0Pt->Fill(pt , lambda0, GetEventWeight());
1151  fhLam1E ->Fill(energy, lambda1, GetEventWeight());
1152  fhLam1Pt->Fill(pt , lambda1, GetEventWeight());
1153  fhDispE ->Fill(energy, disp , GetEventWeight());
1154  fhDispPt->Fill(pt , disp , GetEventWeight());
1155 
1157  {
1158  if(nlm==1)
1159  {
1160  fhLam0PtNLM1->Fill(pt, lambda0, GetEventWeight());
1161  fhLam1PtNLM1->Fill(pt, lambda1, GetEventWeight());
1162  }
1163  else if(nlm==2)
1164  {
1165  fhLam0PtNLM2->Fill(pt, lambda0, GetEventWeight());
1166  fhLam1PtNLM2->Fill(pt, lambda1, GetEventWeight());
1167  }
1168  }
1169 
1170  if(GetCalorimeter() == kEMCAL && GetFirstSMCoveredByTRD() >= 0 &&
1171  GetModuleNumber(cluster) >= GetFirstSMCoveredByTRD() )
1172  {
1173  fhLam0ETRD ->Fill(energy, lambda0, GetEventWeight());
1174  fhLam0PtTRD->Fill(pt , lambda0, GetEventWeight());
1175  fhLam1ETRD ->Fill(energy, lambda1, GetEventWeight());
1176  fhDispETRD ->Fill(energy, disp, GetEventWeight());
1177  }
1178 
1179  //
1180  // EMCal SM regions
1181  //
1182  if(cluster->IsEMCAL() && fFillEMCALRegionSSHistograms)
1183  {
1184  Int_t sm = GetModuleNumber(cluster);
1185 
1186 // Bool_t shared = GetCaloUtils()->IsClusterSharedByTwoSuperModules(GetEMCALGeometry(),cluster);
1187 
1188  fhLam0PerSM[sm]->Fill(pt, lambda0, GetEventWeight());
1189  fhLam1PerSM[sm]->Fill(pt, lambda1, GetEventWeight());
1190 
1191 // if(GetReader()->IsPileUpFromSPD())
1192 // {
1193 // fhLam0PerSMSPDPileUp[sm]->Fill(pt, lambda0, GetEventWeight());
1194 // fhLam1PerSMSPDPileUp[sm]->Fill(pt, lambda1, GetEventWeight());
1195 // }
1196 
1197  Int_t etaRegion = -1, phiRegion = -1;
1198  GetCaloUtils()->GetEMCALSubregion(cluster,GetReader()->GetEMCALCells(),etaRegion,phiRegion);
1199  if(etaRegion >= 0 && etaRegion < 4 && phiRegion >=0 && phiRegion < 3)
1200  {
1201 
1202  // fhLam0EMCALRegion[etaRegion][phiRegion]->Fill(pt,lambda0, GetEventWeight());
1203  //
1204  // if(GetFirstSMCoveredByTRD() >= 0 && sm >= GetFirstSMCoveredByTRD() )
1205  // {
1206  // fhLam0EMCALRegionTRD[etaRegion][phiRegion]->Fill(pt, lambda0, GetEventWeight());
1207  // }
1208 
1209  fhLam0EMCALRegionPerSM[etaRegion][phiRegion][sm]->Fill(pt,lambda0, GetEventWeight());
1210  fhLam1EMCALRegionPerSM[etaRegion][phiRegion][sm]->Fill(pt,lambda1, GetEventWeight());
1211 
1212 // if(shared)
1213 // {
1214 // fhLam0PerSMShared[sm]->Fill(pt,lambda0, GetEventWeight());
1215 // fhLam1PerSMShared[sm]->Fill(pt,lambda1, GetEventWeight());
1216 // }
1217  }
1218 
1219  // Define l0 bin and cluster pt bin
1220  Int_t l0bin = -1;
1221  if ( lambda0 >=0.30 && lambda0 <= 0.40 ) l0bin = 1;
1222  else if( lambda0 >=0.23 && lambda0 <= 0.26 ) l0bin = 0;
1223 
1224  Float_t ptLimit[] = {2,3,4,5,6,8,10,12};
1225  Int_t ptbin = -1;
1226  for(Int_t ipt = 0; ipt < 7; ipt++)
1227  {
1228  if( pt >= ptLimit[ipt] && pt < ptLimit[ipt+1] )
1229  {
1230  ptbin = ipt;
1231  break;
1232  }
1233  }
1234 
1235  if(l0bin!=-1)
1236  {
1237  fhLam1Lam0BinPerSM[l0bin][sm]->Fill(pt,lambda1);
1238 
1239  if(ptbin >= 0)
1240  {
1241  fhEtaPhiLam0BinPtBin[l0bin][ptbin]->Fill(eta, phi, GetEventWeight());
1242 // if(shared) fhEtaPhiLam0BinPtBinSMShared[l0bin][ptbin]->Fill(eta, phi, GetEventWeight());
1243  }
1244  }
1245 
1246  //
1247  // Histograms with time of secondary cells, only if they contribute with a non null
1248  // weight to the shower shape
1249  //
1250  AliVCaloCells* cells = GetReader()->GetEMCALCells();
1251 
1252  Int_t bc = GetReader()->GetInputEvent()->GetBunchCrossNumber();
1253 
1254  //
1255  // Sort the cluster cells in energy, find maximum energy cell
1256  //
1257  Int_t ncell = cluster->GetNCells();
1258  Int_t * sortList = new Int_t [ncell];
1259  Float_t * enerList = new Float_t [ncell];
1260  Double_t * timeList = new Double_t[ncell];
1261 
1262  for(Int_t icell = 0; icell < ncell; icell++)
1263  {
1264  Int_t absId = cluster->GetCellAbsId(icell);
1265  Float_t cellE = cells->GetCellAmplitude(absId);
1266  Double_t cellTime = cells->GetCellTime(absId);
1267 
1268  GetCaloUtils()->GetEMCALRecoUtils()->AcceptCalibrateCell(absId,bc,cellE,cellTime,cells);
1269  cellTime*=1e9;
1270  cellTime-=fConstantTimeShift;
1271 
1272  enerList[icell] = cellE;
1273  timeList[icell] = cellTime;
1274  }
1275 
1276  TMath::Sort(ncell,enerList,sortList);
1277 
1278  //
1279  // Location of max cell
1280  // First in sorted list.
1281  //
1282  Int_t absIdMax = cluster->GetCellAbsId(sortList[0]);
1283  // Float_t cellEMax = enerList[sortList[0]];
1284  Double_t cellTimeMax = timeList[sortList[0]];
1285 
1286  //
1287  // Cluster cell loop, select only secondary cells with enough contribution to cluster
1288  // Start from second highest energy cell
1289  //
1290  largeTime = 0;
1291  Float_t largeTimeE = 0;
1292  for(Int_t icell = 1; icell < ncell; icell++)
1293  {
1294  Int_t absId = cluster->GetCellAbsId(sortList[icell]);
1295  Float_t cellE = enerList[sortList[icell]];
1296  Double_t cellTime = timeList[sortList[icell]];
1297 
1298  //printf("\t icell %d, absId %d, ener %2.2f, time %2.2f, sort %d\n",icell,absId,cellE,cellTime,sortList[icell]);
1299 
1300  if ( absId == absIdMax ) continue;
1301 
1302  Float_t weight = GetCaloUtils()->GetEMCALRecoUtils()->GetCellWeight(cellE,cluster->E());
1303 
1304  //printf("Cluster E %2.2f, cell E %2.2f, time %2.2f, weight %2.3f\n",cluster->E(),cellE,cellTime,weight);
1305  if(weight < 0.01) continue;
1306 
1307  if ( TMath::Abs(cellTime) > 50 )
1308  {
1309  largeTime++;
1310  largeTimeE += cellE;
1311  //printf(" cluster pT %2.2f, cell E %2.2f, large time %d; ",pt, cellE, largeTime);
1312  }
1313  if ( l0bin == -1 ) continue;
1314 
1315  fhTimeLam0BinPerSM [l0bin][sm]->Fill(pt, cellTime, GetEventWeight());
1316  fhTimeLam0BinPerSMWeighted[l0bin][sm]->Fill(pt, cellTime, weight*GetEventWeight());
1317 
1318  fhDTimeLam0BinPerSM [l0bin][sm]->Fill(pt, cellTimeMax-cellTime, GetEventWeight());
1319  fhDTimeLam0BinPerSMWeighted[l0bin][sm]->Fill(pt, cellTimeMax-cellTime, weight*GetEventWeight());
1320 
1321  fhCellClusterEFracLam0BinPerSM [l0bin][sm]->Fill(pt, cellE/cluster->E(), GetEventWeight());
1322 // fhCellClusterEFracLam0BinPerSMWeighted[l0bin][sm]->Fill(pt, cellE/cluster->E(), weight*GetEventWeight());
1323 
1324  fhCellClusterELam0BinPerSM [l0bin][sm]->Fill(pt, cellE, GetEventWeight());
1325  fhCellClusterELam0BinPerSMWeighted[l0bin][sm]->Fill(pt, cellE, weight*GetEventWeight());
1326 
1327  if(ptbin >= 0)
1328  {
1329  fhCellClusterIndexEAndTime[l0bin][ptbin] ->Fill(cellTime, icell , GetEventWeight());
1330  fhCellClusterEAndTime [l0bin][ptbin] ->Fill(cellTime, cellE , GetEventWeight());
1331  fhCellClusterEFracAndTime [l0bin][ptbin] ->Fill(cellTime, cellE/cluster->E(), GetEventWeight());
1332  }
1333 
1334  if ( TMath::Abs(cellTime) > 50 )
1335  {
1336  fhCellClusterEFracLam0BinPerSMLargeTime [l0bin][sm]->Fill(pt, cellE/cluster->E(), GetEventWeight());
1337  fhCellClusterELam0BinPerSMLargeTime [l0bin][sm]->Fill(pt, cellE , GetEventWeight());
1338  fhCellClusterIndexELam0BinPerSMLargeTime[l0bin][sm]->Fill(pt, icell , GetEventWeight());
1339  }
1340 // if(shared)
1341 // fhTimeLam0BinPerSMShared[l0bin][sm]->Fill(pt, cellTime, GetEventWeight());
1342 
1343  if(ptbin >= 0)
1344  {
1345  Int_t icol = -1, icolAbs = -1;
1346  Int_t irow = -1, irowAbs = -1;
1347  Int_t iRCU = -1;
1348  GetModuleNumberCellIndexesAbsCaloMap(absId,GetCalorimeter(), icol, irow, iRCU, icolAbs, irowAbs);
1349 
1350  fhColRowLam0BinPtBin [l0bin][ptbin]->Fill(icolAbs, irowAbs, GetEventWeight());
1351  fhColRowLam0BinPtBinWeighted[l0bin][ptbin]->Fill(icolAbs, irowAbs, weight*GetEventWeight());
1352 
1353  if ( TMath::Abs(cellTime) > 50 )
1354  fhColRowLam0BinPtBinLargeTime[l0bin][ptbin]->Fill(icolAbs, irowAbs, GetEventWeight());
1355 
1356 // if(shared) fhColRowLam0BinPtBinSMShared[l0bin][ptbin]->Fill(icolAbs, irowAbs, GetEventWeight());
1357  }
1358 
1359 
1360  } // cluster cell loop
1361  //if(largeTime > 0) printf("\n");
1362 
1363  if ( largeTime > 0 )
1364  {
1365  if ( l0bin != -1 )
1366  {
1367  if ( ptbin > 0 )
1368  fhEtaPhiLargeTimeInClusterCell[l0bin][ptbin]->Fill(eta, phi, GetEventWeight());
1369 
1370  fhCellClusterEFracLam0BinPerSMLargeTimeTotal[l0bin][sm]->Fill(pt, largeTimeE/cluster->E(), GetEventWeight());
1371  fhNCellsWithLargeTimeInCluster [l0bin][sm]->Fill(pt, largeTime , GetEventWeight());
1372  }
1373 
1374  fhLam0PerSMLargeTimeInClusterCell[sm]->Fill(pt, lambda0, GetEventWeight());
1375  fhLam1PerSMLargeTimeInClusterCell[sm]->Fill(pt, lambda1, GetEventWeight());
1376 
1377  if(largeTime < 6)
1378  {
1379  fhLam0PerNLargeTimeInClusterCell[largeTime-1]->Fill(pt, lambda0, GetEventWeight());
1380  fhLam1PerNLargeTimeInClusterCell[largeTime-1]->Fill(pt, lambda1, GetEventWeight());
1381  }
1382  }
1383 
1384  delete [] sortList ;
1385  delete [] enerList ;
1386  delete [] timeList ;
1387  //printf("Cluster %d, E %2.2f, sm %d, eta %2.2f, phi %2.2f ---> region %d %d\n",cluster->GetID(),cluster->E(),GetModuleNumber(cluster),eta,RadToDeg(phi),etaRegion,phiRegion);
1388  } // Check EMCal SM and regions
1389 
1390  Float_t l0 = 0., l1 = 0.;
1391  Float_t dispp= 0., dEta = 0., dPhi = 0.;
1392  Float_t sEta = 0., sPhi = 0., sEtaPhi = 0.;
1394  {
1395  GetCaloUtils()->GetEMCALRecoUtils()->RecalculateClusterShowerShapeParameters(GetEMCALGeometry(), GetReader()->GetInputEvent()->GetEMCALCells(), cluster,
1396  l0, l1, dispp, dEta, dPhi, sEta, sPhi, sEtaPhi);
1397  //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",
1398  // l0, l1, dispp, dEta, dPhi, sEta, sPhi, sEtaPhi );
1399  //printf("AliAnaPhoton::FillShowerShapeHistogram - dispersion %f, dispersion eta+phi %f \n",
1400  // disp, dPhi+dEta );
1401  fhDispEtaE -> Fill(energy, dEta , GetEventWeight());
1402  fhDispPhiE -> Fill(energy, dPhi , GetEventWeight());
1403  fhSumEtaE -> Fill(energy, sEta , GetEventWeight());
1404  fhSumPhiE -> Fill(energy, sPhi , GetEventWeight());
1405  fhSumEtaPhiE -> Fill(energy, sEtaPhi , GetEventWeight());
1406  fhDispEtaPhiDiffE -> Fill(energy, dPhi-dEta, GetEventWeight());
1407 
1408  if(dEta+dPhi>0)fhSphericityE -> Fill(energy,(dPhi-dEta)/(dEta+dPhi) , GetEventWeight());
1409  if(dEta+sEta>0)fhDispSumEtaDiffE -> Fill(energy,(dEta-sEta)/((dEta+sEta)/2.), GetEventWeight());
1410  if(dPhi+sPhi>0)fhDispSumPhiDiffE -> Fill(energy,(dPhi-sPhi)/((dPhi+sPhi)/2.), GetEventWeight());
1411 
1412  Int_t ebin = -1;
1413  if (energy < 2 ) ebin = 0;
1414  else if (energy < 4 ) ebin = 1;
1415  else if (energy < 6 ) ebin = 2;
1416  else if (energy < 10) ebin = 3;
1417  else if (energy < 15) ebin = 4;
1418  else if (energy < 20) ebin = 5;
1419  else ebin = 6;
1420 
1421  fhDispEtaDispPhi[ebin]->Fill(dEta , dPhi, GetEventWeight());
1422  fhLambda0DispEta[ebin]->Fill(lambda0, dEta, GetEventWeight());
1423  fhLambda0DispPhi[ebin]->Fill(lambda0, dPhi, GetEventWeight());
1424  }
1425 
1426  // If track-matching was off, check effect of track-matching residual cut
1427 
1428  if(!fRejectTrackMatch)
1429  {
1430  Float_t dZ = cluster->GetTrackDz();
1431  Float_t dR = cluster->GetTrackDx();
1432 // if(cluster->IsEMCAL() && GetCaloUtils()->IsRecalculationOfClusterTrackMatchingOn())
1433 // {
1434 // dR = 2000., dZ = 2000.;
1435 // GetCaloUtils()->GetEMCALRecoUtils()->GetMatchedResiduals(cluster->GetID(),dZ,dR);
1436 // }
1437 
1438  if(TMath::Abs(dZ) < 0.05 && TMath::Abs(dR) < 0.05)
1439  {
1440  fhLam0ETM ->Fill(energy, lambda0, GetEventWeight());
1441  fhLam0PtTM->Fill(pt , lambda0, GetEventWeight());
1442  fhLam1ETM ->Fill(energy, lambda1, GetEventWeight());
1443  fhDispETM ->Fill(energy, disp , GetEventWeight());
1444 
1445  if(GetCalorimeter() == kEMCAL && GetFirstSMCoveredByTRD() >= 0 &&
1446  GetModuleNumber(cluster) >= GetFirstSMCoveredByTRD() )
1447  {
1448  fhLam0ETMTRD ->Fill(energy, lambda0, GetEventWeight());
1449  fhLam0PtTMTRD->Fill(pt , lambda0, GetEventWeight());
1450  fhLam1ETMTRD ->Fill(energy, lambda1, GetEventWeight());
1451  fhDispETMTRD ->Fill(energy, disp , GetEventWeight());
1452  }
1453  }
1454  } // If track-matching was off, check effect of matching residual cut
1455 
1457  {
1458  if(energy < 2)
1459  {
1460  fhNCellsLam0LowE ->Fill(ncells, lambda0, GetEventWeight());
1461  fhNCellsLam1LowE ->Fill(ncells, lambda1, GetEventWeight());
1462  fhNCellsDispLowE ->Fill(ncells, disp , GetEventWeight());
1463 
1464  fhLam1Lam0LowE ->Fill(lambda1, lambda0, GetEventWeight());
1465  fhLam0DispLowE ->Fill(lambda0, disp , GetEventWeight());
1466  fhDispLam1LowE ->Fill(disp , lambda1, GetEventWeight());
1467  fhEtaLam0LowE ->Fill(eta , lambda0, GetEventWeight());
1468  fhPhiLam0LowE ->Fill(phi , lambda0, GetEventWeight());
1469  }
1470  else
1471  {
1472  fhNCellsLam0HighE ->Fill(ncells, lambda0, GetEventWeight());
1473  fhNCellsLam1HighE ->Fill(ncells, lambda1, GetEventWeight());
1474  fhNCellsDispHighE ->Fill(ncells, disp , GetEventWeight());
1475 
1476  fhLam1Lam0HighE ->Fill(lambda1, lambda0, GetEventWeight());
1477  fhLam0DispHighE ->Fill(lambda0, disp , GetEventWeight());
1478  fhDispLam1HighE ->Fill(disp , lambda1, GetEventWeight());
1479  fhEtaLam0HighE ->Fill(eta , lambda0, GetEventWeight());
1480  fhPhiLam0HighE ->Fill(phi , lambda0, GetEventWeight());
1481  }
1482  }
1483 
1484  if(IsDataMC())
1485  {
1486  AliVCaloCells* cells = 0;
1487  if(GetCalorimeter() == kEMCAL) cells = GetEMCALCells();
1488  else cells = GetPHOSCells();
1489 
1490  // Fill histograms to check shape of embedded clusters
1491  Float_t fraction = 0;
1492  // printf("check embedding %i\n",GetReader()->IsEmbeddedClusterSelectionOn());
1493 
1494  if(GetReader()->IsEmbeddedClusterSelectionOn())
1495  {
1496  // Only working for EMCAL
1497  // printf("embedded\n");
1498 
1499  Float_t clusterE = 0; // recalculate in case corrections applied.
1500  Float_t cellE = 0;
1501  for(Int_t icell = 0; icell < cluster->GetNCells(); icell++)
1502  {
1503  cellE = cells->GetCellAmplitude(cluster->GetCellAbsId(icell));
1504  clusterE+=cellE;
1505  fraction+=cellE*cluster->GetCellAmplitudeFraction(icell);
1506  }
1507 
1508  //Fraction of total energy due to the embedded signal
1509  fraction/=clusterE;
1510 
1511  AliDebug(1,Form("Energy fraction of embedded signal %2.3f, Energy %2.3f",fraction, clusterE));
1512 
1513  fhEmbeddedSignalFractionEnergy->Fill(clusterE, fraction, GetEventWeight());
1514  } // embedded fraction
1515 
1516  // Check the origin and fill histograms
1517 
1518  Int_t mcIndex = -1;
1519 
1520  if( GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPhoton) &&
1521  !GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCConversion) &&
1522  !GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPi0) &&
1523  !GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCEta))
1524  {
1525  mcIndex = kmcssPhoton ;
1526 
1527  if(!GetReader()->IsEmbeddedClusterSelectionOn())
1528  {
1529  //Check particle overlaps in cluster
1530 
1531  // Compare the primary depositing more energy with the rest,
1532  // if no photon/electron as comon ancestor (conversions), count as other particle
1533  const UInt_t nlabels = cluster->GetNLabels();
1534  Int_t overpdg[nlabels];
1535  Int_t overlab[nlabels];
1536  Int_t noverlaps = GetMCAnalysisUtils()->GetNOverlaps(cluster->GetLabels(), nlabels,mcTag,-1,GetReader(),overpdg,overlab);
1537 
1538  //printf("N overlaps %d \n",noverlaps);
1539 
1540  if(noverlaps == 0)
1541  {
1542  fhMCPhotonELambda0NoOverlap ->Fill(energy, lambda0, GetEventWeight());
1543  }
1544  else if(noverlaps == 1)
1545  {
1546  fhMCPhotonELambda0TwoOverlap ->Fill(energy, lambda0, GetEventWeight());
1547  }
1548  else if(noverlaps > 1)
1549  {
1550  fhMCPhotonELambda0NOverlap ->Fill(energy, lambda0, GetEventWeight());
1551  }
1552  else
1553  {
1554  AliWarning(Form("n overlaps = %d!!", noverlaps));
1555  }
1556  } // No embedding
1557 
1558  // Fill histograms to check shape of embedded clusters
1559  if(GetReader()->IsEmbeddedClusterSelectionOn())
1560  {
1561  if (fraction > 0.9)
1562  {
1563  fhEmbedPhotonELambda0FullSignal ->Fill(energy, lambda0, GetEventWeight());
1564  }
1565  else if(fraction > 0.5)
1566  {
1567  fhEmbedPhotonELambda0MostlySignal ->Fill(energy, lambda0, GetEventWeight());
1568  }
1569  else if(fraction > 0.1)
1570  {
1571  fhEmbedPhotonELambda0MostlyBkg ->Fill(energy, lambda0, GetEventWeight());
1572  }
1573  else
1574  {
1575  fhEmbedPhotonELambda0FullBkg ->Fill(energy, lambda0, GetEventWeight());
1576  }
1577  } // embedded
1578  } // photon no conversion
1579  else if ( GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPhoton) &&
1580  GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCConversion) &&
1581  !GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPi0) &&
1582  !GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCEta))
1583  {
1584  mcIndex = kmcssConversion ;
1585  }//conversion photon
1586 
1587  else if ( GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCElectron))
1588  {
1589  mcIndex = kmcssElectron ;
1590  }//electron
1591  else if ( GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPi0) )
1592  {
1593  mcIndex = kmcssPi0 ;
1594 
1595  //Fill histograms to check shape of embedded clusters
1596  if(GetReader()->IsEmbeddedClusterSelectionOn())
1597  {
1598  if (fraction > 0.9)
1599  {
1600  fhEmbedPi0ELambda0FullSignal ->Fill(energy, lambda0, GetEventWeight());
1601  }
1602  else if(fraction > 0.5)
1603  {
1604  fhEmbedPi0ELambda0MostlySignal ->Fill(energy, lambda0, GetEventWeight());
1605  }
1606  else if(fraction > 0.1)
1607  {
1608  fhEmbedPi0ELambda0MostlyBkg ->Fill(energy, lambda0, GetEventWeight());
1609  }
1610  else
1611  {
1612  fhEmbedPi0ELambda0FullBkg ->Fill(energy, lambda0, GetEventWeight());
1613  }
1614  } // embedded
1615 
1616  }//pi0
1617  else if ( GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCEta) )
1618  {
1619  mcIndex = kmcssEta ;
1620  }//eta
1621  else
1622  {
1623  mcIndex = kmcssOther ;
1624  }//other particles
1625 
1626  fhMCELambda0 [mcIndex]->Fill(energy, lambda0, GetEventWeight());
1627  fhMCPtLambda0 [mcIndex]->Fill(pt , lambda0, GetEventWeight());
1628  fhMCELambda1 [mcIndex]->Fill(energy, lambda1, GetEventWeight());
1629  fhMCEDispersion [mcIndex]->Fill(energy, disp , GetEventWeight());
1630  fhMCNCellsE [mcIndex]->Fill(energy, ncells , GetEventWeight());
1631  fhMCMaxCellDiffClusterE[mcIndex]->Fill(energy, maxCellFraction, GetEventWeight());
1632 
1634  {
1635  if (energy < 2.)
1636  {
1637  fhMCLambda0vsClusterMaxCellDiffE0[mcIndex]->Fill(lambda0, maxCellFraction, GetEventWeight());
1638  fhMCNCellsvsClusterMaxCellDiffE0 [mcIndex]->Fill(ncells, maxCellFraction, GetEventWeight());
1639  }
1640  else if(energy < 6.)
1641  {
1642  fhMCLambda0vsClusterMaxCellDiffE2[mcIndex]->Fill(lambda0, maxCellFraction, GetEventWeight());
1643  fhMCNCellsvsClusterMaxCellDiffE2 [mcIndex]->Fill(ncells, maxCellFraction, GetEventWeight());
1644  }
1645  else
1646  {
1647  fhMCLambda0vsClusterMaxCellDiffE6[mcIndex]->Fill(lambda0, maxCellFraction, GetEventWeight());
1648  fhMCNCellsvsClusterMaxCellDiffE6 [mcIndex]->Fill(ncells, maxCellFraction, GetEventWeight());
1649  }
1650 
1651  if(GetCalorimeter() == kEMCAL)
1652  {
1653  fhMCEDispEta [mcIndex]-> Fill(energy, dEta , GetEventWeight());
1654  fhMCEDispPhi [mcIndex]-> Fill(energy, dPhi , GetEventWeight());
1655  fhMCESumEtaPhi [mcIndex]-> Fill(energy, sEtaPhi , GetEventWeight());
1656  fhMCEDispEtaPhiDiff [mcIndex]-> Fill(energy, dPhi-dEta, GetEventWeight());
1657 
1658  if( dEta+dPhi > 0 ) fhMCESphericity[mcIndex]-> Fill(energy, (dPhi-dEta)/(dEta+dPhi), GetEventWeight());
1659 
1660  Int_t ebin = -1;
1661  if (energy < 2 ) ebin = 0;
1662  else if (energy < 4 ) ebin = 1;
1663  else if (energy < 6 ) ebin = 2;
1664  else if (energy < 10) ebin = 3;
1665  else if (energy < 15) ebin = 4;
1666  else if (energy < 20) ebin = 5;
1667  else ebin = 6;
1668 
1669  fhMCDispEtaDispPhi[ebin][mcIndex]->Fill(dEta , dPhi, GetEventWeight());
1670  fhMCLambda0DispEta[ebin][mcIndex]->Fill(lambda0, dEta, GetEventWeight());
1671  fhMCLambda0DispPhi[ebin][mcIndex]->Fill(lambda0, dPhi, GetEventWeight());
1672  }
1673  }
1674  }// MC data
1675 }
1676 
1677 //__________________________________________________________________________
1683 //__________________________________________________________________________
1685  Int_t cut)
1686 {
1687  Float_t dZ = cluster->GetTrackDz();
1688  Float_t dR = cluster->GetTrackDx();
1689 
1690 // if(cluster->IsEMCAL() && GetCaloUtils()->IsRecalculationOfClusterTrackMatchingOn())
1691 // {
1692 // dR = 2000., dZ = 2000.;
1693 // GetCaloUtils()->GetEMCALRecoUtils()->GetMatchedResiduals(cluster->GetID(),dZ,dR);
1694 // }
1695 
1696  AliVTrack *track = GetCaloUtils()->GetMatchedTrack(cluster, GetReader()->GetInputEvent());
1697 
1698  Bool_t positive = kFALSE;
1699  if(track) positive = (track->Charge()>0);
1700 
1701  if(fhTrackMatchedDEta[cut] && TMath::Abs(dR) < 999)
1702  {
1703  fhTrackMatchedDEta[cut]->Fill(cluster->E(), dZ, GetEventWeight());
1704  fhTrackMatchedDPhi[cut]->Fill(cluster->E(), dR, GetEventWeight());
1705  if(cluster->E() > 0.5) fhTrackMatchedDEtaDPhi[cut]->Fill(dZ, dR, GetEventWeight());
1706 
1707  if(track)
1708  {
1709  if(positive)
1710  {
1711  fhTrackMatchedDEtaPos[cut]->Fill(cluster->E(), dZ, GetEventWeight());
1712  fhTrackMatchedDPhiPos[cut]->Fill(cluster->E(), dR, GetEventWeight());
1713  if(cluster->E() > 0.5) fhTrackMatchedDEtaDPhiPos[cut]->Fill(dZ,dR, GetEventWeight());
1714  }
1715  else
1716  {
1717  fhTrackMatchedDEtaNeg[cut]->Fill(cluster->E(), dZ, GetEventWeight());
1718  fhTrackMatchedDPhiNeg[cut]->Fill(cluster->E(), dR, GetEventWeight());
1719  if(cluster->E() > 0.5) fhTrackMatchedDEtaDPhiNeg[cut]->Fill(dZ, dR, GetEventWeight());
1720  }
1721  }
1722 
1723  Int_t nSMod = GetModuleNumber(cluster);
1724 
1725  if( GetCalorimeter() == kEMCAL && GetFirstSMCoveredByTRD() >= 0 &&
1726  nSMod >= GetFirstSMCoveredByTRD() )
1727  {
1728  fhTrackMatchedDEtaTRD[cut]->Fill(cluster->E(), dZ, GetEventWeight());
1729  fhTrackMatchedDPhiTRD[cut]->Fill(cluster->E(), dR, GetEventWeight());
1730  }
1731 
1732  // Check dEdx and E/p of matched clusters
1733 
1734  if(TMath::Abs(dZ) < 0.05 && TMath::Abs(dR) < 0.05)
1735  {
1736  if(track)
1737  {
1738  Float_t dEdx = track->GetTPCsignal();
1739  Float_t eOverp = cluster->E()/track->P();
1740 
1741  fhdEdx[cut] ->Fill(cluster->E(), dEdx , GetEventWeight());
1742  fhEOverP[cut]->Fill(cluster->E(), eOverp, GetEventWeight());
1743 
1744  if(GetCalorimeter()==kEMCAL && GetFirstSMCoveredByTRD() >= 0 &&
1745  nSMod >= GetFirstSMCoveredByTRD() )
1746  fhEOverPTRD[cut]->Fill(cluster->E(), eOverp, GetEventWeight());
1747 
1748 
1749  }
1750  else
1751  AliWarning(Form("Residual OK but (dR, dZ)= (%2.4f,%2.4f) no track associated WHAT?", dR,dZ));
1752 
1753  if(IsDataMC())
1754  {
1755 
1756  Int_t tag = GetMCAnalysisUtils()->CheckOrigin(cluster->GetLabels(),cluster->GetNLabels(),GetReader(),GetCalorimeter());
1757 
1759  {
1762  fhTrackMatchedMCParticle[cut]->Fill(cluster->E(), 2.5, GetEventWeight());
1764  fhTrackMatchedMCParticle[cut]->Fill(cluster->E(), 0.5, GetEventWeight());
1766  fhTrackMatchedMCParticle[cut]->Fill(cluster->E(), 1.5, GetEventWeight());
1767  else
1768  fhTrackMatchedMCParticle[cut]->Fill(cluster->E(), 3.5, GetEventWeight());
1769 
1770  // Check if several particles contributed to cluster and discard overlapped mesons
1773  {
1774  if(cluster->GetNLabels()==1)
1775  {
1776  fhTrackMatchedDEtaMCNoOverlap[cut]->Fill(cluster->E(), dZ, GetEventWeight());
1777  fhTrackMatchedDPhiMCNoOverlap[cut]->Fill(cluster->E(), dR, GetEventWeight());
1778  }
1779  else
1780  {
1781  fhTrackMatchedDEtaMCOverlap[cut]->Fill(cluster->E(), dZ, GetEventWeight());
1782  fhTrackMatchedDPhiMCOverlap[cut]->Fill(cluster->E(), dR, GetEventWeight());
1783  }
1784 
1785  }// Check overlaps
1786 
1787  }
1788  else
1789  {
1792  fhTrackMatchedMCParticle[cut]->Fill(cluster->E(), 6.5, GetEventWeight());
1794  fhTrackMatchedMCParticle[cut]->Fill(cluster->E(), 4.5, GetEventWeight());
1796  fhTrackMatchedMCParticle[cut]->Fill(cluster->E(), 5.5, GetEventWeight());
1797  else
1798  fhTrackMatchedMCParticle[cut]->Fill(cluster->E(), 7.5, GetEventWeight());
1799 
1800  // Check if several particles contributed to cluster
1803  {
1804  fhTrackMatchedDEtaMCConversion[cut]->Fill(cluster->E(), dZ, GetEventWeight());
1805  fhTrackMatchedDPhiMCConversion[cut]->Fill(cluster->E(), dR, GetEventWeight());
1806 
1807  }// Check overlaps
1808 
1809  }
1810 
1811  } // MC
1812 
1813  } // residuals window
1814 
1815  } // Small residual
1816 }
1817 
1818 //___________________________________________
1820 //___________________________________________
1822 {
1823  TString parList ; //this will be list of parameters used for this analysis.
1824  const Int_t buffersize = 255;
1825  char onePar[buffersize] ;
1826 
1827  snprintf(onePar,buffersize,"--- AliAnaPhoton ---:") ;
1828  parList+=onePar ;
1829  snprintf(onePar,buffersize,"Calorimeter: %s;",GetCalorimeterString().Data()) ;
1830  parList+=onePar ;
1831  snprintf(onePar,buffersize,"fMinDist =%2.2f (Minimal distance to bad channel to accept cluster);",fMinDist) ;
1832  parList+=onePar ;
1833  snprintf(onePar,buffersize,"fMinDist2=%2.2f (Cuts on Minimal distance to study acceptance evaluation);",fMinDist2) ;
1834  parList+=onePar ;
1835  snprintf(onePar,buffersize,"fMinDist3=%2.2f (One more cut on distance used for acceptance-efficiency study);",fMinDist3) ;
1836  parList+=onePar ;
1837  snprintf(onePar,buffersize,"fRejectTrackMatch: %d",fRejectTrackMatch) ;
1838  parList+=onePar ;
1839 
1840  // Get parameters set in base class.
1841  parList += GetBaseParametersList() ;
1842 
1843  // Get parameters set in PID class.
1844  parList += GetCaloPID()->GetPIDParametersList() ;
1845 
1846  // Get parameters set in FiducialCut class (not available yet)
1847  //parlist += GetFidCut()->GetFidCutParametersList()
1848 
1849  return new TObjString(parList) ;
1850 }
1851 
1852 //_____________________________________________
1855 //_____________________________________________
1857 {
1858  TList * outputContainer = new TList() ;
1859  outputContainer->SetName("PhotonHistos") ;
1860 
1867 
1874 
1875  Int_t ndedxbins = GetHistogramRanges()->GetHistodEdxBins();
1878  Int_t nPoverEbins = GetHistogramRanges()->GetHistoPOverEBins();
1879  Float_t pOverEmax = GetHistogramRanges()->GetHistoPOverEMax();
1880  Float_t pOverEmin = GetHistogramRanges()->GetHistoPOverEMin();
1881 
1882  Int_t nratbins = GetHistogramRanges()->GetHistoRatioBins();
1885  Int_t ndifbins = GetHistogramRanges()->GetHistoEDiffBins();
1888 
1889 
1890  Int_t bin[] = {0,2,4,6,10,15,20,100}; // energy bins for SS studies
1891 
1892  TString cut[] = {"Open","Reader","E","Time","NCells","NLM","Fidutial","Matching","Bad","PID"};
1893  for (Int_t i = 0; i < 10 ; i++)
1894  {
1895  fhClusterCutsE[i] = new TH1F(Form("hE_Cut_%d_%s", i, cut[i].Data()),
1896  Form("Number of clusters that pass cuts <= %d, %s", i, cut[i].Data()),
1897  nptbins,ptmin,ptmax);
1898  fhClusterCutsE[i]->SetYTitle("d#it{N}/d#it{E} ");
1899  fhClusterCutsE[i]->SetXTitle("#it{E} (GeV)");
1900  outputContainer->Add(fhClusterCutsE[i]) ;
1901 
1902  fhClusterCutsPt[i] = new TH1F(Form("hPt_Cut_%d_%s", i, cut[i].Data()),
1903  Form("Number of clusters that pass cuts <= %d, %s", i, cut[i].Data()),
1904  nptbins,ptmin,ptmax);
1905  fhClusterCutsPt[i]->SetYTitle("d#it{N}/d#it{E} ");
1906  fhClusterCutsPt[i]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
1907  outputContainer->Add(fhClusterCutsPt[i]) ;
1908  }
1909 
1910  fhEClusterSM = new TH2F("hEClusterSM","Raw clusters E and super-module number",
1911  nptbins,ptmin,ptmax,
1912  GetCaloUtils()->GetNumberOfSuperModulesUsed(),0,GetCaloUtils()->GetNumberOfSuperModulesUsed());
1913  fhEClusterSM->SetYTitle("SuperModule ");
1914  fhEClusterSM->SetXTitle("#it{E} (GeV)");
1915  outputContainer->Add(fhEClusterSM) ;
1916 
1917  fhPtClusterSM = new TH2F("hPtClusterSM","Raw clusters #it{p}_{T} and super-module number",
1918  nptbins,ptmin,ptmax,
1919  GetCaloUtils()->GetNumberOfSuperModulesUsed(),0,GetCaloUtils()->GetNumberOfSuperModulesUsed());
1920  fhPtClusterSM->SetYTitle("SuperModule ");
1921  fhPtClusterSM->SetXTitle("#it{E} (GeV)");
1922  outputContainer->Add(fhPtClusterSM) ;
1923 
1924  fhEPhotonSM = new TH2F("hEPhotonSM","Selected clusters E and super-module number",
1925  nptbins,ptmin,ptmax,
1926  GetCaloUtils()->GetNumberOfSuperModulesUsed(),0,GetCaloUtils()->GetNumberOfSuperModulesUsed());
1927  fhEPhotonSM->SetYTitle("SuperModule ");
1928  fhEPhotonSM->SetXTitle("#it{E} (GeV)");
1929  outputContainer->Add(fhEPhotonSM) ;
1930 
1931  fhPtPhotonSM = new TH2F("hPtPhotonSM","Selected clusters #it{p}_{T} and super-module number",
1932  nptbins,ptmin,ptmax,
1933  GetCaloUtils()->GetNumberOfSuperModulesUsed(),0,GetCaloUtils()->GetNumberOfSuperModulesUsed());
1934  fhPtPhotonSM->SetYTitle("SuperModule ");
1935  fhPtPhotonSM->SetXTitle("#it{p}_{T} (GeV/#it{c})");
1936  outputContainer->Add(fhPtPhotonSM) ;
1937 
1938  fhNCellsE = new TH2F ("hNCellsE","# of cells in cluster vs E of clusters", nptbins,ptmin,ptmax, nbins,nmin,nmax);
1939  fhNCellsE->SetXTitle("#it{E} (GeV)");
1940  fhNCellsE->SetYTitle("# of cells in cluster");
1941  outputContainer->Add(fhNCellsE);
1942 
1943  fhCellsE = new TH2F ("hCellsE","energy of cells in cluster vs E of clusters", nptbins,ptmin,ptmax, nptbins*2,ptmin,ptmax);
1944  fhCellsE->SetXTitle("#it{E}_{cluster} (GeV)");
1945  fhCellsE->SetYTitle("#it{E}_{cell} (GeV)");
1946  outputContainer->Add(fhCellsE);
1947 
1948  fhTimePt = new TH2F ("hTimePt","time of cluster vs pT of clusters", nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
1949  fhTimePt->SetXTitle("#it{p}_{T} (GeV/#it{c})");
1950  fhTimePt->SetYTitle("#it{time} (ns)");
1951  outputContainer->Add(fhTimePt);
1952 
1953  fhMaxCellDiffClusterE = new TH2F ("hMaxCellDiffClusterE","energy vs difference of cluster energy - max cell energy / cluster energy, good clusters",
1954  nptbins,ptmin,ptmax, 500,0,1.);
1955  fhMaxCellDiffClusterE->SetXTitle("#it{E}_{cluster} (GeV) ");
1956  fhMaxCellDiffClusterE->SetYTitle("(#it{E}_{cluster} - #it{E}_{cell max})/ #it{E}_{cluster}");
1957  outputContainer->Add(fhMaxCellDiffClusterE);
1958 
1959  fhEPhoton = new TH1F("hEPhoton","Number of #gamma over calorimeter vs energy",nptbins,ptmin,ptmax);
1960  fhEPhoton->SetYTitle("#it{counts}");
1961  fhEPhoton->SetXTitle("#it{E}_{#gamma}(GeV)");
1962  outputContainer->Add(fhEPhoton) ;
1963 
1964  fhPtPhoton = new TH1F("hPtPhoton","Number of #gamma over calorimeter vs #it{p}_{T}",nptbins,ptmin,ptmax);
1965  fhPtPhoton->SetYTitle("#it{counts}");
1966  fhPtPhoton->SetXTitle("p_{T #gamma}(GeV/#it{c})");
1967  outputContainer->Add(fhPtPhoton) ;
1968 
1970  {
1971  fhPtCentralityPhoton = new TH2F("hPtCentralityPhoton","centrality vs #it{p}_{T}",nptbins,ptmin,ptmax, 100,0,100);
1972  fhPtCentralityPhoton->SetYTitle("Centrality");
1973  fhPtCentralityPhoton->SetXTitle("#it{p}_{T}(GeV/#it{c})");
1974  outputContainer->Add(fhPtCentralityPhoton) ;
1975 
1976  fhPtEventPlanePhoton = new TH2F("hPtEventPlanePhoton","centrality vs #it{p}_{T}",nptbins,ptmin,ptmax, 100,0,TMath::Pi());
1977  fhPtEventPlanePhoton->SetYTitle("Event plane angle (rad)");
1978  fhPtEventPlanePhoton->SetXTitle("#it{p}_{T} (GeV/#it{c})");
1979  outputContainer->Add(fhPtEventPlanePhoton) ;
1980  }
1981 
1982  fhEtaPhi = new TH2F
1983  ("hEtaPhi","cluster,#it{E} > 0.5 GeV, #eta vs #phi",netabins,etamin,etamax,nphibins,phimin,phimax);
1984  fhEtaPhi->SetYTitle("#phi (rad)");
1985  fhEtaPhi->SetXTitle("#eta");
1986  outputContainer->Add(fhEtaPhi) ;
1987 
1988  fhPhiPhoton = new TH2F
1989  ("hPhiPhoton","#phi_{#gamma} vs #it{p}_{T}",nptbins,ptmin,ptmax,nphibins,phimin,phimax);
1990  fhPhiPhoton->SetYTitle("#phi (rad)");
1991  fhPhiPhoton->SetXTitle("p_{T #gamma} (GeV/#it{c})");
1992  outputContainer->Add(fhPhiPhoton) ;
1993 
1994  fhEtaPhoton = new TH2F
1995  ("hEtaPhoton","#eta_{#gamma} vs #it{p}_{T}",nptbins,ptmin,ptmax,netabins,etamin,etamax);
1996  fhEtaPhoton->SetYTitle("#eta");
1997  fhEtaPhoton->SetXTitle("p_{T #gamma} (GeV/#it{c})");
1998  outputContainer->Add(fhEtaPhoton) ;
1999 
2000  fhEtaPhiPhoton = new TH2F
2001  ("hEtaPhiPhoton","#eta vs #phi",netabins,etamin,etamax,nphibins,phimin,phimax);
2002  fhEtaPhiPhoton->SetYTitle("#phi (rad)");
2003  fhEtaPhiPhoton->SetXTitle("#eta");
2004  outputContainer->Add(fhEtaPhiPhoton) ;
2005  if(GetMinPt() < 0.5)
2006  {
2007  fhEtaPhi05Photon = new TH2F
2008  ("hEtaPhi05Photon","#eta vs #phi, E < 0.5",netabins,etamin,etamax,nphibins,phimin,phimax);
2009  fhEtaPhi05Photon->SetYTitle("#phi (rad)");
2010  fhEtaPhi05Photon->SetXTitle("#eta");
2011  outputContainer->Add(fhEtaPhi05Photon) ;
2012  }
2013 
2014  fhNLocMax = new TH2F("hNLocMax","Number of local maxima in cluster",
2015  nptbins,ptmin,ptmax,10,0,10);
2016  fhNLocMax ->SetYTitle("N maxima");
2017  fhNLocMax ->SetXTitle("#it{E} (GeV)");
2018  outputContainer->Add(fhNLocMax) ;
2019 
2020  //Shower shape
2021  if(fFillSSHistograms)
2022  {
2023  fhLam0E = new TH2F ("hLam0E","#lambda_{0}^{2} vs E", nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
2024  fhLam0E->SetYTitle("#lambda_{0}^{2}");
2025  fhLam0E->SetXTitle("#it{E} (GeV)");
2026  outputContainer->Add(fhLam0E);
2027 
2028  fhLam0Pt = new TH2F ("hLam0Pt","#lambda_{0}^{2} vs #it{p}_{T}", nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
2029  fhLam0Pt->SetYTitle("#lambda_{0}^{2}");
2030  fhLam0Pt->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2031  outputContainer->Add(fhLam0Pt);
2032 
2033  fhLam1E = new TH2F ("hLam1E","#lambda_{1}^{2} vs E", nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
2034  fhLam1E->SetYTitle("#lambda_{1}^{2}");
2035  fhLam1E->SetXTitle("#it{E} (GeV)");
2036  outputContainer->Add(fhLam1E);
2037 
2038  fhLam1Pt = new TH2F ("hLam1Pt","#lambda_{1}^{2} vs E", nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
2039  fhLam1Pt->SetYTitle("#lambda_{1}^{2}");
2040  fhLam1Pt->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2041  outputContainer->Add(fhLam1Pt);
2042 
2043  fhDispE = new TH2F ("hDispE"," dispersion^{2} vs E", nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
2044  fhDispE->SetYTitle("D^{2}");
2045  fhDispE->SetXTitle("#it{E} (GeV) ");
2046  outputContainer->Add(fhDispE);
2047 
2048  fhDispPt = new TH2F ("hDispPt"," dispersion^{2} vs #it{p}_{T}", nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
2049  fhDispPt->SetYTitle("D^{2}");
2050  fhDispPt->SetXTitle("#it{p}_{T} (GeV/#it{c}) ");
2051  outputContainer->Add(fhDispPt);
2052 
2054  {
2055  fhLam0PtNLM1 = new TH2F ("hLam0PtNLM1","#lambda_{0}^{2} vs #it{p}_{T}, #it{n}_{LM}=1", nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
2056  fhLam0PtNLM1->SetYTitle("#lambda_{0}^{2}");
2057  fhLam0PtNLM1->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2058  outputContainer->Add(fhLam0PtNLM1);
2059 
2060  fhLam0PtNLM2 = new TH2F ("hLam0PtNLM2","#lambda_{0}^{2} vs #it{p}_{T}, #it{n}_{LM}=2", nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
2061  fhLam0PtNLM2->SetYTitle("#lambda_{0}^{2}");
2062  fhLam0PtNLM2->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2063  outputContainer->Add(fhLam0PtNLM2);
2064 
2065  fhLam1PtNLM1 = new TH2F ("hLam1PtNLM1","#lambda_{1}^{2} vs #it{p}_{T}, #it{n}_{LM}=1", nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
2066  fhLam1PtNLM1->SetYTitle("#lambda_{1}^{2}");
2067  fhLam1PtNLM1->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2068  outputContainer->Add(fhLam1PtNLM1);
2069 
2070  fhLam1PtNLM2 = new TH2F ("hLam1PtNLM2","#lambda_{1}^{2} vs #it{p}_{T}, #it{n}_{LM}=2", nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
2071  fhLam1PtNLM2->SetYTitle("#lambda_{1}^{2}");
2072  fhLam1PtNLM2->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2073  outputContainer->Add(fhLam1PtNLM2);
2074  }
2075 
2076  if(!fRejectTrackMatch)
2077  {
2078  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);
2079  fhLam0ETM->SetYTitle("#lambda_{0}^{2}");
2080  fhLam0ETM->SetXTitle("#it{E} (GeV)");
2081  outputContainer->Add(fhLam0ETM);
2082 
2083  fhLam0PtTM = new TH2F ("hLam0PtTM","#lambda_{0}^{2} vs #it{p}_{T}, cut on track-matching residual |#Delta #eta| < 0.05, |#Delta #phi| < 0.05", nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
2084  fhLam0PtTM->SetYTitle("#lambda_{0}^{2}");
2085  fhLam0PtTM->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2086  outputContainer->Add(fhLam0PtTM);
2087 
2088  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);
2089  fhLam1ETM->SetYTitle("#lambda_{1}^{2}");
2090  fhLam1ETM->SetXTitle("#it{E} (GeV)");
2091  outputContainer->Add(fhLam1ETM);
2092 
2093  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);
2094  fhDispETM->SetYTitle("D^{2}");
2095  fhDispETM->SetXTitle("#it{E} (GeV) ");
2096  outputContainer->Add(fhDispETM);
2097  }
2098 
2099  if(GetCalorimeter() == kEMCAL && GetFirstSMCoveredByTRD() >= 0)
2100  {
2101  fhLam0ETRD = new TH2F ("hLam0ETRD","#lambda_{0}^{2} vs E, EMCAL SM covered by TRD", nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
2102  fhLam0ETRD->SetYTitle("#lambda_{0}^{2}");
2103  fhLam0ETRD->SetXTitle("#it{E} (GeV)");
2104  outputContainer->Add(fhLam0ETRD);
2105 
2106  fhLam0PtTRD = new TH2F ("hLam0PtTRD","#lambda_{0}^{2} vs #it{p}_{T}, EMCAL SM covered by TRD", nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
2107  fhLam0PtTRD->SetYTitle("#lambda_{0}^{2}");
2108  fhLam0PtTRD->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2109  outputContainer->Add(fhLam0PtTRD);
2110 
2111  fhLam1ETRD = new TH2F ("hLam1ETRD","#lambda_{1}^{2} vs E, EMCAL SM covered by TRD", nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
2112  fhLam1ETRD->SetYTitle("#lambda_{1}^{2}");
2113  fhLam1ETRD->SetXTitle("#it{E} (GeV)");
2114  outputContainer->Add(fhLam1ETRD);
2115 
2116  fhDispETRD = new TH2F ("hDispETRD"," dispersion^{2} vs E, EMCAL SM covered by TRD", nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
2117  fhDispETRD->SetYTitle("Dispersion^{2}");
2118  fhDispETRD->SetXTitle("#it{E} (GeV) ");
2119  outputContainer->Add(fhDispETRD);
2120 
2122  {
2123  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);
2124  fhLam0ETMTRD->SetYTitle("#lambda_{0}^{2}");
2125  fhLam0ETMTRD->SetXTitle("#it{E} (GeV)");
2126  outputContainer->Add(fhLam0ETMTRD);
2127 
2128  fhLam0PtTMTRD = new TH2F ("hLam0PtTMTRD","#lambda_{0}^{2} vs #it{p}_{T}, EMCAL SM covered by TRD, cut on track-matching residual |#Delta #eta| < 0.05, |#Delta #phi| < 0.05", nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
2129  fhLam0PtTMTRD->SetYTitle("#lambda_{0}^{2}");
2130  fhLam0PtTMTRD->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2131  outputContainer->Add(fhLam0PtTMTRD);
2132 
2133  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);
2134  fhLam1ETMTRD->SetYTitle("#lambda_{1}^{2}");
2135  fhLam1ETMTRD->SetXTitle("#it{E} (GeV)");
2136  outputContainer->Add(fhLam1ETMTRD);
2137 
2138  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);
2139  fhDispETMTRD->SetYTitle("Dispersion^{2}");
2140  fhDispETMTRD->SetXTitle("#it{E} (GeV) ");
2141  outputContainer->Add(fhDispETMTRD);
2142  }
2143  }
2144 
2146  {
2147  fhNCellsLam0LowE = new TH2F ("hNCellsLam0LowE","N_{cells} in cluster vs #lambda_{0}^{2}, E < 2 GeV", nbins,nmin, nmax, ssbins,ssmin,ssmax);
2148  fhNCellsLam0LowE->SetXTitle("N_{cells}");
2149  fhNCellsLam0LowE->SetYTitle("#lambda_{0}^{2}");
2150  outputContainer->Add(fhNCellsLam0LowE);
2151 
2152  fhNCellsLam0HighE = new TH2F ("hNCellsLam0HighE","N_{cells} in cluster vs #lambda_{0}^{2}, #it{E} > 2 GeV", nbins,nmin, nmax, ssbins,ssmin,ssmax);
2153  fhNCellsLam0HighE->SetXTitle("N_{cells}");
2154  fhNCellsLam0HighE->SetYTitle("#lambda_{0}^{2}");
2155  outputContainer->Add(fhNCellsLam0HighE);
2156 
2157  fhNCellsLam1LowE = new TH2F ("hNCellsLam1LowE","N_{cells} in cluster vs #lambda_{1}^{2}, E < 2 GeV", nbins,nmin, nmax, ssbins,ssmin,ssmax);
2158  fhNCellsLam1LowE->SetXTitle("N_{cells}");
2159  fhNCellsLam1LowE->SetYTitle("#lambda_{0}^{2}");
2160  outputContainer->Add(fhNCellsLam1LowE);
2161 
2162  fhNCellsLam1HighE = new TH2F ("hNCellsLam1HighE","N_{cells} in cluster vs #lambda_{1}^{2}, #it{E} > 2 GeV", nbins,nmin, nmax, ssbins,ssmin,ssmax);
2163  fhNCellsLam1HighE->SetXTitle("N_{cells}");
2164  fhNCellsLam1HighE->SetYTitle("#lambda_{0}^{2}");
2165  outputContainer->Add(fhNCellsLam1HighE);
2166 
2167  fhNCellsDispLowE = new TH2F ("hNCellsDispLowE","N_{cells} in cluster vs dispersion^{2}, E < 2 GeV", nbins,nmin, nmax, ssbins,ssmin,ssmax);
2168  fhNCellsDispLowE->SetXTitle("N_{cells}");
2169  fhNCellsDispLowE->SetYTitle("D^{2}");
2170  outputContainer->Add(fhNCellsDispLowE);
2171 
2172  fhNCellsDispHighE = new TH2F ("hNCellsDispHighE","N_{cells} in cluster vs dispersion^{2}, E < 2 GeV", nbins,nmin, nmax, ssbins,ssmin,ssmax);
2173  fhNCellsDispHighE->SetXTitle("N_{cells}");
2174  fhNCellsDispHighE->SetYTitle("D^{2}");
2175  outputContainer->Add(fhNCellsDispHighE);
2176 
2177  fhEtaLam0LowE = new TH2F ("hEtaLam0LowE","#eta vs #lambda_{0}^{2}, E < 2 GeV", netabins,etamin,etamax, ssbins,ssmin,ssmax);
2178  fhEtaLam0LowE->SetYTitle("#lambda_{0}^{2}");
2179  fhEtaLam0LowE->SetXTitle("#eta");
2180  outputContainer->Add(fhEtaLam0LowE);
2181 
2182  fhPhiLam0LowE = new TH2F ("hPhiLam0LowE","#phi vs #lambda_{0}^{2}, E < 2 GeV", nphibins,phimin,phimax, ssbins,ssmin,ssmax);
2183  fhPhiLam0LowE->SetYTitle("#lambda_{0}^{2}");
2184  fhPhiLam0LowE->SetXTitle("#phi");
2185  outputContainer->Add(fhPhiLam0LowE);
2186 
2187  fhEtaLam0HighE = new TH2F ("hEtaLam0HighE","#eta vs #lambda_{0}^{2}, #it{E} > 2 GeV", netabins,etamin,etamax, ssbins,ssmin,ssmax);
2188  fhEtaLam0HighE->SetYTitle("#lambda_{0}^{2}");
2189  fhEtaLam0HighE->SetXTitle("#eta");
2190  outputContainer->Add(fhEtaLam0HighE);
2191 
2192  fhPhiLam0HighE = new TH2F ("hPhiLam0HighE","#phi vs #lambda_{0}^{2}, #it{E} > 2 GeV", nphibins,phimin,phimax, ssbins,ssmin,ssmax);
2193  fhPhiLam0HighE->SetYTitle("#lambda_{0}^{2}");
2194  fhPhiLam0HighE->SetXTitle("#phi");
2195  outputContainer->Add(fhPhiLam0HighE);
2196 
2197  fhLam1Lam0LowE = new TH2F ("hLam1Lam0LowE","#lambda_{0}^{2} vs #lambda_{1}^{2} in cluster of E < 2 GeV", ssbins,ssmin,ssmax, ssbins,ssmin,ssmax);
2198  fhLam1Lam0LowE->SetYTitle("#lambda_{0}^{2}");
2199  fhLam1Lam0LowE->SetXTitle("#lambda_{1}^{2}");
2200  outputContainer->Add(fhLam1Lam0LowE);
2201 
2202  fhLam1Lam0HighE = new TH2F ("hLam1Lam0HighE","#lambda_{0}^{2} vs #lambda_{1}^{2} in cluster of #it{E} > 2 GeV", ssbins,ssmin,ssmax, ssbins,ssmin,ssmax);
2203  fhLam1Lam0HighE->SetYTitle("#lambda_{0}^{2}");
2204  fhLam1Lam0HighE->SetXTitle("#lambda_{1}^{2}");
2205  outputContainer->Add(fhLam1Lam0HighE);
2206 
2207  fhLam0DispLowE = new TH2F ("hLam0DispLowE","#lambda_{0}^{2} vs dispersion^{2} in cluster of E < 2 GeV", ssbins,ssmin,ssmax, ssbins,ssmin,ssmax);
2208  fhLam0DispLowE->SetXTitle("#lambda_{0}^{2}");
2209  fhLam0DispLowE->SetYTitle("D^{2}");
2210  outputContainer->Add(fhLam0DispLowE);
2211 
2212  fhLam0DispHighE = new TH2F ("hLam0DispHighE","#lambda_{0}^{2} vs dispersion^{2} in cluster of #it{E} > 2 GeV", ssbins,ssmin,ssmax, ssbins,ssmin,ssmax);
2213  fhLam0DispHighE->SetXTitle("#lambda_{0}^{2}");
2214  fhLam0DispHighE->SetYTitle("D^{2}");
2215  outputContainer->Add(fhLam0DispHighE);
2216 
2217  fhDispLam1LowE = new TH2F ("hDispLam1LowE","Dispersion^{2} vs #lambda_{1}^{2} in cluster of E < 2 GeV", ssbins,ssmin,ssmax, ssbins,ssmin,ssmax);
2218  fhDispLam1LowE->SetXTitle("D^{2}");
2219  fhDispLam1LowE->SetYTitle("#lambda_{1}^{2}");
2220  outputContainer->Add(fhDispLam1LowE);
2221 
2222  fhDispLam1HighE = new TH2F ("hDispLam1HighE","Dispersion^{2} vs #lambda_{1^{2}} in cluster of #it{E} > 2 GeV", ssbins,ssmin,ssmax, ssbins,ssmin,ssmax);
2223  fhDispLam1HighE->SetXTitle("D^{2}");
2224  fhDispLam1HighE->SetYTitle("#lambda_{1}^{2}");
2225  outputContainer->Add(fhDispLam1HighE);
2226 
2227  if(GetCalorimeter() == kEMCAL)
2228  {
2229  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);
2230  fhDispEtaE->SetXTitle("#it{E} (GeV)");
2231  fhDispEtaE->SetYTitle("#sigma^{2}_{#eta #eta}");
2232  outputContainer->Add(fhDispEtaE);
2233 
2234  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);
2235  fhDispPhiE->SetXTitle("#it{E} (GeV)");
2236  fhDispPhiE->SetYTitle("#sigma^{2}_{#phi #phi}");
2237  outputContainer->Add(fhDispPhiE);
2238 
2239  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);
2240  fhSumEtaE->SetXTitle("#it{E} (GeV)");
2241  fhSumEtaE->SetYTitle("#delta^{2}_{#eta #eta}");
2242  outputContainer->Add(fhSumEtaE);
2243 
2244  fhSumPhiE = new TH2F ("hSumPhiE","#delta^{2}_{#phi #phi} = #Sigma w_{i}(#phi_{i})^{2}/ #Sigma w_{i} - <#phi>^{2} vs E",
2245  nptbins,ptmin,ptmax, ssbins,ssmin,ssmax);
2246  fhSumPhiE->SetXTitle("#it{E} (GeV)");
2247  fhSumPhiE->SetYTitle("#delta^{2}_{#phi #phi}");
2248  outputContainer->Add(fhSumPhiE);
2249 
2250  fhSumEtaPhiE = new TH2F ("hSumEtaPhiE","#delta^{2}_{#eta #phi} = #Sigma w_{i}(#phi_{i} #eta_{i} ) / #Sigma w_{i} - <#phi><#eta> vs E",
2251  nptbins,ptmin,ptmax, 2*ssbins,-ssmax,ssmax);
2252  fhSumEtaPhiE->SetXTitle("#it{E} (GeV)");
2253  fhSumEtaPhiE->SetYTitle("#delta^{2}_{#eta #phi}");
2254  outputContainer->Add(fhSumEtaPhiE);
2255 
2256  fhDispEtaPhiDiffE = new TH2F ("hDispEtaPhiDiffE","#sigma^{2}_{#phi #phi} - #sigma^{2}_{#eta #eta} vs E",
2257  nptbins,ptmin,ptmax,200, -10,10);
2258  fhDispEtaPhiDiffE->SetXTitle("#it{E} (GeV)");
2259  fhDispEtaPhiDiffE->SetYTitle("#sigma^{2}_{#phi #phi}-#sigma^{2}_{#eta #eta}");
2260  outputContainer->Add(fhDispEtaPhiDiffE);
2261 
2262  fhSphericityE = new TH2F ("hSphericityE","(#sigma^{2}_{#phi #phi} - #sigma^{2}_{#eta #eta}) / (#sigma^{2}_{#eta #eta} + #sigma^{2}_{#phi #phi}) vs E",
2263  nptbins,ptmin,ptmax, 200, -1,1);
2264  fhSphericityE->SetXTitle("#it{E} (GeV)");
2265  fhSphericityE->SetYTitle("s = (#sigma^{2}_{#phi #phi} - #sigma^{2}_{#eta #eta}) / (#sigma^{2}_{#eta #eta} + #sigma^{2}_{#phi #phi})");
2266  outputContainer->Add(fhSphericityE);
2267 
2268  fhDispSumEtaDiffE = new TH2F ("hDispSumEtaDiffE","#sigma^{2}_{#eta #eta} - #delta^{2}_{#eta #eta} / average vs E", nptbins,ptmin,ptmax, 200,-0.01,0.01);
2269  fhDispSumEtaDiffE->SetXTitle("#it{E} (GeV)");
2270  fhDispSumEtaDiffE->SetYTitle("#sigma^{2}_{#eta #eta} - #delta^{2}_{#eta #eta} / average");
2271  outputContainer->Add(fhDispSumEtaDiffE);
2272 
2273  fhDispSumPhiDiffE = new TH2F ("hDispSumPhiDiffE","#sigma^{2}_{#phi #phi} - #delta^{2}_{#phi #phi} / average vs E", nptbins,ptmin,ptmax, 200,-0.01,0.01);
2274  fhDispSumPhiDiffE->SetXTitle("#it{E} (GeV)");
2275  fhDispSumPhiDiffE->SetYTitle("#sigma^{2}_{#phi #phi} - #delta^{2}_{#phi #phi} / average");
2276  outputContainer->Add(fhDispSumPhiDiffE);
2277 
2278  for(Int_t i = 0; i < 7; i++)
2279  {
2280  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]),
2281  ssbins,ssmin,ssmax , ssbins,ssmin,ssmax);
2282  fhDispEtaDispPhi[i]->SetXTitle("#sigma^{2}_{#eta #eta}");
2283  fhDispEtaDispPhi[i]->SetYTitle("#sigma^{2}_{#phi #phi}");
2284  outputContainer->Add(fhDispEtaDispPhi[i]);
2285 
2286  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]),
2287  ssbins,ssmin,ssmax , ssbins,ssmin,ssmax);
2288  fhLambda0DispEta[i]->SetXTitle("#lambda^{2}_{0}");
2289  fhLambda0DispEta[i]->SetYTitle("#sigma^{2}_{#eta #eta}");
2290  outputContainer->Add(fhLambda0DispEta[i]);
2291 
2292  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]),
2293  ssbins,ssmin,ssmax , ssbins,ssmin,ssmax);
2294  fhLambda0DispPhi[i]->SetXTitle("#lambda^{2}_{0}");
2295  fhLambda0DispPhi[i]->SetYTitle("#sigma^{2}_{#phi #phi}");
2296  outputContainer->Add(fhLambda0DispPhi[i]);
2297  }
2298  }
2299  }
2300  } // Shower shape
2301 
2302  // Track Matching
2303 
2304  if(fFillTMHisto)
2305  {
2306  TString cutTM [] = {"NoCut",""};
2307 
2308  for(Int_t i = 0; i < 2; i++)
2309  {
2310  fhTrackMatchedDEta[i] = new TH2F
2311  (Form("hTrackMatchedDEta%s",cutTM[i].Data()),
2312  Form("d#eta of cluster-track vs cluster energy, %s",cutTM[i].Data()),
2313  nptbins,ptmin,ptmax,nresetabins,resetamin,resetamax);
2314  fhTrackMatchedDEta[i]->SetYTitle("d#eta");
2315  fhTrackMatchedDEta[i]->SetXTitle("#it{E}_{cluster} (GeV)");
2316 
2317  fhTrackMatchedDPhi[i] = new TH2F
2318  (Form("hTrackMatchedDPhi%s",cutTM[i].Data()),
2319  Form("d#phi of cluster-track vs cluster energy, %s",cutTM[i].Data()),
2320  nptbins,ptmin,ptmax,nresphibins,resphimin,resphimax);
2321  fhTrackMatchedDPhi[i]->SetYTitle("d#phi (rad)");
2322  fhTrackMatchedDPhi[i]->SetXTitle("#it{E}_{cluster} (GeV)");
2323 
2324  fhTrackMatchedDEtaDPhi[i] = new TH2F
2325  (Form("hTrackMatchedDEtaDPhi%s",cutTM[i].Data()),
2326  Form("d#eta vs d#phi of cluster-track vs cluster energy, %s",cutTM[i].Data()),
2327  nresetabins,resetamin,resetamax,nresphibins,resphimin,resphimax);
2328  fhTrackMatchedDEtaDPhi[i]->SetYTitle("d#phi (rad)");
2329  fhTrackMatchedDEtaDPhi[i]->SetXTitle("d#eta");
2330 
2331  fhTrackMatchedDEtaPos[i] = new TH2F
2332  (Form("hTrackMatchedDEtaPos%s",cutTM[i].Data()),
2333  Form("d#eta of cluster-track vs cluster energy, %s",cutTM[i].Data()),
2334  nptbins,ptmin,ptmax,nresetabins,resetamin,resetamax);
2335  fhTrackMatchedDEtaPos[i]->SetYTitle("d#eta");
2336  fhTrackMatchedDEtaPos[i]->SetXTitle("#it{E}_{cluster} (GeV)");
2337 
2338  fhTrackMatchedDPhiPos[i] = new TH2F
2339  (Form("hTrackMatchedDPhiPos%s",cutTM[i].Data()),
2340  Form("d#phi of cluster-track vs cluster energy, %s",cutTM[i].Data()),
2341  nptbins,ptmin,ptmax,nresphibins,resphimin,resphimax);
2342  fhTrackMatchedDPhiPos[i]->SetYTitle("d#phi (rad)");
2343  fhTrackMatchedDPhiPos[i]->SetXTitle("#it{E}_{cluster} (GeV)");
2344 
2346  (Form("hTrackMatchedDEtaDPhiPos%s",cutTM[i].Data()),
2347  Form("d#eta vs d#phi of cluster-track vs cluster energy, %s",cutTM[i].Data()),
2348  nresetabins,resetamin,resetamax,nresphibins,resphimin,resphimax);
2349  fhTrackMatchedDEtaDPhiPos[i]->SetYTitle("d#phi (rad)");
2350  fhTrackMatchedDEtaDPhiPos[i]->SetXTitle("d#eta");
2351 
2352  fhTrackMatchedDEtaNeg[i] = new TH2F
2353  (Form("hTrackMatchedDEtaNeg%s",cutTM[i].Data()),
2354  Form("d#eta of cluster-track vs cluster energy, %s",cutTM[i].Data()),
2355  nptbins,ptmin,ptmax,nresetabins,resetamin,resetamax);
2356  fhTrackMatchedDEtaNeg[i]->SetYTitle("d#eta");
2357  fhTrackMatchedDEtaNeg[i]->SetXTitle("#it{E}_{cluster} (GeV)");
2358 
2359  fhTrackMatchedDPhiNeg[i] = new TH2F
2360  (Form("hTrackMatchedDPhiNeg%s",cutTM[i].Data()),
2361  Form("d#phi of cluster-track vs cluster energy, %s",cutTM[i].Data()),
2362  nptbins,ptmin,ptmax,nresphibins,resphimin,resphimax);
2363  fhTrackMatchedDPhiNeg[i]->SetYTitle("d#phi (rad)");
2364  fhTrackMatchedDPhiNeg[i]->SetXTitle("#it{E}_{cluster} (GeV)");
2365 
2367  (Form("hTrackMatchedDEtaDPhiNeg%s",cutTM[i].Data()),
2368  Form("d#eta vs d#phi of cluster-track vs cluster energy, %s",cutTM[i].Data()),
2369  nresetabins,resetamin,resetamax,nresphibins,resphimin,resphimax);
2370  fhTrackMatchedDEtaDPhiNeg[i]->SetYTitle("d#phi (rad)");
2371  fhTrackMatchedDEtaDPhiNeg[i]->SetXTitle("d#eta");
2372 
2373  fhdEdx[i] = new TH2F (Form("hdEdx%s",cutTM[i].Data()),Form("matched track <dE/dx> vs cluster E, %s",cutTM[i].Data()),
2374  nptbins,ptmin,ptmax,ndedxbins, dedxmin, dedxmax);
2375  fhdEdx[i]->SetXTitle("#it{E} (GeV)");
2376  fhdEdx[i]->SetYTitle("<dE/dx>");
2377 
2378  fhEOverP[i] = new TH2F (Form("hEOverP%s",cutTM[i].Data()),Form("matched track E/p vs cluster E, %s",cutTM[i].Data()),
2379  nptbins,ptmin,ptmax,nPoverEbins,pOverEmin,pOverEmax);
2380  fhEOverP[i]->SetXTitle("#it{E} (GeV)");
2381  fhEOverP[i]->SetYTitle("E/p");
2382 
2383  outputContainer->Add(fhTrackMatchedDEta[i]) ;
2384  outputContainer->Add(fhTrackMatchedDPhi[i]) ;
2385  outputContainer->Add(fhTrackMatchedDEtaDPhi[i]) ;
2386  outputContainer->Add(fhTrackMatchedDEtaPos[i]) ;
2387  outputContainer->Add(fhTrackMatchedDPhiPos[i]) ;
2388  outputContainer->Add(fhTrackMatchedDEtaDPhiPos[i]) ;
2389  outputContainer->Add(fhTrackMatchedDEtaNeg[i]) ;
2390  outputContainer->Add(fhTrackMatchedDPhiNeg[i]) ;
2391  outputContainer->Add(fhTrackMatchedDEtaDPhiNeg[i]) ;
2392  outputContainer->Add(fhdEdx[i]);
2393  outputContainer->Add(fhEOverP[i]);
2394 
2396  {
2397  fhTrackMatchedDEtaTRD[i] = new TH2F
2398  (Form("hTrackMatchedDEtaTRD%s",cutTM[i].Data()),
2399  Form("d#eta of cluster-track vs cluster energy, SM behind TRD, %s",cutTM[i].Data()),
2400  nptbins,ptmin,ptmax,nresetabins,resetamin,resetamax);
2401  fhTrackMatchedDEtaTRD[i]->SetYTitle("d#eta");
2402  fhTrackMatchedDEtaTRD[i]->SetXTitle("#it{E}_{cluster} (GeV)");
2403 
2404  fhTrackMatchedDPhiTRD[i] = new TH2F
2405  (Form("hTrackMatchedDPhiTRD%s",cutTM[i].Data()),
2406  Form("d#phi of cluster-track vs cluster energy, SM behing TRD, %s",cutTM[i].Data()),
2407  nptbins,ptmin,ptmax,nresphibins,resphimin,resphimax);
2408  fhTrackMatchedDPhiTRD[i]->SetYTitle("d#phi (rad)");
2409  fhTrackMatchedDPhiTRD[i]->SetXTitle("#it{E}_{cluster} (GeV)");
2410 
2411  fhEOverPTRD[i] = new TH2F
2412  (Form("hEOverPTRD%s",cutTM[i].Data()),
2413  Form("matched track E/p vs cluster E, behind TRD, %s",cutTM[i].Data()),
2414  nptbins,ptmin,ptmax,nPoverEbins,pOverEmin,pOverEmax);
2415  fhEOverPTRD[i]->SetXTitle("#it{E} (GeV)");
2416  fhEOverPTRD[i]->SetYTitle("E/p");
2417 
2418  outputContainer->Add(fhTrackMatchedDEtaTRD[i]) ;
2419  outputContainer->Add(fhTrackMatchedDPhiTRD[i]) ;
2420  outputContainer->Add(fhEOverPTRD[i]);
2421  }
2422 
2423  if(IsDataMC())
2424  {
2426  (Form("hTrackMatchedDEtaMCNoOverlap%s",cutTM[i].Data()),
2427  Form("d#eta of cluster-track vs cluster energy, no other MC particles overlap %s",cutTM[i].Data()),
2428  nptbins,ptmin,ptmax,nresetabins,resetamin,resetamax);
2429  fhTrackMatchedDEtaMCNoOverlap[i]->SetYTitle("d#eta");
2430  fhTrackMatchedDEtaMCNoOverlap[i]->SetXTitle("#it{E}_{cluster} (GeV)");
2431 
2433  (Form("hTrackMatchedDPhiMCNoOverlap%s",cutTM[i].Data()),
2434  Form("d#phi of cluster-track vs cluster energy, no other MC particles overlap %s",cutTM[i].Data()),
2435  nptbins,ptmin,ptmax,nresphibins,resphimin,resphimax);
2436  fhTrackMatchedDPhiMCNoOverlap[i]->SetYTitle("d#phi (rad)");
2437  fhTrackMatchedDPhiMCNoOverlap[i]->SetXTitle("#it{E}_{cluster} (GeV)");
2438 
2439  outputContainer->Add(fhTrackMatchedDEtaMCNoOverlap[i]) ;
2440  outputContainer->Add(fhTrackMatchedDPhiMCNoOverlap[i]) ;
2442  (Form("hTrackMatchedDEtaMCOverlap%s",cutTM[i].Data()),
2443  Form("d#eta of cluster-track vs cluster energy, several MC particles overlap %s",cutTM[i].Data()),
2444  nptbins,ptmin,ptmax,nresetabins,resetamin,resetamax);
2445  fhTrackMatchedDEtaMCOverlap[i]->SetYTitle("d#eta");
2446  fhTrackMatchedDEtaMCOverlap[i]->SetXTitle("#it{E}_{cluster} (GeV)");
2447 
2449  (Form("hTrackMatchedDPhiMCOverlap%s",cutTM[i].Data()),
2450  Form("d#phi of cluster-track vs cluster energy, several MC particles overlap %s",cutTM[i].Data()),
2451  nptbins,ptmin,ptmax,nresphibins,resphimin,resphimax);
2452  fhTrackMatchedDPhiMCOverlap[i]->SetYTitle("d#phi (rad)");
2453  fhTrackMatchedDPhiMCOverlap[i]->SetXTitle("#it{E}_{cluster} (GeV)");
2454 
2455  outputContainer->Add(fhTrackMatchedDEtaMCOverlap[i]) ;
2456  outputContainer->Add(fhTrackMatchedDPhiMCOverlap[i]) ;
2457 
2459  (Form("hTrackMatchedDEtaMCConversion%s",cutTM[i].Data()),
2460  Form("d#eta of cluster-track vs cluster energy, no other MC particles overlap appart from conversions %s",cutTM[i].Data()),
2461  nptbins,ptmin,ptmax,nresetabins,resetamin,resetamax);
2462  fhTrackMatchedDEtaMCConversion[i]->SetYTitle("d#eta");
2463  fhTrackMatchedDEtaMCConversion[i]->SetXTitle("#it{E}_{cluster} (GeV)");
2464 
2466  (Form("hTrackMatchedDPhiMCConversion%s",cutTM[i].Data()),
2467  Form("d#phi of cluster-track vs cluster energy, no other MC particles overlap appart from conversions %s",cutTM[i].Data()),
2468  nptbins,ptmin,ptmax,nresphibins,resphimin,resphimax);
2469  fhTrackMatchedDPhiMCConversion[i]->SetYTitle("d#phi (rad)");
2470  fhTrackMatchedDPhiMCConversion[i]->SetXTitle("#it{E}_{cluster} (GeV)");
2471 
2472  outputContainer->Add(fhTrackMatchedDEtaMCConversion[i]) ;
2473  outputContainer->Add(fhTrackMatchedDPhiMCConversion[i]) ;
2474 
2476  (Form("hTrackMatchedMCParticle%s",cutTM[i].Data()),
2477  Form("Origin of particle vs energy %s",cutTM[i].Data()),
2478  nptbins,ptmin,ptmax,8,0,8);
2479  fhTrackMatchedMCParticle[i]->SetXTitle("#it{E} (GeV)");
2480  //fhTrackMatchedMCParticle[i]->SetYTitle("Particle type");
2481 
2482  fhTrackMatchedMCParticle[i]->GetYaxis()->SetBinLabel(1 ,"Photon");
2483  fhTrackMatchedMCParticle[i]->GetYaxis()->SetBinLabel(2 ,"Electron");
2484  fhTrackMatchedMCParticle[i]->GetYaxis()->SetBinLabel(3 ,"Meson Merged");
2485  fhTrackMatchedMCParticle[i]->GetYaxis()->SetBinLabel(4 ,"Rest");
2486  fhTrackMatchedMCParticle[i]->GetYaxis()->SetBinLabel(5 ,"Conv. Photon");
2487  fhTrackMatchedMCParticle[i]->GetYaxis()->SetBinLabel(6 ,"Conv. Electron");
2488  fhTrackMatchedMCParticle[i]->GetYaxis()->SetBinLabel(7 ,"Conv. Merged");
2489  fhTrackMatchedMCParticle[i]->GetYaxis()->SetBinLabel(8 ,"Conv. Rest");
2490 
2491  outputContainer->Add(fhTrackMatchedMCParticle[i]);
2492  }
2493  }
2494  }
2495 
2496  if(IsPileUpAnalysisOn())
2497  {
2498  TString pileUpName[] = {"SPD","EMCAL","SPDOrEMCAL","SPDAndEMCAL","SPDAndNotEMCAL","EMCALAndNotSPD","NotSPDAndNotEMCAL"} ;
2499 
2500  for(Int_t i = 0 ; i < 7 ; i++)
2501  {
2502  fhPtPhotonPileUp[i] = new TH1F(Form("hPtPhotonPileUp%s",pileUpName[i].Data()),
2503  Form("Selected photon #it{p}_{T} distribution, %s Pile-Up event",pileUpName[i].Data()), nptbins,ptmin,ptmax);
2504  fhPtPhotonPileUp[i]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2505  outputContainer->Add(fhPtPhotonPileUp[i]);
2506 
2507  fhClusterTimeDiffPhotonPileUp[i] = new TH2F(Form("hClusterTimeDiffPhotonPileUp%s",pileUpName[i].Data()),
2508  Form("Photon cluster E vs #it{t}_{max}-#it{t}_{cell} in cluster, %s Pile-Up event",pileUpName[i].Data()),
2509  nptbins,ptmin,ptmax,400,-200,200);
2510  fhClusterTimeDiffPhotonPileUp[i]->SetXTitle("#it{E} (GeV)");
2511  fhClusterTimeDiffPhotonPileUp[i]->SetYTitle("#it{t}_{max}-#it{t}_{cell} (ns)");
2512  outputContainer->Add(fhClusterTimeDiffPhotonPileUp[i]);
2513  }
2514 
2515  fhTimePtPhotonNoCut = new TH2F ("hTimePtPhoton_NoCut","time of photon cluster vs pT of clusters, no cut", nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
2516  fhTimePtPhotonNoCut->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2517  fhTimePtPhotonNoCut->SetYTitle("#it{time} (ns)");
2518  outputContainer->Add(fhTimePtPhotonNoCut);
2519 
2520  fhTimePtPhotonSPD = new TH2F ("hTimePtPhoton_SPD","time of photon cluster vs pT of clusters, SPD cut", nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
2521  fhTimePtPhotonSPD->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2522  fhTimePtPhotonSPD->SetYTitle("#it{time} (ns)");
2523  outputContainer->Add(fhTimePtPhotonSPD);
2524 
2525  fhTimeNPileUpVertSPD = new TH2F ("hTime_NPileUpVertSPD","time of cluster vs N pile-up SPD vertex", ntimebins,timemin,timemax,20,0,20);
2526  fhTimeNPileUpVertSPD->SetYTitle("# vertex ");
2527  fhTimeNPileUpVertSPD->SetXTitle("#it{time} (ns)");
2528  outputContainer->Add(fhTimeNPileUpVertSPD);
2529 
2530  fhTimeNPileUpVertTrack = new TH2F ("hTime_NPileUpVertTracks","time of cluster vs N pile-up Tracks vertex", ntimebins,timemin,timemax, 20,0,20 );
2531  fhTimeNPileUpVertTrack->SetYTitle("# vertex ");
2532  fhTimeNPileUpVertTrack->SetXTitle("#it{time} (ns)");
2533  outputContainer->Add(fhTimeNPileUpVertTrack);
2534 
2535  fhPtPhotonNPileUpSPDVtx = new TH2F ("hPtPhoton_NPileUpVertSPD","pT of cluster vs N pile-up SPD vertex",
2536  nptbins,ptmin,ptmax,20,0,20);
2537  fhPtPhotonNPileUpSPDVtx->SetYTitle("# vertex ");
2538  fhPtPhotonNPileUpSPDVtx->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2539  outputContainer->Add(fhPtPhotonNPileUpSPDVtx);
2540 
2541  fhPtPhotonNPileUpTrkVtx = new TH2F ("hPtPhoton_NPileUpVertTracks","pT of cluster vs N pile-up Tracks vertex",
2542  nptbins,ptmin,ptmax, 20,0,20 );
2543  fhPtPhotonNPileUpTrkVtx->SetYTitle("# vertex ");
2544  fhPtPhotonNPileUpTrkVtx->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2545  outputContainer->Add(fhPtPhotonNPileUpTrkVtx);
2546 
2547  fhPtPhotonNPileUpSPDVtxTimeCut = new TH2F ("hPtPhoton_NPileUpVertSPD_TimeCut","pT of cluster vs N pile-up SPD vertex, |tof| < 25 ns",
2548  nptbins,ptmin,ptmax,20,0,20);
2549  fhPtPhotonNPileUpSPDVtxTimeCut->SetYTitle("# vertex ");
2550  fhPtPhotonNPileUpSPDVtxTimeCut->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2551  outputContainer->Add(fhPtPhotonNPileUpSPDVtxTimeCut);
2552 
2553  fhPtPhotonNPileUpTrkVtxTimeCut = new TH2F ("hPtPhoton_NPileUpVertTracks_TimeCut","pT of cluster vs N pile-up Tracks vertex, |tof| < 25 ns",
2554  nptbins,ptmin,ptmax, 20,0,20 );
2555  fhPtPhotonNPileUpTrkVtxTimeCut->SetYTitle("# vertex ");
2556  fhPtPhotonNPileUpTrkVtxTimeCut->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2557  outputContainer->Add(fhPtPhotonNPileUpTrkVtxTimeCut);
2558 
2559  fhPtPhotonNPileUpSPDVtxTimeCut2 = new TH2F ("hPtPhoton_NPileUpVertSPD_TimeCut2","pT of cluster vs N pile-up SPD vertex, -25 < tof < 75 ns",
2560  nptbins,ptmin,ptmax,20,0,20);
2561  fhPtPhotonNPileUpSPDVtxTimeCut2->SetYTitle("# vertex ");
2562  fhPtPhotonNPileUpSPDVtxTimeCut2->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2563  outputContainer->Add(fhPtPhotonNPileUpSPDVtxTimeCut2);
2564 
2565  fhPtPhotonNPileUpTrkVtxTimeCut2 = new TH2F ("hPtPhoton_NPileUpVertTracks_TimeCut2","pT of cluster vs N pile-up Tracks vertex, -25 < tof < 75 ns",
2566  nptbins,ptmin,ptmax, 20,0,20 );
2567  fhPtPhotonNPileUpTrkVtxTimeCut2->SetYTitle("# vertex ");
2568  fhPtPhotonNPileUpTrkVtxTimeCut2->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2569  outputContainer->Add(fhPtPhotonNPileUpTrkVtxTimeCut2);
2570 
2571  }
2572 
2574  {
2575  for(Int_t ieta = 0; ieta < 4; ieta++)
2576  {
2577  for(Int_t iphi = 0; iphi < 3; iphi++)
2578  {
2579 // fhLam0EMCALRegion[ieta][iphi] =
2580 // new TH2F(Form("hLam0_eta%d_phi%d",ieta,iphi),
2581 // Form("cluster from converted photon, #it{p}_{T} vs #lambda_{0}^{2}, region eta %d, phi %d",ieta,iphi),
2582 // nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
2583 // fhLam0EMCALRegion[ieta][iphi]->SetYTitle("#lambda_{0}^{2}");
2584 // fhLam0EMCALRegion[ieta][iphi]->SetXTitle("#it{p}_{T} (GeV)");
2585 // outputContainer->Add(fhLam0EMCALRegion[ieta][iphi]) ;
2586 //
2587 // if(GetFirstSMCoveredByTRD() >= 0)
2588 // {
2589 // fhLam0EMCALRegionTRD[ieta][iphi] =
2590 // new TH2F(Form("hLam0TRD_eta%d_phi%d",ieta,iphi),
2591 // Form("cluster from converted photon, #it{p}_{T} vs #lambda_{0}^{2}, region eta %d, phi %d, SM covered by TRD",ieta,iphi),
2592 // nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
2593 // fhLam0EMCALRegionTRD[ieta][iphi]->SetYTitle("#lambda_{0}^{2}");
2594 // fhLam0EMCALRegionTRD[ieta][iphi]->SetXTitle("#it{p}_{T} (GeV)");
2595 // outputContainer->Add(fhLam0EMCALRegionTRD[ieta][iphi]) ;
2596 // } // TRD
2597 
2598  for(Int_t ism = 0; ism < GetCaloUtils()->GetNumberOfSuperModulesUsed(); ism++)
2599  {
2600  fhLam0EMCALRegionPerSM[ieta][iphi][ism] =
2601  new TH2F(Form("hLam0_eta%d_phi%d_sm%d",ieta,iphi,ism),
2602  Form("cluster from converted photon, #it{p}_{T} vs #lambda_{0}^{2}, sm %d, region eta %d, phi %d",ism,ieta,iphi),
2603  nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
2604  fhLam0EMCALRegionPerSM[ieta][iphi][ism]->SetYTitle("#lambda_{0}^{2}");
2605  fhLam0EMCALRegionPerSM[ieta][iphi][ism]->SetXTitle("#it{p}_{T} (GeV)");
2606  outputContainer->Add(fhLam0EMCALRegionPerSM[ieta][iphi][ism]) ;
2607 
2608  fhLam1EMCALRegionPerSM[ieta][iphi][ism] =
2609  new TH2F(Form("hLam1_eta%d_phi%d_sm%d",ieta,iphi,ism),
2610  Form("cluster from converted photon, #it{p}_{T} vs #lambda_{1}^{2}, sm %d, region eta %d, phi %d",ism,ieta,iphi),
2611  nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
2612  fhLam1EMCALRegionPerSM[ieta][iphi][ism]->SetYTitle("#lambda_{1}^{2}");
2613  fhLam1EMCALRegionPerSM[ieta][iphi][ism]->SetXTitle("#it{p}_{T} (GeV)");
2614  outputContainer->Add(fhLam1EMCALRegionPerSM[ieta][iphi][ism]) ;
2615  }
2616  } // iphi
2617  } // ieta
2618 
2619  Float_t ptLimit[] = {2,3,4,5,6,8,10,12};
2620  TString l0bin [] = {"0.23<#lambda^{2}_{0}<0.26","0.3<#lambda^{2}_{0}<0.4"};
2621 
2622  for(Int_t il0 = 0; il0 < 2; il0++)
2623  {
2624  for(Int_t ipt = 0; ipt < 7; ipt++)
2625  {
2626  fhEtaPhiLam0BinPtBin[il0][ipt] = new TH2F
2627  (Form("hEtaPhiLam0Bin%d_PtBin%d",il0,ipt),
2628  Form("#eta vs #phi in #it{p}_{T}=[%2.1f,%2.1f] GeV/#it{c}, %s",
2629  ptLimit[ipt],ptLimit[ipt+1],l0bin[il0].Data()),
2630  netabins,etamin,etamax,nphibins,phimin,phimax);
2631  fhEtaPhiLam0BinPtBin[il0][ipt]->SetYTitle("#phi (rad)");
2632  fhEtaPhiLam0BinPtBin[il0][ipt]->SetXTitle("#eta");
2633  outputContainer->Add(fhEtaPhiLam0BinPtBin[il0][ipt]) ;
2634 
2635  fhColRowLam0BinPtBin[il0][ipt] = new TH2F
2636  (Form("hColRowLam0Bin%d_PtBin%d",il0,ipt),
2637  Form("row vs column in #it{p}_{T}=[%2.1f,%2.1f] GeV/#it{c}, %s, w > 0",
2638  ptLimit[ipt],ptLimit[ipt+1],l0bin[il0].Data()),
2639  96+2,-1.5,96+0.5,(8*24+2*8)+2,-1.5,(8*24+2*8)+0.5); // fix to generalize to other periods
2640  fhColRowLam0BinPtBin[il0][ipt]->SetYTitle("row");
2641  fhColRowLam0BinPtBin[il0][ipt]->SetXTitle("column");
2642  outputContainer->Add(fhColRowLam0BinPtBin[il0][ipt]) ;
2643 
2644  fhColRowLam0BinPtBinWeighted[il0][ipt] = new TH2F
2645  (Form("hColRowLam0Bin%d_PtBin%dWeighted",il0,ipt),
2646  Form("cluster cell row vs column weighted in #it{p}_{T}=[%2.1f,%2.1f] GeV/#it{c}, %s",
2647  ptLimit[ipt],ptLimit[ipt+1],l0bin[il0].Data()),
2648  96+2,-1.5,96+0.5,(8*24+2*8)+2,-1.5,(8*24+2*8)+0.5); // fix to generalize to other periods
2649  fhColRowLam0BinPtBinWeighted[il0][ipt]->SetYTitle("row");
2650  fhColRowLam0BinPtBinWeighted[il0][ipt]->SetXTitle("column");
2651  outputContainer->Add(fhColRowLam0BinPtBinWeighted[il0][ipt]) ;
2652 
2653  fhColRowLam0BinPtBinLargeTime[il0][ipt] = new TH2F
2654  (Form("hColRowLam0Bin%d_PtBin%d_LargeTimeInClusterCell",il0,ipt),
2655  Form("row vs column in #it{p}_{T}=[%2.1f,%2.1f] GeV/#it{c}, %s, |t| > 50 ns, w > 0",
2656  ptLimit[ipt],ptLimit[ipt+1],l0bin[il0].Data()),
2657  96+2,-1.5,96+0.5,(8*24+2*8)+2,-1.5,(8*24+2*8)+0.5); // fix to generalize to other periods
2658  fhColRowLam0BinPtBinLargeTime[il0][ipt]->SetYTitle("row");
2659  fhColRowLam0BinPtBinLargeTime[il0][ipt]->SetXTitle("column");
2660  outputContainer->Add(fhColRowLam0BinPtBinLargeTime[il0][ipt]) ;
2661 
2662 // fhEtaPhiLam0BinPtBinSMShared[il0][ipt] = new TH2F
2663 // (Form("hEtaPhiLam0Bin%d_PtBin%d_SMShared",il0,ipt),
2664 // Form("#eta vs #phi in #it{p}_{T}=[%2.1f,%2.1f] GeV/#it{c}, %s",
2665 // ptLimit[ipt],ptLimit[ipt+1],l0bin[il0].Data()),
2666 // netabins,etamin,etamax,nphibins,phimin,phimax);
2667 // fhEtaPhiLam0BinPtBinSMShared[il0][ipt]->SetYTitle("#phi (rad)");
2668 // fhEtaPhiLam0BinPtBinSMShared[il0][ipt]->SetXTitle("#eta");
2669 // outputContainer->Add(fhEtaPhiLam0BinPtBinSMShared[il0][ipt]) ;
2670 //
2671 // fhColRowLam0BinPtBinSMShared[il0][ipt] = new TH2F
2672 // (Form("hColRowLam0Bin%d_PtBin%d_SMShared",il0,ipt),
2673 // Form("row vs column in #it{p}_{T}=[%2.1f,%2.1f] GeV/#it{c}, %s, w > 0",
2674 // ptLimit[ipt],ptLimit[ipt+1],l0bin[il0].Data()),
2675 // 96,0,96,5*24,0,5*24); // fix to generalize to other periods
2676 // fhColRowLam0BinPtBinSMShared[il0][ipt]->SetYTitle("row");
2677 // fhColRowLam0BinPtBinSMShared[il0][ipt]->SetXTitle("column");
2678 // outputContainer->Add(fhColRowLam0BinPtBinSMShared[il0][ipt]) ;
2679 
2680  fhEtaPhiLargeTimeInClusterCell[il0][ipt] = new TH2F
2681  (Form("hEtaPhiLam0Bin%d_PtBin%d_LargeTimeInClusterCell",il0,ipt),
2682  Form("#eta vs #phi in #it{p}_{T}=[%2.1f,%2.1f] GeV/#it{c}, w > 0, t > 50 ns, %s",
2683  ptLimit[ipt],ptLimit[ipt+1],l0bin[il0].Data()),
2684  netabins,etamin,etamax,nphibins,phimin,phimax);
2685  fhEtaPhiLargeTimeInClusterCell[il0][ipt]->SetYTitle("#phi (rad)");
2686  fhEtaPhiLargeTimeInClusterCell[il0][ipt]->SetXTitle("#eta");
2687  outputContainer->Add(fhEtaPhiLargeTimeInClusterCell[il0][ipt]) ;
2688 
2689  fhCellClusterIndexEAndTime[il0][ipt] = new TH2F
2690  (Form("hCellClusterIndexEAndTimeLam0Bin%d_PtBin%d",il0,ipt),
2691  Form("#it{t}_{cell} vs cell index (E sorted) in #it{p}_{T}=[%2.1f,%2.1f] GeV/#it{c}, w > 0, %s",
2692  ptLimit[ipt],ptLimit[ipt+1],l0bin[il0].Data()),
2693  ntimebins,timemin,timemax,30,0,30);
2694  fhCellClusterIndexEAndTime[il0][ipt] ->SetYTitle("Index (sorted #it{E})");
2695  fhCellClusterIndexEAndTime[il0][ipt] ->SetXTitle("#it{time}_{cell} (ns)");
2696  outputContainer->Add(fhCellClusterIndexEAndTime[il0][ipt] ) ;
2697 
2698  fhCellClusterEAndTime[il0][ipt] = new TH2F
2699  (Form("hCellClusterEAndTimeLam0Bin%d_PtBin%d",il0,ipt),
2700  Form("#it{E}_{cell} vs #it{t}_{cell} in #it{p}_{T}=[%2.1f,%2.1f] GeV/#it{c}, w > 0, %s",
2701  ptLimit[ipt],ptLimit[ipt+1],l0bin[il0].Data()),
2702  ntimebins,timemin,timemax,100,0,5);
2703  fhCellClusterEAndTime[il0][ipt] ->SetYTitle("#it{E}_{cell} (GeV)");
2704  fhCellClusterEAndTime[il0][ipt] ->SetXTitle("#it{time}_{cell} (ns)");
2705  outputContainer->Add(fhCellClusterEAndTime[il0][ipt] ) ;
2706 
2707  fhCellClusterEFracAndTime[il0][ipt] = new TH2F
2708  (Form("hCellClusterEFracAndTimeLam0Bin%d_PtBin%d",il0,ipt),
2709  Form("#it{E}_{cell} vs #it{t}_{cell} in #it{p}_{T}=[%2.1f,%2.1f] GeV/#it{c}, w > 0, %s",
2710  ptLimit[ipt],ptLimit[ipt+1],l0bin[il0].Data()),
2711  ntimebins,timemin,timemax,100,0,1);
2712  fhCellClusterEFracAndTime[il0][ipt] ->SetYTitle("#it{E}_{cell} (GeV)");
2713  fhCellClusterEFracAndTime[il0][ipt] ->SetXTitle("#it{time}_{cell} (ns)");
2714  outputContainer->Add(fhCellClusterEFracAndTime[il0][ipt] ) ;
2715  } // pt bin
2716 
2717  for(Int_t ism = 0; ism < GetCaloUtils()->GetNumberOfSuperModulesUsed(); ism++)
2718  {
2719  fhLam1Lam0BinPerSM[il0][ism] = new TH2F
2720  (Form("hLam1Lam0Bin%d_sm%d",il0,ism),
2721  Form("#it{p}_{T} vs #lambda^{2}_{1} in sm %d, %s",ism,l0bin[il0].Data()),
2722  nptbins,ptmin,ptmax,40,0,0.4);
2723  fhLam1Lam0BinPerSM[il0][ism]->SetYTitle("#lambda^{2}_{1}");
2724  fhLam1Lam0BinPerSM[il0][ism]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2725  outputContainer->Add(fhLam1Lam0BinPerSM[il0][ism]) ;
2726 
2727  fhTimeLam0BinPerSM[il0][ism] = new TH2F
2728  (Form("hTimeLam0Bin%d_sm%d",il0,ism),
2729  Form("#it{p}_{T} vs cluster cell time in sm %d, %s, w > 0",ism,l0bin[il0].Data()),
2730  nptbins,ptmin,ptmax,ntimebins,timemin,timemax);
2731  fhTimeLam0BinPerSM[il0][ism]->SetYTitle("cell time (ns)");
2732  fhTimeLam0BinPerSM[il0][ism]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2733  outputContainer->Add(fhTimeLam0BinPerSM[il0][ism]) ;
2734 
2735  fhTimeLam0BinPerSMWeighted[il0][ism] = new TH2F
2736  (Form("hTimeLam0Bin%d_sm%d_Weighted",il0,ism),
2737  Form("#it{p}_{T} vs cluster cell time weighted in sm %d, %s",ism,l0bin[il0].Data()),
2738  nptbins,ptmin,ptmax,ntimebins,timemin,timemax);
2739  fhTimeLam0BinPerSMWeighted[il0][ism]->SetYTitle("cell time (ns)");
2740  fhTimeLam0BinPerSMWeighted[il0][ism]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2741  outputContainer->Add(fhTimeLam0BinPerSMWeighted[il0][ism]) ;
2742 
2743  fhDTimeLam0BinPerSM[il0][ism] = new TH2F
2744  (Form("hDTimeLam0Bin%d_sm%d",il0,ism),
2745  Form("#it{p}_{T} vs t_{cluster}-t_{cell} in sm %d, %s, w > 0",ism,l0bin[il0].Data()),
2746  nptbins,ptmin,ptmax,ntimebins,timemin,timemax);
2747  fhDTimeLam0BinPerSM[il0][ism]->SetYTitle("t_{cluster}-t_{cell} (ns)");
2748  fhDTimeLam0BinPerSM[il0][ism]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2749  outputContainer->Add(fhDTimeLam0BinPerSM[il0][ism]) ;
2750 
2751  fhDTimeLam0BinPerSMWeighted[il0][ism] = new TH2F
2752  (Form("hDTimeLam0Bin%d_sm%d_Weighted",il0,ism),
2753  Form("#it{p}_{T} vs t_{cluster}-t_{cell} weighted in sm %d, %s",ism,l0bin[il0].Data()),
2754  nptbins,ptmin,ptmax,ntimebins,timemin,timemax);
2755  fhDTimeLam0BinPerSMWeighted[il0][ism]->SetYTitle("t_{cluster}-t_{cell} (ns)");
2756  fhDTimeLam0BinPerSMWeighted[il0][ism]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2757  outputContainer->Add(fhDTimeLam0BinPerSMWeighted[il0][ism]) ;
2758 
2759  fhCellClusterEFracLam0BinPerSM[il0][ism] = new TH2F
2760  (Form("hCellClusterEFracLam0Bin%d_sm%d",il0,ism),
2761  Form("#it{p}_{T} vs cell E / cluster E in sm %d, %s, w > 0",ism,l0bin[il0].Data()),
2762  nptbins,ptmin,ptmax,100,0,1);
2763  fhCellClusterEFracLam0BinPerSM[il0][ism]->SetYTitle("cell E / cluster E ");
2764  fhCellClusterEFracLam0BinPerSM[il0][ism]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2765  outputContainer->Add(fhCellClusterEFracLam0BinPerSM[il0][ism]) ;
2766 
2767 // fhCellClusterEFracLam0BinPerSMWeighted[il0][ism] = new TH2F
2768 // (Form("hCellClusterEFracLam0Bin%d_sm%d_Weighted",il0,ism),
2769 // Form("#it{p}_{T} vs cell E / cluster E weighted in sm %d, %s",ism,l0bin[il0].Data()),
2770 // nptbins,ptmin,ptmax,100,0,1);
2771 // fhCellClusterEFracLam0BinPerSMWeighted[il0][ism]->SetYTitle("cell E / cluster E ");
2772 // fhCellClusterEFracLam0BinPerSMWeighted[il0][ism]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2773 // outputContainer->Add(fhCellClusterEFracLam0BinPerSMWeighted[il0][ism]) ;
2774 
2776  (Form("hCellClusterEFracLam0Bin%d_sm%d_LargeTimeInClusterCell",il0,ism),
2777  Form("#it{p}_{T} vs cell E / cluster E in sm %d, %s, w > 0, |t_{cell}| > 50 ns",ism,l0bin[il0].Data()),
2778  nptbins,ptmin,ptmax,100,0,1);
2779  fhCellClusterEFracLam0BinPerSMLargeTime[il0][ism]->SetYTitle("cell E / cluster E ");
2780  fhCellClusterEFracLam0BinPerSMLargeTime[il0][ism]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2781  outputContainer->Add(fhCellClusterEFracLam0BinPerSMLargeTime[il0][ism]) ;
2782 
2784  (Form("hCellClusterEFracLam0Bin%d_sm%d_LargeTimeInClusterCell_Total",il0,ism),
2785  Form("#it{p}_{T} vs cell E / cluster E in sm %d, %s, w > 0, all |t_{cell}| > 50 ns",ism,l0bin[il0].Data()),
2786  nptbins,ptmin,ptmax,100,0,1);
2787  fhCellClusterEFracLam0BinPerSMLargeTimeTotal[il0][ism]->SetYTitle("cell E / cluster E ");
2788  fhCellClusterEFracLam0BinPerSMLargeTimeTotal[il0][ism]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2789  outputContainer->Add(fhCellClusterEFracLam0BinPerSMLargeTimeTotal[il0][ism]) ;
2790 
2791  fhCellClusterELam0BinPerSM[il0][ism] = new TH2F
2792  (Form("hCellClusterELam0Bin%d_sm%d",il0,ism),
2793  Form("#it{p}_{T} vs cell E in sm %d, %s, w > 0",ism,l0bin[il0].Data()),
2794  nptbins,ptmin,ptmax,500,0,10);
2795  fhCellClusterELam0BinPerSM[il0][ism]->SetYTitle("cell E (GeV)");
2796  fhCellClusterELam0BinPerSM[il0][ism]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2797  outputContainer->Add(fhCellClusterELam0BinPerSM[il0][ism]) ;
2798 
2800  (Form("hCellClusterELam0Bin%d_sm%d_Weighted",il0,ism),
2801  Form("#it{p}_{T} vs cell E weighted in sm %d, %s",ism,l0bin[il0].Data()),
2802  nptbins,ptmin,ptmax,500,0,10);
2803  fhCellClusterELam0BinPerSMWeighted[il0][ism]->SetYTitle("cell E (GeV)");
2804  fhCellClusterELam0BinPerSMWeighted[il0][ism]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2805  outputContainer->Add(fhCellClusterELam0BinPerSMWeighted[il0][ism]) ;
2806 
2808  (Form("hCellClusterELam0Bin%d_sm%d_LargeTimeInClusterCell",il0,ism),
2809  Form("#it{p}_{T} vs cell E in sm %d, %s, w > 0, |t_{cell}| > 50 ns",ism,l0bin[il0].Data()),
2810  nptbins,ptmin,ptmax,500,0,10);
2811  fhCellClusterELam0BinPerSMLargeTime[il0][ism]->SetYTitle("cell E (GeV)");
2812  fhCellClusterELam0BinPerSMLargeTime[il0][ism]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2813  outputContainer->Add(fhCellClusterELam0BinPerSMLargeTime[il0][ism]) ;
2814 
2816  (Form("hCellClusterIndexELam0Bin%d_sm%d_LargeTimeInClusterCell",il0,ism),
2817  Form("#it{p}_{T} vs cell index (E sorted) in sm %d, %s, w > 0, |t_{cell}| > 50 ns",ism,l0bin[il0].Data()),
2818  nptbins,ptmin,ptmax,30,0,30);
2819  fhCellClusterIndexELam0BinPerSMLargeTime[il0][ism]->SetYTitle("Index");
2820  fhCellClusterIndexELam0BinPerSMLargeTime[il0][ism]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2821  outputContainer->Add(fhCellClusterIndexELam0BinPerSMLargeTime[il0][ism]) ;
2822 
2823  fhNCellsWithLargeTimeInCluster[il0][ism] = new TH2F
2824  (Form("hNCellsWithLargeTimeInClusterLam0Bin%d_sm%d",il0,ism),
2825  Form("#it{p}_{T} vs number of cells in sm %d, %s, w > 0, |t_{cell}| > 50 ns",ism,l0bin[il0].Data()),
2826  nptbins,ptmin,ptmax,30,0,30);
2827  fhNCellsWithLargeTimeInCluster[il0][ism]->SetYTitle("#it{n}_{cells}");
2828  fhNCellsWithLargeTimeInCluster[il0][ism]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2829  outputContainer->Add(fhNCellsWithLargeTimeInCluster[il0][ism]) ;
2830 
2831 // if(ism < 12)
2832 // {
2833 // fhTimeLam0BinPerSMShared[il0][ism] = new TH2F
2834 // (Form("hTimeLam0Bin%d_sm%d_SMShared",il0,ism),
2835 // Form("#it{p}_{T} vs cluster cell time in sm %d, %s, w > 0",ism,l0bin[il0].Data()),
2836 // nptbins,ptmin,ptmax,ntimebins,timemin,timemax);
2837 // fhTimeLam0BinPerSMShared[il0][ism]->SetYTitle("cell time (ns)");
2838 // fhTimeLam0BinPerSMShared[il0][ism]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2839 // outputContainer->Add(fhTimeLam0BinPerSMShared[il0][ism]) ;
2840 // } // Run 1 SM
2841  } // sm
2842  } // l0 bin
2843 
2844  for(Int_t ism = 0; ism < GetCaloUtils()->GetNumberOfSuperModulesUsed(); ism++)
2845  {
2846  fhLam0PerSM[ism] = new TH2F
2847  (Form("hLam0_sm%d",ism),
2848  Form("#it{p}_{T} vs #lambda^{2}_{0} in sm %d",ism),
2849  nptbins,ptmin,ptmax,40,0,0.4);
2850  fhLam0PerSM[ism]->SetYTitle("#lambda^{2}_{0}");
2851  fhLam0PerSM[ism]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2852  outputContainer->Add(fhLam0PerSM[ism]) ;
2853 
2854  fhLam1PerSM[ism] = new TH2F
2855  (Form("hLam1_sm%d",ism),
2856  Form("#it{p}_{T} vs #lambda^{2}_{1} in sm %d",ism),
2857  nptbins,ptmin,ptmax,40,0,0.4);
2858  fhLam1PerSM[ism]->SetYTitle("#lambda^{2}_{1}");
2859  fhLam1PerSM[ism]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2860  outputContainer->Add(fhLam1PerSM[ism]) ;
2861 
2863  (Form("hLam0_sm%d_LargeTimeInClusterCell",ism),
2864  Form("#it{p}_{T} vs #lambda^{2}_{0} in sm %d,|t_{secondary cell}| > 50 ns",ism),
2865  nptbins,ptmin,ptmax,40,0,0.4);
2866  fhLam0PerSMLargeTimeInClusterCell[ism]->SetYTitle("#lambda^{2}_{0}");
2867  fhLam0PerSMLargeTimeInClusterCell[ism]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2868  outputContainer->Add(fhLam0PerSMLargeTimeInClusterCell[ism]) ;
2869 
2871  (Form("hLam1_sm%d_LargeTimeInClusterCell",ism),
2872  Form("#it{p}_{T} vs #lambda^{2}_{1} in sm %d, |t_{secondary cell}| > 50 ns",ism),
2873  nptbins,ptmin,ptmax,40,0,0.4);
2874  fhLam1PerSMLargeTimeInClusterCell[ism]->SetYTitle("#lambda^{2}_{1}");
2875  fhLam1PerSMLargeTimeInClusterCell[ism]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2876  outputContainer->Add(fhLam1PerSMLargeTimeInClusterCell[ism]) ;
2877 
2878 // fhLam0PerSMSPDPileUp[ism] = new TH2F
2879 // (Form("hLam0_sm%d_SPDPileUp",ism),
2880 // Form("#it{p}_{T} vs #lambda^{2}_{0} in sm %d, Pile-up event SPD",ism),
2881 // nptbins,ptmin,ptmax,40,0,0.4);
2882 // fhLam0PerSMSPDPileUp[ism]->SetYTitle("#lambda^{2}_{0}");
2883 // fhLam0PerSMSPDPileUp[ism]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2884 // outputContainer->Add(fhLam0PerSMSPDPileUp[ism]) ;
2885 //
2886 // fhLam1PerSMSPDPileUp[ism] = new TH2F
2887 // (Form("hLam1_sm%d_SPDPileUp",ism),
2888 // Form("#it{p}_{T} vs #lambda^{2}_{1} in sm %d, Pile-up event SPD",ism),
2889 // nptbins,ptmin,ptmax,40,0,0.4);
2890 // fhLam1PerSMSPDPileUp[ism]->SetYTitle("#lambda^{2}_{1}");
2891 // fhLam1PerSMSPDPileUp[ism]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2892 // outputContainer->Add(fhLam1PerSMSPDPileUp[ism]) ;
2893 
2894 // if(ism < 12)
2895 // {
2896 // fhLam0PerSMShared[ism] = new TH2F
2897 // (Form("hLam0_sm%d_SMShared",ism),
2898 // Form("#it{p}_{T} vs #lambda^{2}_{0} in sm %d, SM shared",ism),
2899 // nptbins,ptmin,ptmax,40,0,0.4);
2900 // fhLam0PerSMShared[ism]->SetYTitle("#lambda^{2}_{0}");
2901 // fhLam0PerSMShared[ism]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2902 // outputContainer->Add(fhLam0PerSMShared[ism]) ;
2903 //
2904 // fhLam1PerSMShared[ism] = new TH2F
2905 // (Form("hLam1_sm%d_SMShared",ism),
2906 // Form("#it{p}_{T} vs #lambda^{2}_{1} in sm %d, SM shared",ism),
2907 // nptbins,ptmin,ptmax,40,0,0.4);
2908 // fhLam1PerSMShared[ism]->SetYTitle("#lambda^{2}_{1}");
2909 // fhLam1PerSMShared[ism]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2910 // outputContainer->Add(fhLam1PerSMShared[ism]) ;
2911 // } // run1
2912  } // sm
2913 
2914  for(Int_t ilarge = 0; ilarge < 5; ilarge++)
2915  {
2917  (Form("hLam0_NLargeTimeInClusterCell%d",ilarge),
2918  Form("#it{p}_{T} vs #lambda^{2}_{0} in sm %d,|t_{secondary cell}| > 50 ns",ilarge),
2919  nptbins,ptmin,ptmax,40,0,0.4);
2920  fhLam0PerNLargeTimeInClusterCell[ilarge]->SetYTitle("#lambda^{2}_{0}");
2921  fhLam0PerNLargeTimeInClusterCell[ilarge]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2922  outputContainer->Add(fhLam0PerNLargeTimeInClusterCell[ilarge]) ;
2923 
2925  (Form("hLam1_NLargeTimeInClusterCell%d",ilarge),
2926  Form("#it{p}_{T} vs #lambda^{2}_{1} in sm %d, |t_{secondary cell}| > 50 ns",ilarge),
2927  nptbins,ptmin,ptmax,40,0,0.4);
2928  fhLam1PerNLargeTimeInClusterCell[ilarge]->SetYTitle("#lambda^{2}_{1}");
2929  fhLam1PerNLargeTimeInClusterCell[ilarge]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2930  outputContainer->Add(fhLam1PerNLargeTimeInClusterCell[ilarge]) ;
2931  }
2932 
2933  } // regions in EMCal
2934 
2935 
2936  if(IsDataMC())
2937  {
2938  TString ptype[] = { "#gamma" , "#gamma_{#pi decay}" , "#gamma_{#eta decay}", "#gamma_{other decay}",
2939  "#pi^{0}" , "#eta" , "e^{#pm}" , "#gamma->e^{#pm}" ,
2940  "hadron?" , "Anti-N" , "Anti-P" ,
2941  "Neutron" , "Proton" , "#pi^{#pm}" ,
2942  "#gamma_{prompt}", "#gamma_{fragmentation}", "#gamma_{ISR}" , "String" } ;
2943 
2944  TString pname[] = { "Photon" , "PhotonPi0Decay" , "PhotonEtaDecay", "PhotonOtherDecay",
2945  "Pi0" , "Eta" , "Electron" , "Conversion" ,
2946  "Hadron" , "AntiNeutron" , "AntiProton" ,
2947  "Neutron" , "Proton" , "ChPion" ,
2948  "PhotonPrompt", "PhotonFragmentation", "PhotonISR" , "String" } ;
2949 
2950  for(Int_t i = 0; i < fNOriginHistograms; i++)
2951  {
2952  fhMCE[i] = new TH1F(Form("hE_MC%s",pname[i].Data()),
2953  Form("cluster from %s : E ",ptype[i].Data()),
2954  nptbins,ptmin,ptmax);
2955  fhMCE[i]->SetXTitle("#it{E} (GeV)");
2956  outputContainer->Add(fhMCE[i]) ;
2957 
2958  fhMCPt[i] = new TH1F(Form("hPt_MC%s",pname[i].Data()),
2959  Form("cluster from %s : #it{p}_{T} ",ptype[i].Data()),
2960  nptbins,ptmin,ptmax);
2961  fhMCPt[i]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2962  outputContainer->Add(fhMCPt[i]) ;
2963 
2964  fhMCEta[i] = new TH2F(Form("hEta_MC%s",pname[i].Data()),
2965  Form("cluster from %s : #eta ",ptype[i].Data()),
2966  nptbins,ptmin,ptmax,netabins,etamin,etamax);
2967  fhMCEta[i]->SetYTitle("#eta");
2968  fhMCEta[i]->SetXTitle("#it{E} (GeV)");
2969  outputContainer->Add(fhMCEta[i]) ;
2970 
2971  fhMCPhi[i] = new TH2F(Form("hPhi_MC%s",pname[i].Data()),
2972  Form("cluster from %s : #phi ",ptype[i].Data()),
2973  nptbins,ptmin,ptmax,nphibins,phimin,phimax);
2974  fhMCPhi[i]->SetYTitle("#phi (rad)");
2975  fhMCPhi[i]->SetXTitle("#it{E} (GeV)");
2976  outputContainer->Add(fhMCPhi[i]) ;
2977 
2978 
2979  fhMCDeltaE[i] = new TH2F (Form("hDeltaE_MC%s",pname[i].Data()),
2980  Form("MC - Reco E from %s",pname[i].Data()),
2981  nptbins,ptmin,ptmax, 200,-50,50);
2982  fhMCDeltaE[i]->SetYTitle("#Delta #it{E} (GeV)");
2983  fhMCDeltaE[i]->SetXTitle("#it{E} (GeV)");
2984  outputContainer->Add(fhMCDeltaE[i]);
2985 
2986  fhMCDeltaPt[i] = new TH2F (Form("hDeltaPt_MC%s",pname[i].Data()),
2987  Form("MC - Reco #it{p}_{T} from %s",pname[i].Data()),
2988  nptbins,ptmin,ptmax, 200,-50,50);
2989  fhMCDeltaPt[i]->SetXTitle("p_{T,rec} (GeV/#it{c})");
2990  fhMCDeltaPt[i]->SetYTitle("#Delta #it{p}_{T} (GeV/#it{c})");
2991  outputContainer->Add(fhMCDeltaPt[i]);
2992 
2993  fhMC2E[i] = new TH2F (Form("h2E_MC%s",pname[i].Data()),
2994  Form("E distribution, reconstructed vs generated from %s",pname[i].Data()),
2995  nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
2996  fhMC2E[i]->SetXTitle("#it{E}_{rec} (GeV)");
2997  fhMC2E[i]->SetYTitle("#it{E}_{gen} (GeV)");
2998  outputContainer->Add(fhMC2E[i]);
2999 
3000  fhMC2Pt[i] = new TH2F (Form("h2Pt_MC%s",pname[i].Data()),
3001  Form("p_T distribution, reconstructed vs generated from %s",pname[i].Data()),
3002  nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
3003  fhMC2Pt[i]->SetXTitle("p_{T,rec} (GeV/#it{c})");
3004  fhMC2Pt[i]->SetYTitle("p_{T,gen} (GeV/#it{c})");
3005  outputContainer->Add(fhMC2Pt[i]);
3006  }
3007 
3008  TString pptype[] = { "#gamma" , "#gamma_{#pi decay}" ,
3009  "#gamma_{#eta decay}", "#gamma_{other decay}" ,
3010  "#gamma_{prompt}" , "#gamma_{fragmentation}", "#gamma_{ISR}" } ;
3011 
3012  TString ppname[] = { "Photon" , "PhotonPi0Decay" ,
3013  "PhotonEtaDecay", "PhotonOtherDecay" ,
3014  "PhotonPrompt" , "PhotonFragmentation", "PhotonISR" } ;
3015 
3016  for(Int_t i = 0; i < fNPrimaryHistograms; i++)
3017  {
3018  fhEPrimMC[i] = new TH1F(Form("hEPrim_MC%s",ppname[i].Data()),
3019  Form("primary photon %s : E ",pptype[i].Data()),
3020  nptbins,ptmin,ptmax);
3021  fhEPrimMC[i]->SetXTitle("#it{E} (GeV)");
3022  outputContainer->Add(fhEPrimMC[i]) ;
3023 
3024  fhPtPrimMC[i] = new TH1F(Form("hPtPrim_MC%s",ppname[i].Data()),
3025  Form("primary photon %s : #it{p}_{T} ",pptype[i].Data()),
3026  nptbins,ptmin,ptmax);
3027  fhPtPrimMC[i]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
3028  outputContainer->Add(fhPtPrimMC[i]) ;
3029 
3030  fhYPrimMC[i] = new TH2F(Form("hYPrim_MC%s",ppname[i].Data()),
3031  Form("primary photon %s : Rapidity ",pptype[i].Data()),
3032  nptbins,ptmin,ptmax,200,-2,2);
3033  fhYPrimMC[i]->SetYTitle("Rapidity");
3034  fhYPrimMC[i]->SetXTitle("#it{E} (GeV)");
3035  outputContainer->Add(fhYPrimMC[i]) ;
3036 
3037  fhEtaPrimMC[i] = new TH2F(Form("hEtaPrim_MC%s",ppname[i].Data()),
3038  Form("primary photon %s : #eta",pptype[i].Data()),
3039  nptbins,ptmin,ptmax,200,-2,2);
3040  fhEtaPrimMC[i]->SetYTitle("#eta");
3041  fhEtaPrimMC[i]->SetXTitle("#it{E} (GeV)");
3042  outputContainer->Add(fhEtaPrimMC[i]) ;
3043 
3044  fhPhiPrimMC[i] = new TH2F(Form("hPhiPrim_MC%s",ppname[i].Data()),
3045  Form("primary photon %s : #phi ",pptype[i].Data()),
3046  nptbins,ptmin,ptmax,nphibins,0,TMath::TwoPi());
3047  fhPhiPrimMC[i]->SetYTitle("#phi (rad)");
3048  fhPhiPrimMC[i]->SetXTitle("#it{E} (GeV)");
3049  outputContainer->Add(fhPhiPrimMC[i]) ;
3050 
3051 
3052  fhEPrimMCAcc[i] = new TH1F(Form("hEPrimAcc_MC%s",ppname[i].Data()),
3053  Form("primary photon %s in acceptance: E ",pptype[i].Data()),
3054  nptbins,ptmin,ptmax);
3055  fhEPrimMCAcc[i]->SetXTitle("#it{E} (GeV)");
3056  outputContainer->Add(fhEPrimMCAcc[i]) ;
3057 
3058  fhPtPrimMCAcc[i] = new TH1F(Form("hPtPrimAcc_MC%s",ppname[i].Data()),
3059  Form("primary photon %s in acceptance: #it{p}_{T} ",pptype[i].Data()),
3060  nptbins,ptmin,ptmax);
3061  fhPtPrimMCAcc[i]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
3062  outputContainer->Add(fhPtPrimMCAcc[i]) ;
3063 
3064  fhYPrimMCAcc[i] = new TH2F(Form("hYPrimAcc_MC%s",ppname[i].Data()),
3065  Form("primary photon %s in acceptance: Rapidity ",pptype[i].Data()),
3066  nptbins,ptmin,ptmax,100,-1,1);
3067  fhYPrimMCAcc[i]->SetYTitle("Rapidity");
3068  fhYPrimMCAcc[i]->SetXTitle("#it{E} (GeV)");
3069  outputContainer->Add(fhYPrimMCAcc[i]) ;
3070 
3071  fhEtaPrimMCAcc[i] = new TH2F(Form("hEtaPrimAcc_MC%s",ppname[i].Data()),
3072  Form("primary photon %s in acceptance: #eta ",pptype[i].Data()),
3073  nptbins,ptmin,ptmax,netabins,etamin,etamax);
3074  fhEtaPrimMCAcc[i]->SetYTitle("#eta");
3075  fhEtaPrimMCAcc[i]->SetXTitle("#it{E} (GeV)");
3076  outputContainer->Add(fhEtaPrimMCAcc[i]) ;
3077 
3078  fhPhiPrimMCAcc[i] = new TH2F(Form("hPhiPrimAcc_MC%s",ppname[i].Data()),
3079  Form("primary photon %s in acceptance: #phi ",pptype[i].Data()),
3080  nptbins,ptmin,ptmax,nphibins,phimin,phimax);
3081  fhPhiPrimMCAcc[i]->SetYTitle("#phi (rad)");
3082  fhPhiPrimMCAcc[i]->SetXTitle("#it{E} (GeV)");
3083  outputContainer->Add(fhPhiPrimMCAcc[i]) ;
3084  }
3085 
3086  if(fFillSSHistograms)
3087  {
3088  TString ptypess[] = { "#gamma","hadron?","#pi^{0}","#eta","#gamma->e^{#pm}","e^{#pm}"} ;
3089 
3090  TString pnamess[] = { "Photon","Hadron","Pi0","Eta","Conversion","Electron"} ;
3091 
3092  for(Int_t i = 0; i < 6; i++)
3093  {
3094  fhMCELambda0[i] = new TH2F(Form("hELambda0_MC%s",pnamess[i].Data()),
3095  Form("cluster from %s : E vs #lambda_{0}^{2}",ptypess[i].Data()),
3096  nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
3097  fhMCELambda0[i]->SetYTitle("#lambda_{0}^{2}");
3098  fhMCELambda0[i]->SetXTitle("#it{E} (GeV)");
3099  outputContainer->Add(fhMCELambda0[i]) ;
3100 
3101  fhMCPtLambda0[i] = new TH2F(Form("hPtLambda0_MC%s",pnamess[i].Data()),
3102  Form("cluster from %s : #it{p}_{T} vs #lambda_{0}^{2}",ptypess[i].Data()),
3103  nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
3104  fhMCPtLambda0[i]->SetYTitle("#lambda_{0}^{2}");
3105  fhMCPtLambda0[i]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
3106  outputContainer->Add(fhMCPtLambda0[i]) ;
3107 
3108  fhMCELambda1[i] = new TH2F(Form("hELambda1_MC%s",pnamess[i].Data()),
3109  Form("cluster from %s : E vs #lambda_{1}^{2}",ptypess[i].Data()),
3110  nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
3111  fhMCELambda1[i]->SetYTitle("#lambda_{1}^{2}");
3112  fhMCELambda1[i]->SetXTitle("#it{E} (GeV)");
3113  outputContainer->Add(fhMCELambda1[i]) ;
3114 
3115  fhMCEDispersion[i] = new TH2F(Form("hEDispersion_MC%s",pnamess[i].Data()),
3116  Form("cluster from %s : E vs dispersion^{2}",ptypess[i].Data()),
3117  nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
3118  fhMCEDispersion[i]->SetYTitle("D^{2}");
3119  fhMCEDispersion[i]->SetXTitle("#it{E} (GeV)");
3120  outputContainer->Add(fhMCEDispersion[i]) ;
3121 
3122  fhMCNCellsE[i] = new TH2F (Form("hNCellsE_MC%s",pnamess[i].Data()),
3123  Form("# of cells in cluster from %s vs E of clusters",ptypess[i].Data()),
3124  nptbins,ptmin,ptmax, nbins,nmin,nmax);
3125  fhMCNCellsE[i]->SetXTitle("#it{E} (GeV)");
3126  fhMCNCellsE[i]->SetYTitle("# of cells in cluster");
3127  outputContainer->Add(fhMCNCellsE[i]);
3128 
3129  fhMCMaxCellDiffClusterE[i] = new TH2F (Form("hMaxCellDiffClusterE_MC%s",pnamess[i].Data()),
3130  Form("energy vs difference of cluster energy from %s - max cell energy / cluster energy, good clusters",ptypess[i].Data()),
3131  nptbins,ptmin,ptmax, 500,0,1.);
3132  fhMCMaxCellDiffClusterE[i]->SetXTitle("#it{E}_{cluster} (GeV) ");
3133  fhMCMaxCellDiffClusterE[i]->SetYTitle("(#it{E}_{cluster} - #it{E}_{cell max})/ #it{E}_{cluster}");
3134  outputContainer->Add(fhMCMaxCellDiffClusterE[i]);
3135 
3137  {
3138  fhMCLambda0vsClusterMaxCellDiffE0[i] = new TH2F(Form("hLambda0vsClusterMaxCellDiffE0_MC%s",pnamess[i].Data()),
3139  Form("cluster from %s : #lambda^{2}_{0} vs fraction of energy carried by max cell, E < 2 GeV",ptypess[i].Data()),
3140  ssbins,ssmin,ssmax,500,0,1.);
3141  fhMCLambda0vsClusterMaxCellDiffE0[i]->SetXTitle("#lambda_{0}^{2}");
3142  fhMCLambda0vsClusterMaxCellDiffE0[i]->SetYTitle("(#it{E}_{cluster} - #it{E}_{cell max})/ #it{E}_{cluster}");
3143  outputContainer->Add(fhMCLambda0vsClusterMaxCellDiffE0[i]) ;
3144 
3145  fhMCLambda0vsClusterMaxCellDiffE2[i] = new TH2F(Form("hLambda0vsClusterMaxCellDiffE2_MC%s",pnamess[i].Data()),
3146  Form("cluster from %s : #lambda^{2}_{0} vs fraction of energy carried by max cell, 2< E < 6 GeV",ptypess[i].Data()),
3147  ssbins,ssmin,ssmax,500,0,1.);
3148  fhMCLambda0vsClusterMaxCellDiffE2[i]->SetXTitle("#lambda_{0}^{2}");
3149  fhMCLambda0vsClusterMaxCellDiffE2[i]->SetYTitle("(#it{E}_{cluster} - #it{E}_{cell max})/ #it{E}_{cluster}");
3150  outputContainer->Add(fhMCLambda0vsClusterMaxCellDiffE2[i]) ;
3151 
3152  fhMCLambda0vsClusterMaxCellDiffE6[i] = new TH2F(Form("hLambda0vsClusterMaxCellDiffE6_MC%s",pnamess[i].Data()),
3153  Form("cluster from %s : #lambda^{2}_{0} vs fraction of energy carried by max cell, #it{E} > 6 GeV",ptypess[i].Data()),
3154  ssbins,ssmin,ssmax,500,0,1.);
3155  fhMCLambda0vsClusterMaxCellDiffE6[i]->SetXTitle("#lambda_{0}^{2}");
3156  fhMCLambda0vsClusterMaxCellDiffE6[i]->SetYTitle("(#it{E}_{cluster} - #it{E}_{cell max})/ #it{E}_{cluster}");
3157  outputContainer->Add(fhMCLambda0vsClusterMaxCellDiffE6[i]) ;
3158 
3159  fhMCNCellsvsClusterMaxCellDiffE0[i] = new TH2F(Form("hNCellsvsClusterMaxCellDiffE0_MC%s",pnamess[i].Data()),
3160  Form("cluster from %s : N cells in cluster vs fraction of energy carried by max cell, E < 2 GeV",ptypess[i].Data()),
3161  nbins/5,nmin,nmax/5,500,0,1.);
3162  fhMCNCellsvsClusterMaxCellDiffE0[i]->SetXTitle("N cells in cluster");
3163  fhMCNCellsvsClusterMaxCellDiffE0[i]->SetYTitle("(#it{E}_{cluster} - #it{E}_{cell max})/ #it{E}_{cluster}");
3164  outputContainer->Add(fhMCNCellsvsClusterMaxCellDiffE0[i]) ;
3165 
3166  fhMCNCellsvsClusterMaxCellDiffE2[i] = new TH2F(Form("hNCellsvsClusterMaxCellDiffE2_MC%s",pnamess[i].Data()),
3167  Form("cluster from %s : N cells in cluster vs fraction of energy carried by max cell, 2< E < 6 GeV",ptypess[i].Data()),
3168  nbins/5,nmin,nmax/5,500,0,1.);
3169  fhMCNCellsvsClusterMaxCellDiffE2[i]->SetXTitle("N cells in cluster");
3170  fhMCNCellsvsClusterMaxCellDiffE2[i]->SetYTitle("(#it{E}_{cluster} - #it{E}_{cell max})/ #it{E}_{cluster}");
3171  outputContainer->Add(fhMCNCellsvsClusterMaxCellDiffE2[i]) ;
3172 
3173  fhMCNCellsvsClusterMaxCellDiffE6[i] = new TH2F(Form("hNCellsvsClusterMaxCellDiffE6_MC%s",pnamess[i].Data()),
3174  Form("cluster from %s : N cells in cluster vs fraction of energy carried by max cell, #it{E} > 6 GeV",ptypess[i].Data()),
3175  nbins/5,nmin,nmax/5,500,0,1.);
3176  fhMCNCellsvsClusterMaxCellDiffE6[i]->SetXTitle("N cells in cluster");
3177  fhMCNCellsvsClusterMaxCellDiffE6[i]->SetYTitle("#it{E} (GeV)");
3178  outputContainer->Add(fhMCNCellsvsClusterMaxCellDiffE6[i]) ;
3179 
3180  if(GetCalorimeter()==kEMCAL)
3181  {
3182  fhMCEDispEta[i] = new TH2F (Form("hEDispEtaE_MC%s",pnamess[i].Data()),
3183  Form("cluster from %s : #sigma^{2}_{#eta #eta} = #Sigma w_{i}(#eta_{i} - <#eta>)^{2}/ #Sigma w_{i} vs E",ptypess[i].Data()),
3184  nptbins,ptmin,ptmax, ssbins,ssmin,ssmax);
3185  fhMCEDispEta[i]->SetXTitle("#it{E} (GeV)");
3186  fhMCEDispEta[i]->SetYTitle("#sigma^{2}_{#eta #eta}");
3187  outputContainer->Add(fhMCEDispEta[i]);
3188 
3189  fhMCEDispPhi[i] = new TH2F (Form("hEDispPhiE_MC%s",pnamess[i].Data()),
3190  Form("cluster from %s : #sigma^{2}_{#phi #phi} = #Sigma w_{i}(#phi_{i} - <#phi>)^{2} / #Sigma w_{i} vs E",ptypess[i].Data()),
3191  nptbins,ptmin,ptmax, ssbins,ssmin,ssmax);
3192  fhMCEDispPhi[i]->SetXTitle("#it{E} (GeV)");
3193  fhMCEDispPhi[i]->SetYTitle("#sigma^{2}_{#phi #phi}");
3194  outputContainer->Add(fhMCEDispPhi[i]);
3195 
3196  fhMCESumEtaPhi[i] = new TH2F (Form("hESumEtaPhiE_MC%s",pnamess[i].Data()),
3197  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()),
3198  nptbins,ptmin,ptmax, 2*ssbins,-ssmax,ssmax);
3199  fhMCESumEtaPhi[i]->SetXTitle("#it{E} (GeV)");
3200  fhMCESumEtaPhi[i]->SetYTitle("#delta^{2}_{#eta #phi}");
3201  outputContainer->Add(fhMCESumEtaPhi[i]);
3202 
3203  fhMCEDispEtaPhiDiff[i] = new TH2F (Form("hEDispEtaPhiDiffE_MC%s",pnamess[i].Data()),
3204  Form("cluster from %s : #sigma^{2}_{#phi #phi} - #sigma^{2}_{#eta #eta} vs E",ptypess[i].Data()),
3205  nptbins,ptmin,ptmax,200,-10,10);
3206  fhMCEDispEtaPhiDiff[i]->SetXTitle("#it{E} (GeV)");
3207  fhMCEDispEtaPhiDiff[i]->SetYTitle("#sigma^{2}_{#phi #phi}-#sigma^{2}_{#eta #eta}");
3208  outputContainer->Add(fhMCEDispEtaPhiDiff[i]);
3209 
3210  fhMCESphericity[i] = new TH2F (Form("hESphericity_MC%s",pnamess[i].Data()),
3211  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()),
3212  nptbins,ptmin,ptmax, 200,-1,1);
3213  fhMCESphericity[i]->SetXTitle("#it{E} (GeV)");
3214  fhMCESphericity[i]->SetYTitle("s = (#sigma^{2}_{#phi #phi} - #sigma^{2}_{#eta #eta}) / (#sigma^{2}_{#eta #eta} + #sigma^{2}_{#phi #phi})");
3215  outputContainer->Add(fhMCESphericity[i]);
3216 
3217  for(Int_t ie = 0; ie < 7; ie++)
3218  {
3219  fhMCDispEtaDispPhi[ie][i] = new TH2F (Form("hMCDispEtaDispPhi_EBin%d_MC%s",ie,pnamess[i].Data()),
3220  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]),
3221  ssbins,ssmin,ssmax , ssbins,ssmin,ssmax);
3222  fhMCDispEtaDispPhi[ie][i]->SetXTitle("#sigma^{2}_{#eta #eta}");
3223  fhMCDispEtaDispPhi[ie][i]->SetYTitle("#sigma^{2}_{#phi #phi}");
3224  outputContainer->Add(fhMCDispEtaDispPhi[ie][i]);
3225 
3226  fhMCLambda0DispEta[ie][i] = new TH2F (Form("hMCLambda0DispEta_EBin%d_MC%s",ie,pnamess[i].Data()),
3227  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]),
3228  ssbins,ssmin,ssmax , ssbins,ssmin,ssmax);
3229  fhMCLambda0DispEta[ie][i]->SetXTitle("#lambda^{2}_{0}");
3230  fhMCLambda0DispEta[ie][i]->SetYTitle("#sigma^{2}_{#phi #phi}");
3231  outputContainer->Add(fhMCLambda0DispEta[ie][i]);
3232 
3233  fhMCLambda0DispPhi[ie][i] = new TH2F (Form("hMCLambda0DispPhi_EBin%d_MC%s",ie,pnamess[i].Data()),
3234  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]),
3235  ssbins,ssmin,ssmax , ssbins,ssmin,ssmax);
3236  fhMCLambda0DispPhi[ie][i]->SetXTitle("#lambda^{2}_{0}");
3237  fhMCLambda0DispPhi[ie][i]->SetYTitle("#sigma^{2}_{#phi #phi}");
3238  outputContainer->Add(fhMCLambda0DispPhi[ie][i]);
3239  }
3240  }
3241  }
3242  }// loop
3243 
3244  if(!GetReader()->IsEmbeddedClusterSelectionOn())
3245  {
3246  fhMCPhotonELambda0NoOverlap = new TH2F("hELambda0_MCPhoton_NoOverlap",
3247  "cluster from Photon : E vs #lambda_{0}^{2}",
3248  nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
3249  fhMCPhotonELambda0NoOverlap->SetYTitle("#lambda_{0}^{2}");
3250  fhMCPhotonELambda0NoOverlap->SetXTitle("#it{E} (GeV)");
3251  outputContainer->Add(fhMCPhotonELambda0NoOverlap) ;
3252 
3253  fhMCPhotonELambda0TwoOverlap = new TH2F("hELambda0_MCPhoton_TwoOverlap",
3254  "cluster from Photon : E vs #lambda_{0}^{2}",
3255  nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
3256  fhMCPhotonELambda0TwoOverlap->SetYTitle("#lambda_{0}^{2}");
3257  fhMCPhotonELambda0TwoOverlap->SetXTitle("#it{E} (GeV)");
3258  outputContainer->Add(fhMCPhotonELambda0TwoOverlap) ;
3259 
3260  fhMCPhotonELambda0NOverlap = new TH2F("hELambda0_MCPhoton_NOverlap",
3261  "cluster from Photon : E vs #lambda_{0}^{2}",
3262  nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
3263  fhMCPhotonELambda0NOverlap->SetYTitle("#lambda_{0}^{2}");
3264  fhMCPhotonELambda0NOverlap->SetXTitle("#it{E} (GeV)");
3265  outputContainer->Add(fhMCPhotonELambda0NOverlap) ;
3266  } // No embedding
3267 
3268  if(GetReader()->IsEmbeddedClusterSelectionOn())
3269  {
3270  fhEmbeddedSignalFractionEnergy = new TH2F("hEmbeddedSignal_FractionEnergy",
3271  "Energy Fraction of embedded signal versus cluster energy",
3272  nptbins,ptmin,ptmax,100,0.,1.);
3273  fhEmbeddedSignalFractionEnergy->SetYTitle("Fraction");
3274  fhEmbeddedSignalFractionEnergy->SetXTitle("#it{E} (GeV)");
3275  outputContainer->Add(fhEmbeddedSignalFractionEnergy) ;
3276 
3277  fhEmbedPhotonELambda0FullSignal = new TH2F("hELambda0_EmbedPhoton_FullSignal",
3278  "cluster from Photon embedded with more than 90% energy in cluster : E vs #lambda_{0}^{2}",
3279  nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
3280  fhEmbedPhotonELambda0FullSignal->SetYTitle("#lambda_{0}^{2}");
3281  fhEmbedPhotonELambda0FullSignal->SetXTitle("#it{E} (GeV)");
3282  outputContainer->Add(fhEmbedPhotonELambda0FullSignal) ;
3283 
3284  fhEmbedPhotonELambda0MostlySignal = new TH2F("hELambda0_EmbedPhoton_MostlySignal",
3285  "cluster from Photon embedded with 50% to 90% energy in cluster : E vs #lambda_{0}^{2}",
3286  nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
3287  fhEmbedPhotonELambda0MostlySignal->SetYTitle("#lambda_{0}^{2}");
3288  fhEmbedPhotonELambda0MostlySignal->SetXTitle("#it{E} (GeV)");
3289  outputContainer->Add(fhEmbedPhotonELambda0MostlySignal) ;
3290 
3291  fhEmbedPhotonELambda0MostlyBkg = new TH2F("hELambda0_EmbedPhoton_MostlyBkg",
3292  "cluster from Photon embedded with 10% to 50% energy in cluster : E vs #lambda_{0}^{2}",
3293  nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
3294  fhEmbedPhotonELambda0MostlyBkg->SetYTitle("#lambda_{0}^{2}");
3295  fhEmbedPhotonELambda0MostlyBkg->SetXTitle("#it{E} (GeV)");
3296  outputContainer->Add(fhEmbedPhotonELambda0MostlyBkg) ;
3297 
3298  fhEmbedPhotonELambda0FullBkg = new TH2F("hELambda0_EmbedPhoton_FullBkg",
3299  "cluster from Photonm embedded with 0% to 10% energy in cluster : E vs #lambda_{0}^{2}",
3300  nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
3301  fhEmbedPhotonELambda0FullBkg->SetYTitle("#lambda_{0}^{2}");
3302  fhEmbedPhotonELambda0FullBkg->SetXTitle("#it{E} (GeV)");
3303  outputContainer->Add(fhEmbedPhotonELambda0FullBkg) ;
3304 
3305  fhEmbedPi0ELambda0FullSignal = new TH2F("hELambda0_EmbedPi0_FullSignal",
3306  "cluster from Pi0 embedded with more than 90% energy in cluster : E vs #lambda_{0}^{2}",
3307  nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
3308  fhEmbedPi0ELambda0FullSignal->SetYTitle("#lambda_{0}^{2}");
3309  fhEmbedPi0ELambda0FullSignal->SetXTitle("#it{E} (GeV)");
3310  outputContainer->Add(fhEmbedPi0ELambda0FullSignal) ;
3311 
3312  fhEmbedPi0ELambda0MostlySignal = new TH2F("hELambda0_EmbedPi0_MostlySignal",
3313  "cluster from Pi0 embedded with 50% to 90% energy in cluster : E vs #lambda_{0}^{2}",
3314  nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
3315  fhEmbedPi0ELambda0MostlySignal->SetYTitle("#lambda_{0}^{2}");
3316  fhEmbedPi0ELambda0MostlySignal->SetXTitle("#it{E} (GeV)");
3317  outputContainer->Add(fhEmbedPi0ELambda0MostlySignal) ;
3318 
3319  fhEmbedPi0ELambda0MostlyBkg = new TH2F("hELambda0_EmbedPi0_MostlyBkg",
3320  "cluster from Pi0 embedded with 10% to 50% energy in cluster : E vs #lambda_{0}^{2}",
3321  nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
3322  fhEmbedPi0ELambda0MostlyBkg->SetYTitle("#lambda_{0}^{2}");
3323  fhEmbedPi0ELambda0MostlyBkg->SetXTitle("#it{E} (GeV)");
3324  outputContainer->Add(fhEmbedPi0ELambda0MostlyBkg) ;
3325 
3326  fhEmbedPi0ELambda0FullBkg = new TH2F("hELambda0_EmbedPi0_FullBkg",
3327  "cluster from Pi0 embedded with 0% to 10% energy in cluster : E vs #lambda_{0}^{2}",
3328  nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
3329  fhEmbedPi0ELambda0FullBkg->SetYTitle("#lambda_{0}^{2}");
3330  fhEmbedPi0ELambda0FullBkg->SetXTitle("#it{E} (GeV)");
3331  outputContainer->Add(fhEmbedPi0ELambda0FullBkg) ;
3332  }// embedded histograms
3333  }// Fill SS MC histograms
3334 
3335 
3337  {
3338  fhMCConversionVertex = new TH2F("hMCPhotonConversionVertex","cluster from converted photon, #it{p}_{T} vs vertex distance",
3339  nptbins,ptmin,ptmax,500,0,500);
3340  fhMCConversionVertex->SetYTitle("#it{R} (cm)");
3341  fhMCConversionVertex->SetXTitle("#it{p}_{T} (GeV/#it{c})");
3342  outputContainer->Add(fhMCConversionVertex) ;
3343 
3344  if(GetCalorimeter() == kEMCAL && GetFirstSMCoveredByTRD() >= 0)
3345  {
3346  fhMCConversionVertexTRD = new TH2F("hMCPhotonConversionVertexTRD","cluster from converted photon, #it{p}_{T} vs vertex distance, SM covered by TRD",
3347  nptbins,ptmin,ptmax,500,0,500);
3348  fhMCConversionVertexTRD->SetYTitle("#it{R} (cm)");
3349  fhMCConversionVertexTRD->SetXTitle("#it{p}_{T} (GeV/#it{c})");
3350  outputContainer->Add(fhMCConversionVertexTRD) ;
3351  }
3352 
3353  if(fFillSSHistograms)
3354  {
3355  TString region[] = {"ITS","TPC","TRD","TOF","Top EMCal","In EMCal"};
3356  for(Int_t iR = 0; iR < 6; iR++)
3357  {
3358  fhMCConversionLambda0Rcut[iR] = new TH2F(Form("hMCPhotonConversionLambda0_R%d",iR),
3359  Form("cluster from converted photon, #it{p}_{T} vs #lambda_{0}^{2}, conversion in %s",region[iR].Data()),
3360  nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
3361  fhMCConversionLambda0Rcut[iR]->SetYTitle("#lambda_{0}^{2}");
3362  fhMCConversionLambda0Rcut[iR]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
3363  outputContainer->Add(fhMCConversionLambda0Rcut[iR]) ;
3364 
3365  fhMCConversionLambda1Rcut[iR] = new TH2F(Form("hMCPhotonConversionLambda1_R%d",iR),
3366  Form("cluster from converted photon, #it{p}_{T} vs #lambda_{1}^{2}, conversion in %s",region[iR].Data()),
3367  nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
3368  fhMCConversionLambda1Rcut[iR]->SetYTitle("#lambda_{1}^{2}");
3369  fhMCConversionLambda1Rcut[iR]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
3370  outputContainer->Add(fhMCConversionLambda1Rcut[iR]) ;
3371  } // R cut
3372 
3373 
3374  if(GetCalorimeter() == kEMCAL && GetFirstSMCoveredByTRD() >= 0)
3375  {
3376  for(Int_t iR = 0; iR < 6; iR++)
3377  {
3378  fhMCConversionLambda0RcutTRD[iR] = new TH2F(Form("hMCPhotonConversionLambda0TRD_R%d",iR),
3379  Form("cluster from converted photon, #it{p}_{T} vs #lambda_{0}^{2}, conversion in %s, SM covered by TRD",region[iR].Data()),
3380  nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
3381  fhMCConversionLambda0RcutTRD[iR]->SetYTitle("#lambda_{0}^{2}");
3382  fhMCConversionLambda0RcutTRD[iR]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
3383  outputContainer->Add(fhMCConversionLambda0RcutTRD[iR]) ;
3384 
3385  fhMCConversionLambda1RcutTRD[iR] = new TH2F(Form("hMCPhotonConversionLambda1TRD_R%d",iR),
3386  Form("cluster from converted photon, #it{p}_{T} vs #lambda_{1}^{2}, conversion in %s, SM covered by TRD",region[iR].Data()),
3387  nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
3388  fhMCConversionLambda1RcutTRD[iR]->SetYTitle("#lambda_{1}^{2}");
3389  fhMCConversionLambda1RcutTRD[iR]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
3390  outputContainer->Add(fhMCConversionLambda1RcutTRD[iR]) ;
3391  } // R cut
3392  }
3393 
3394  // if(GetCalorimeter() == kEMCAL && fFillEMCALRegionSSHistograms)
3395  // {
3396  // for(Int_t ieta = 0; ieta < 4; ieta++)
3397  // {
3398  // for(Int_t iphi = 0; iphi < 3; iphi++)
3399  // {
3400  // for(Int_t iR = 0; iR < 6; iR++)
3401  // {
3402  // fhLam0EMCALRegionMCConvRcut[ieta][iphi][iR] =
3403  // new TH2F(Form("hMCPhotonConversionLambda0_R%d_eta%d_phi%d",iR,ieta,iphi),
3404  // Form("cluster from converted photon, #it{p}_{T} vs #lambda_{0}^{2}, conversion in %s, region eta %d, phi %d",region[iR].Data(),ieta,iphi),
3405  // nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
3406  // fhLam0EMCALRegionMCConvRcut[ieta][iphi][iR]->SetYTitle("#lambda_{0}^{2}");
3407  // fhLam0EMCALRegionMCConvRcut[ieta][iphi][iR]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
3408  // outputContainer->Add(fhLam0EMCALRegionMCConvRcut[ieta][iphi][iR]) ;
3409  //
3410  // if(GetFirstSMCoveredByTRD() >= 0)
3411  // {
3412  // fhLam0EMCALRegionTRDMCConvRcut[ieta][iphi][iR] =
3413  // new TH2F(Form("hMCPhotonConversionLambda0TRD_R%d_eta%d_phi%d",iR,ieta,iphi),
3414  // Form("cluster from converted photon, #it{p}_{T} vs #lambda_{0}^{2}, conversion in %s, region eta %d, phi %d, SM covered by TRD",region[iR].Data(),ieta,iphi),
3415  // nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
3416  // fhLam0EMCALRegionTRDMCConvRcut[ieta][iphi][iR]->SetYTitle("#lambda_{0}^{2}");
3417  // fhLam0EMCALRegionTRDMCConvRcut[ieta][iphi][iR]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
3418  // outputContainer->Add(fhLam0EMCALRegionTRDMCConvRcut[ieta][iphi][iR]) ;
3419  // } // TRD
3420  //
3421  // } // iR
3422  // } // iphi
3423  // } // ieta
3424  // } // regions in EMCal
3425  } // shower shape
3426  } // conversion vertex
3427  } // Histos with MC
3428 
3430  {
3431  for(Int_t ie=0; ie<fNEBinCuts; ie++)
3432  {
3433  fhEBinClusterEtaPhi[ie] = new TH2F
3434  (Form("hEBin%d_Cluster_EtaPhi",ie),
3435  Form("#eta vs #phi, cluster, %2.2f<#it{p}_{T}<%2.2f GeV/#it{c}",fEBinCuts[ie],fEBinCuts[ie+1]),
3436  netabins,etamin,etamax,nphibins,phimin,phimax);
3437  fhEBinClusterEtaPhi[ie]->SetYTitle("#phi (rad)");
3438  fhEBinClusterEtaPhi[ie]->SetXTitle("#eta");
3439  outputContainer->Add(fhEBinClusterEtaPhi[ie]) ;
3440 
3441  fhEBinClusterColRow[ie] = new TH2F
3442  (Form("hEBin%d_Cluster_ColRow",ie),
3443  Form("column vs row, cluster max E cell, %2.2f<#it{p}_{T}<%2.2f GeV/#it{c}",fEBinCuts[ie],fEBinCuts[ie+1]),
3444  96+2,-1.5,96+0.5,(8*24+2*8)+2,-1.5,(8*24+2*8)+0.5);
3445  fhEBinClusterColRow[ie]->SetYTitle("row");
3446  fhEBinClusterColRow[ie]->SetXTitle("column");
3447  outputContainer->Add(fhEBinClusterColRow[ie]) ;
3448 
3449  fhEBinClusterEtaPhiPID[ie] = new TH2F
3450  (Form("hEBin%d_Cluster_EtaPhi_PID",ie),
3451  Form("#eta vs #phi, cluster, %2.2f<#it{p}_{T}<%2.2f GeV/#it{c}, PID cut",fEBinCuts[ie],fEBinCuts[ie+1]),
3452  netabins,etamin,etamax,nphibins,phimin,phimax);
3453  fhEBinClusterEtaPhiPID[ie]->SetYTitle("#phi (rad)");
3454  fhEBinClusterEtaPhiPID[ie]->SetXTitle("#eta");
3455  outputContainer->Add(fhEBinClusterEtaPhiPID[ie]) ;
3456 
3457  fhEBinClusterColRowPID[ie] = new TH2F
3458  (Form("hEBin%d_Cluster_ColRow_PID",ie),
3459  Form("column vs row, cluster max E cell, %2.2f<#it{p}_{T}<%2.2f GeV/#it{c}, PID cut",fEBinCuts[ie],fEBinCuts[ie+1]),
3460  96+2,-1.5,96+0.5,(8*24+2*8)+2,-1.5,(8*24+2*8)+0.5);
3461  fhEBinClusterColRowPID[ie]->SetYTitle("row");
3462  fhEBinClusterColRowPID[ie]->SetXTitle("column");
3463  outputContainer->Add(fhEBinClusterColRowPID[ie]) ;
3464  }
3465  }
3466 
3468  {
3469  TString caseTitle[] = {"","CleanCluster","MergedClusterHijingBkg","MergedClusterNotHijingBkg","MergedClusterHijingAndOtherBkg","MergedCluster"};
3470 
3471  for(Int_t icase = 0; icase < 6; icase++)
3472  {
3473  fhLocalRegionClusterEtaPhi[icase] = new TH2F
3474  (Form("hLocalRegionClusterEtaPhi%s",caseTitle[icase].Data()),
3475  "cluster,#it{E} > 0.5 GeV, #eta vs #phi",netabins,etamin,etamax,nphibins,phimin,phimax);
3476  fhLocalRegionClusterEtaPhi[icase]->SetYTitle("#phi (rad)");
3477  fhLocalRegionClusterEtaPhi[icase]->SetXTitle("#eta");
3478  outputContainer->Add(fhLocalRegionClusterEtaPhi[icase]) ;
3479 
3480  fhLocalRegionClusterEnergySum[icase] = new TH2F
3481  (Form("hLocalRegionClusterEnergySum%s",caseTitle[icase].Data()),
3482  "Sum of cluster energy around trigger cluster #it{E} with R=0.2",
3483  nptbins,ptmin,ptmax, 200,0,100);
3484  fhLocalRegionClusterEnergySum[icase]->SetXTitle("#it{E} (GeV)");
3485  fhLocalRegionClusterEnergySum[icase]->SetYTitle("#Sigma #it{E} (GeV)");
3486  outputContainer->Add(fhLocalRegionClusterEnergySum[icase]);
3487 
3489  (Form("hLocalRegionClusterMultiplicity%s",caseTitle[icase].Data()),
3490  "Cluster multiplicity around trigger cluster #it{E} with R=0.2",
3491  nptbins,ptmin,ptmax, 200,0,200);
3492  fhLocalRegionClusterMultiplicity[icase]->SetXTitle("#it{E} (GeV)");
3493  fhLocalRegionClusterMultiplicity[icase]->SetYTitle("Multiplicity");
3494  outputContainer->Add(fhLocalRegionClusterMultiplicity[icase]);
3495 
3497  {
3499  (Form("hLocalRegionClusterEnergySumPerCentrality%s",caseTitle[icase].Data()),
3500  "Sum of cluster energy around trigger cluster vs centrality with R=0.2",
3501  100,0,100, 200,0,100);
3502  fhLocalRegionClusterEnergySumPerCentrality[icase]->SetXTitle("Centrality");
3503  fhLocalRegionClusterEnergySumPerCentrality[icase]->SetYTitle("#Sigma #it{E} (GeV)");
3504  outputContainer->Add(fhLocalRegionClusterEnergySumPerCentrality[icase]);
3505 
3507  (Form("hLocalRegionClusterMultiplicityPerCentrality%s",caseTitle[icase].Data()),
3508  "Cluster multiplicity around trigger cluster vs centrality with R=0.2",
3509  100,0,100, 200,0,200);
3510  fhLocalRegionClusterMultiplicityPerCentrality[icase]->SetXTitle("Centrality");
3511  fhLocalRegionClusterMultiplicityPerCentrality[icase]->SetYTitle("Multiplicity");
3512  outputContainer->Add(fhLocalRegionClusterMultiplicityPerCentrality[icase]);
3513  }
3514 
3516  {
3518  (Form("hLocalRegionClusterEnergySumHijing%s",caseTitle[icase].Data()),
3519  "Sum of cluster energy (HIJING) around trigger cluster #it{E} with R=0.2",
3520  nptbins,ptmin,ptmax, 200,0,100);
3521  fhLocalRegionClusterEnergySumHijing[icase]->SetXTitle("#it{E} (GeV)");
3522  fhLocalRegionClusterEnergySumHijing[icase]->SetYTitle("#Sigma #it{E} (GeV)");
3523  outputContainer->Add(fhLocalRegionClusterEnergySumHijing[icase]);
3524 
3526  (Form("hLocalRegionClusterMultiplicityHijing%s",caseTitle[icase].Data()),
3527  "Cluster multiplicity (HIJING) around trigger cluster #it{E} with R=0.2",
3528  nptbins,ptmin,ptmax, 200,0,200);
3529  fhLocalRegionClusterMultiplicityHijing[icase]->SetXTitle("#it{E} (GeV)");
3530  fhLocalRegionClusterMultiplicityHijing[icase]->SetYTitle("Multiplicity");
3531  outputContainer->Add(fhLocalRegionClusterMultiplicityHijing[icase]);
3532 
3534  (Form("hLocalRegionClusterEnergySumAdded%s",caseTitle[icase].Data()),
3535  "Sum of cluster energy (not HIJING) around trigger cluster #it{E} with R=0.2",
3536  nptbins,ptmin,ptmax, 200,0,100);
3537  fhLocalRegionClusterEnergySumAdded[icase]->SetXTitle("#it{E} (GeV)");
3538  fhLocalRegionClusterEnergySumAdded[icase]->SetYTitle("#Sigma #it{E} (GeV)");
3539  outputContainer->Add(fhLocalRegionClusterEnergySumAdded[icase]);
3540 
3542  (Form("hLocalRegionClusterMultiplicityAdded%s",caseTitle[icase].Data()),
3543  "Cluster multiplicity (not HIJING) around trigger cluster #it{E} with R=0.2",
3544  nptbins,ptmin,ptmax, 200,0,200);
3545  fhLocalRegionClusterMultiplicityAdded[icase]->SetXTitle("#it{E} (GeV)");
3546  fhLocalRegionClusterMultiplicityAdded[icase]->SetYTitle("Multiplicity");
3547  outputContainer->Add(fhLocalRegionClusterMultiplicityAdded[icase]);
3548 
3549 
3551  {
3553  (Form("hLocalRegionClusterEnergySumPerCentralityHijing%s",caseTitle[icase].Data()),
3554  "Sum of cluster energy (HIJING) around trigger cluster vs centrality with R=0.2",
3555  100,0,100, 200,0,100);
3556  fhLocalRegionClusterEnergySumPerCentralityHijing[icase]->SetXTitle("Centrality");
3557  fhLocalRegionClusterEnergySumPerCentralityHijing[icase]->SetYTitle("#Sigma #it{E} (GeV)");
3558  outputContainer->Add(fhLocalRegionClusterEnergySumPerCentralityHijing[icase]);
3559 
3561  (Form("hLocalRegionClusterMultiplicityPerCentralityHijing%s",caseTitle[icase].Data()),
3562  "Cluster multiplicity (HIJING) around trigger cluster vs centrality with R=0.2",
3563  100,0,100, 200,0,200);
3564  fhLocalRegionClusterMultiplicityPerCentralityHijing[icase]->SetXTitle("Centrality");
3565  fhLocalRegionClusterMultiplicityPerCentralityHijing[icase]->SetYTitle("Multiplicity");
3566  outputContainer->Add(fhLocalRegionClusterMultiplicityPerCentralityHijing[icase]);
3567 
3568 
3570  (Form("hLocalRegionClusterEnergySumPerCentralityAdded%s",caseTitle[icase].Data()),
3571  "Sum of cluster energy (not HIJING) around trigger cluster vs centrality with R=0.2",
3572  100,0,100, 200,0,100);
3573  fhLocalRegionClusterEnergySumPerCentralityAdded[icase]->SetXTitle("Centrality");
3574  fhLocalRegionClusterEnergySumPerCentralityAdded[icase]->SetYTitle("#Sigma #it{E} (GeV)");
3575  outputContainer->Add(fhLocalRegionClusterEnergySumPerCentralityAdded[icase]);
3576 
3578  (Form("hLocalRegionClusterMultiplicityPerCentralityAdded%s",caseTitle[icase].Data()),
3579  "Cluster multiplicity (not HIJING) around trigger cluster vs centrality with R=0.2",
3580  100,0,100, 200,0,200);
3581  fhLocalRegionClusterMultiplicityPerCentralityAdded[icase]->SetXTitle("Centrality");
3582  fhLocalRegionClusterMultiplicityPerCentralityAdded[icase]->SetYTitle("Multiplicity");
3583  outputContainer->Add(fhLocalRegionClusterMultiplicityPerCentralityAdded[icase]);
3584  }
3585  }
3586  }
3587  }
3588 
3590  {
3591 
3592  TString mcGenNames[] = {"","_MC_Pi0Merged","_MC_Pi0Decay","_MC_EtaDecay","_MC_PhotonOther","_MC_Electron","_MC_Other"};
3593  TString mcGenTitle[] = {"",",MC Pi0-Merged",",MC Pi0-Decay",", MC Eta-Decay",", MC Photon other sources",", MC Electron",", MC other sources"};
3594  for(Int_t igen = 0; igen < GetNCocktailGenNamesToCheck(); igen++)
3595  {
3596  TString add = "_MainGener_";
3597  if(igen==0) add = "";
3598  for(Int_t imc = 0; imc < fgkNGenTypes; imc++)
3599  {
3600  fhCleanGeneratorCluster[igen][imc] = new TH1F(Form("hCleanGeneratorCluster%s%s%s",add.Data(),GetCocktailGenNameToCheck(igen).Data(), mcGenNames[imc].Data()),
3601  Form("Number of selected clusters with contribution of %s generator%s, no overlap",
3602  GetCocktailGenNameToCheck(igen).Data(), mcGenTitle[imc].Data()),
3603  nptbins,ptmin,ptmax);
3604  fhCleanGeneratorCluster[igen][imc]->SetYTitle("#it{counts}");
3605  fhCleanGeneratorCluster[igen][imc]->SetXTitle("#it{E} (GeV)");
3606  outputContainer->Add(fhCleanGeneratorCluster[igen][imc]) ;
3607 
3608  fhCleanGeneratorClusterEPrimRecoRatio[igen][imc] = new TH2F(Form("hCleanGeneratorClusterEPrimRecoRatio%s%s%s",add.Data(),GetCocktailGenNameToCheck(igen).Data(), mcGenNames[imc].Data()),
3609  Form("#it{E}_{reco}/#it{E}_{gen} clusters with contribution of %s generator%s, no overlap",
3610  GetCocktailGenNameToCheck(igen).Data(), mcGenTitle[imc].Data()),
3611  nptbins,ptmin,ptmax,nratbins,ratmin,ratmax);
3612  fhCleanGeneratorClusterEPrimRecoRatio[igen][imc]->SetYTitle("#it{E}_{reco}/#it{E}_{gen}");
3613  fhCleanGeneratorClusterEPrimRecoRatio[igen][imc]->SetXTitle("#it{E}_{reco} (GeV)");
3614  outputContainer->Add(fhCleanGeneratorClusterEPrimRecoRatio[igen][imc]) ;
3615 
3616  fhCleanGeneratorClusterEPrimRecoDiff[igen][imc] = new TH2F(Form("hCleanGeneratorClusterEPrimRecoDiff%s%s%s",add.Data(),GetCocktailGenNameToCheck(igen).Data(), mcGenNames[imc].Data()),
3617  Form("#it{E}_{reco}-#it{E}_{gen} clusters with contribution of %s generator%s, no overlap",
3618  GetCocktailGenNameToCheck(igen).Data(), mcGenTitle[imc].Data()),
3619  nptbins,ptmin,ptmax,ndifbins,difmin,difmax);
3620  fhCleanGeneratorClusterEPrimRecoDiff[igen][imc]->SetYTitle("#it{E}_{reco}-#it{E}_{gen} (GeV)");
3621  fhCleanGeneratorClusterEPrimRecoDiff[igen][imc]->SetXTitle("#it{E}_{reco} (GeV)");
3622  outputContainer->Add(fhCleanGeneratorClusterEPrimRecoDiff[igen][imc]) ;
3623 
3624  //
3625 
3626  fhMergeGeneratorCluster[igen][imc] = new TH1F(Form("hMergeGeneratorCluster%s%s%s",add.Data(),GetCocktailGenNameToCheck(igen).Data(), mcGenNames[imc].Data()),
3627  Form("Number of selected clusters with contribution of >=2 generators, main %s%s",
3628  GetCocktailGenNameToCheck(igen).Data(), mcGenTitle[imc].Data()),
3629  nptbins,ptmin,ptmax);
3630  fhMergeGeneratorCluster[igen][imc]->SetYTitle("#it{counts}");
3631  fhMergeGeneratorCluster[igen][imc]->SetXTitle("#it{E} (GeV)");
3632  outputContainer->Add(fhMergeGeneratorCluster[igen][imc]) ;
3633 
3634  fhMergeGeneratorClusterEPrimRecoRatio[igen][imc] = new TH2F(Form("hMergeGeneratorClusterEPrimRecoRatio%s%s%s",add.Data(),GetCocktailGenNameToCheck(igen).Data(), mcGenNames[imc].Data()),
3635  Form("#it{E}_{reco}/#it{E}_{gen}clusters with contribution of >=2 generators, main %s%s",
3636  GetCocktailGenNameToCheck(igen).Data(), mcGenTitle[imc].Data()),
3637  nptbins,ptmin,ptmax,nratbins,ratmin,ratmax);
3638  fhMergeGeneratorClusterEPrimRecoRatio[igen][imc]->SetYTitle("#it{E}_{reco}/#it{E}_{gen}");
3639  fhMergeGeneratorClusterEPrimRecoRatio[igen][imc]->SetXTitle("#it{E}_{reco} (GeV)");
3640  outputContainer->Add(fhMergeGeneratorClusterEPrimRecoRatio[igen][imc]) ;
3641 
3642  fhMergeGeneratorClusterEPrimRecoDiff[igen][imc] = new TH2F(Form("hMergeGeneratorClusterEPrimRecoDiff%s%s%s",add.Data(),GetCocktailGenNameToCheck(igen).Data(), mcGenNames[imc].Data()),
3643  Form("#it{E}_{reco}-#it{E}_{gen}clusters with contribution of >=2 generators, main %s%s",
3644  GetCocktailGenNameToCheck(igen).Data(), mcGenTitle[imc].Data()),
3645  nptbins,ptmin,ptmax,ndifbins,difmin,difmax);
3646  fhMergeGeneratorClusterEPrimRecoDiff[igen][imc]->SetYTitle("#it{E}_{reco}-#it{E}_{gen} (GeV)");
3647  fhMergeGeneratorClusterEPrimRecoDiff[igen][imc]->SetXTitle("#it{E}_{reco} (GeV)");
3648  outputContainer->Add(fhMergeGeneratorClusterEPrimRecoDiff[igen][imc]) ;
3649 
3650  if(GetCocktailGenNameToCheck(igen).Contains("ijing")) continue;
3651 
3652  //
3653 
3654  fhMergeGeneratorClusterNotHijingBkg[igen][imc] = new TH1F(Form("hMergeGeneratorClusterNotHijingBkg%s%s%s",add.Data(),GetCocktailGenNameToCheck(igen).Data(), mcGenNames[imc].Data()),
3655  Form("Number of selected clusters with contribution of >=2 generators, , none is HIJING, main %s%s",
3656  GetCocktailGenNameToCheck(igen).Data(), mcGenTitle[imc].Data()),
3657  nptbins,ptmin,ptmax);
3658  fhMergeGeneratorClusterNotHijingBkg[igen][imc]->SetYTitle("#it{counts}");
3659  fhMergeGeneratorClusterNotHijingBkg[igen][imc]->SetXTitle("#it{E} (GeV)");
3660  outputContainer->Add(fhMergeGeneratorClusterNotHijingBkg[igen][imc]) ;
3661 
3662  fhMergeGeneratorClusterNotHijingBkgEPrimRecoRatio[igen][imc] = new TH2F(Form("hMergeGeneratorClusterNotHijingBkgEPrimRecoRatio%s%s%s",add.Data(),GetCocktailGenNameToCheck(igen).Data(), mcGenNames[imc].Data()),
3663  Form("#it{E}_{reco}/#it{E}_{gen} clusters with contribution of >=2 generators, , none is HIJING, main %s%s",
3664  GetCocktailGenNameToCheck(igen).Data(), mcGenTitle[imc].Data()),
3665  nptbins,ptmin,ptmax,nratbins,ratmin,ratmax);
3666  fhMergeGeneratorClusterNotHijingBkgEPrimRecoRatio[igen][imc]->SetYTitle("#it{E}_{reco}/#it{E}_{gen}");
3667  fhMergeGeneratorClusterNotHijingBkgEPrimRecoRatio[igen][imc]->SetXTitle("#it{E}_{reco} (GeV)");
3668  outputContainer->Add(fhMergeGeneratorClusterNotHijingBkgEPrimRecoRatio[igen][imc]) ;
3669 
3670  fhMergeGeneratorClusterNotHijingBkgEPrimRecoDiff[igen][imc] = new TH2F(Form("hMergeGeneratorClusterNotHijingBkgEPrimRecoDiff%s%s%s",add.Data(),GetCocktailGenNameToCheck(igen).Data(), mcGenNames[imc].Data()),
3671  Form("#it{E}_{reco}-#it{E}_{gen} clusters with contribution of >=2 generators, , none is HIJING, main %s%s",
3672  GetCocktailGenNameToCheck(igen).Data(), mcGenTitle[imc].Data()),
3673  nptbins,ptmin,ptmax,ndifbins,difmin,difmax);
3674  fhMergeGeneratorClusterNotHijingBkgEPrimRecoDiff[igen][imc]->SetYTitle("#it{E}_{reco}-#it{E}_{gen} (GeV)");
3675  fhMergeGeneratorClusterNotHijingBkgEPrimRecoDiff[igen][imc]->SetXTitle("#it{E}_{reco} (GeV)");
3676  outputContainer->Add(fhMergeGeneratorClusterNotHijingBkgEPrimRecoDiff[igen][imc]) ;
3677 
3678  //
3679 
3680  fhMergeGeneratorClusterHijingAndOtherBkg[igen][imc] = new TH1F(Form("hMergeGeneratorClusterHijingAndOtherBkg%s%s%s",add.Data(),GetCocktailGenNameToCheck(igen).Data(), mcGenNames[imc].Data()),
3681  Form("Number of selected clusters with contribution of >=3 generators, none is HIJING, main %s%s",
3682  GetCocktailGenNameToCheck(igen).Data(), mcGenTitle[imc].Data()),
3683  nptbins,ptmin,ptmax);
3684  fhMergeGeneratorClusterHijingAndOtherBkg[igen][imc]->SetYTitle("#it{counts}");
3685  fhMergeGeneratorClusterHijingAndOtherBkg[igen][imc]->SetXTitle("#it{E} (GeV)");
3686  outputContainer->Add(fhMergeGeneratorClusterHijingAndOtherBkg[igen][imc]) ;
3687 
3688 
3689  fhMergeGeneratorClusterHijingAndOtherBkgEPrimRecoRatio[igen][imc] = new TH2F(Form("hMergeGeneratorClusterHijingAndOtherBkgEPrimRecoRatio%s%s%s",add.Data(),GetCocktailGenNameToCheck(igen).Data(), mcGenNames[imc].Data()),
3690  Form("#it{E}_{reco}/#it{E}_{gen} clusters with contribution of >=3 generators, none is HIJING, main %s%s",
3691  GetCocktailGenNameToCheck(igen).Data(), mcGenTitle[imc].Data()),
3692  nptbins,ptmin,ptmax,nratbins,ratmin,ratmax);
3693  fhMergeGeneratorClusterHijingAndOtherBkgEPrimRecoRatio[igen][imc]->SetYTitle("#it{E}_{reco}/#it{E}_{gen}");
3694  fhMergeGeneratorClusterHijingAndOtherBkgEPrimRecoRatio[igen][imc]->SetXTitle("#it{E}_{reco} (GeV)");
3695  outputContainer->Add(fhMergeGeneratorClusterHijingAndOtherBkgEPrimRecoRatio[igen][imc]) ;
3696 
3697 
3698  fhMergeGeneratorClusterHijingAndOtherBkgEPrimRecoDiff[igen][imc] = new TH2F(Form("hMergeGeneratorClusterHijingAndOtherBkgEPrimRecoDiff%s%s%s",add.Data(),GetCocktailGenNameToCheck(igen).Data(), mcGenNames[imc].Data()),
3699  Form("#it{E}_{reco}-#it{E}_{gen} clusters with contribution of >=3 generators, none is HIJING, main %s%s",
3700  GetCocktailGenNameToCheck(igen).Data(), mcGenTitle[imc].Data()),
3701  nptbins,ptmin,ptmax,ndifbins,difmin,difmax);
3702  fhMergeGeneratorClusterHijingAndOtherBkgEPrimRecoDiff[igen][imc]->SetYTitle("#it{E}_{reco}-#it{E}_{gen} (GeV)");
3703  fhMergeGeneratorClusterHijingAndOtherBkgEPrimRecoDiff[igen][imc]->SetXTitle("#it{E}_{reco} (GeV)");
3704  outputContainer->Add(fhMergeGeneratorClusterHijingAndOtherBkgEPrimRecoDiff[igen][imc]) ;
3705 
3706  //
3707 
3708  fhMergeGeneratorClusterHijingBkg[igen][imc] = new TH1F(Form("hMergeGeneratorClusterHijingBkg%s%s%s",add.Data(),GetCocktailGenNameToCheck(igen).Data(), mcGenNames[imc].Data()),
3709  Form("Number of selected clusters with contribution of >=3 generators, none is HIJING, main %s%s",
3710  GetCocktailGenNameToCheck(igen).Data(), mcGenTitle[imc].Data()),
3711  nptbins,ptmin,ptmax);
3712  fhMergeGeneratorClusterHijingBkg[igen][imc]->SetYTitle("#it{counts}");
3713  fhMergeGeneratorClusterHijingBkg[igen][imc]->SetXTitle("#it{E} (GeV)");
3714  outputContainer->Add(fhMergeGeneratorClusterHijingBkg[igen][imc]) ;
3715 
3716  fhMergeGeneratorClusterHijingBkgEPrimRecoRatio[igen][imc] = new TH2F(Form("hMergeGeneratorClusterHijingBkgEPrimRecoRatio%s%s%s",add.Data(),GetCocktailGenNameToCheck(igen).Data(), mcGenNames[imc].Data()),
3717  Form("#it{E}_{reco}/#it{E}_{gen} clusters with contribution of >=3 generators, none is HIJING, main %s%s",
3718  GetCocktailGenNameToCheck(igen).Data(), mcGenTitle[imc].Data()),
3719  nptbins,ptmin,ptmax,nratbins,ratmin,ratmax);
3720  fhMergeGeneratorClusterHijingBkgEPrimRecoRatio[igen][imc]->SetYTitle("#it{E}_{reco}/#it{E}_{gen}");
3721  fhMergeGeneratorClusterHijingBkgEPrimRecoRatio[igen][imc]->SetXTitle("#it{E}_{reco} (GeV)");
3722  outputContainer->Add(fhMergeGeneratorClusterHijingBkgEPrimRecoRatio[igen][imc]) ;
3723 
3724  fhMergeGeneratorClusterHijingBkgEPrimRecoDiff[igen][imc] = new TH2F(Form("hMergeGeneratorClusterHijingBkgEPrimRecoDiff%s%s%s",add.Data(),GetCocktailGenNameToCheck(igen).Data(), mcGenNames[imc].Data()),
3725  Form("#it{E}_{reco}-#it{E}_{gen} clusters with contribution of >=3 generators, none is HIJING, main %s%s",
3726  GetCocktailGenNameToCheck(igen).Data(), mcGenTitle[imc].Data()),
3727  nptbins,ptmin,ptmax,ndifbins,difmin,difmax);
3728  fhMergeGeneratorClusterHijingBkgEPrimRecoDiff[igen][imc]->SetYTitle("#it{E}_{reco}-#it{E}_{gen} (GeV)");
3729  fhMergeGeneratorClusterHijingBkgEPrimRecoDiff[igen][imc]->SetXTitle("#it{E}_{reco} (GeV)");
3730  outputContainer->Add(fhMergeGeneratorClusterHijingBkgEPrimRecoDiff[igen][imc]) ;
3731  }
3732  }
3733  }
3734 
3735  return outputContainer ;
3736 }
3737 
3738 //_______________________
3740 //_______________________
3742 {
3743  if ( GetCalorimeter() == kPHOS && !GetReader()->IsPHOSSwitchedOn() && NewOutputAOD() )
3744  AliFatal("!!STOP: You want to use PHOS in analysis but it is not read!! \n!!Check the configuration file!!");
3745  else if( GetCalorimeter() == kEMCAL && !GetReader()->IsEMCALSwitchedOn() && NewOutputAOD() )
3746  AliFatal("!!STOP: You want to use EMCAL in analysis but it is not read!! \n!!Check the configuration file!!");
3747 
3748  // Trick to select primary photon particles in the analysis of pure MC input
3749  if(GetReader()->GetDataType() == AliCaloTrackReader::kMC) GetCaloPID()->SwitchOnBayesian();
3750 }
3751 
3752 //_________________________________
3754 //_________________________________
3756 {
3757  AddToHistogramsName("AnaPhoton_");
3758 
3759  fMinDist = 2.;
3760  fMinDist2 = 4.;
3761  fMinDist3 = 5.;
3762 
3763  fTimeCutMin =-1000000;
3764  fTimeCutMax = 1000000;
3765  fNCellsCut = 0;
3766 
3767  fRejectTrackMatch = kTRUE ;
3768 
3769  fNEBinCuts = 14;
3770  fEBinCuts[0] = 0.; fEBinCuts[1] = 0.3; fEBinCuts[2] = 0.5;
3771  fEBinCuts[3] = 1.; fEBinCuts[4] = 2. ; fEBinCuts[5] = 3. ;
3772  fEBinCuts[6] = 4.; fEBinCuts[7] = 5. ; fEBinCuts[8] = 7. ;
3773  fEBinCuts[9] = 9.; fEBinCuts[10]= 12.; fEBinCuts[11]= 15.;
3774  fEBinCuts[12]= 20.; fEBinCuts[13]= 50.; fEBinCuts[14]= 100.;
3775  for(Int_t i = fNEBinCuts; i < 15; i++) fEBinCuts[i] = 1000.;
3776 }
3777 
3778 //_______________________________________
3782 //_______________________________________
3784 {
3785  // Get the vertex
3786  Double_t v[3] = {0,0,0}; //vertex ;
3787  GetReader()->GetVertex(v);
3788 
3789  // Select the Calorimeter of the photon
3790  TObjArray * pl = 0x0;
3791  AliVCaloCells* cells = 0;
3792  if (GetCalorimeter() == kPHOS )
3793  {
3794  pl = GetPHOSClusters();
3795  cells = GetPHOSCells();
3796  }
3797  else if (GetCalorimeter() == kEMCAL)
3798  {
3799  pl = GetEMCALClusters();
3800  cells = GetEMCALCells();
3801  }
3802 
3803  if(!pl)
3804  {
3805  AliWarning(Form("TObjArray with %s clusters is NULL!",GetCalorimeterString().Data()));
3806  return;
3807  }
3808 
3809  // Loop on raw clusters before filtering in the reader and fill control histogram
3810  if((GetReader()->GetEMCALClusterListName()=="" && GetCalorimeter()==kEMCAL) || GetCalorimeter()==kPHOS)
3811  {
3812  for(Int_t iclus = 0; iclus < GetReader()->GetInputEvent()->GetNumberOfCaloClusters(); iclus++ )
3813  {
3814  AliVCluster * clus = GetReader()->GetInputEvent()->GetCaloCluster(iclus);
3815 
3816  if (GetCalorimeter() == kPHOS && clus->IsPHOS() && clus->E() > GetReader()->GetPHOSPtMin() )
3817  {
3818  clus->GetMomentum(fMomentum,GetVertex(0)) ;
3819 
3820  fhClusterCutsE [0]->Fill(fMomentum.E() , GetEventWeight());
3821  fhClusterCutsPt[0]->Fill(fMomentum.Pt(), GetEventWeight());
3822  }
3823  else if(GetCalorimeter() == kEMCAL && clus->IsEMCAL() && clus->E() > GetReader()->GetEMCALPtMin())
3824  {
3825  clus->GetMomentum(fMomentum,GetVertex(0)) ;
3826 
3827  fhClusterCutsE [0]->Fill(fMomentum.E(), GetEventWeight());
3828  fhClusterCutsPt[0]->Fill(fMomentum.Pt(), GetEventWeight());
3829  }
3830  }
3831  }
3832  else
3833  { // reclusterized
3834  TClonesArray * clusterList = 0;
3835 
3836  if(GetReader()->GetInputEvent()->FindListObject(GetReader()->GetEMCALClusterListName()))
3837  clusterList = dynamic_cast<TClonesArray*> (GetReader()->GetInputEvent()->FindListObject(GetReader()->GetEMCALClusterListName()));
3838  else if(GetReader()->GetOutputEvent())
3839  clusterList = dynamic_cast<TClonesArray*> (GetReader()->GetOutputEvent()->FindListObject(GetReader()->GetEMCALClusterListName()));
3840 
3841  if(clusterList)
3842  {
3843  Int_t nclusters = clusterList->GetEntriesFast();
3844  for (Int_t iclus = 0; iclus < nclusters; iclus++)
3845  {
3846  AliVCluster * clus = dynamic_cast<AliVCluster*> (clusterList->At(iclus));
3847 
3848  if(clus && clus->E() > GetReader()->GetEMCALPtMin())
3849  {
3850  clus->GetMomentum(fMomentum,GetVertex(0)) ;
3851 
3852  fhClusterCutsE [0]->Fill(clus->E() , GetEventWeight());
3853  fhClusterCutsPt[0]->Fill(fMomentum.Pt(), GetEventWeight());
3854  }
3855  }
3856  }
3857  }
3858 
3859  // Init arrays, variables, get number of clusters
3860  Int_t nCaloClusters = pl->GetEntriesFast();
3861 
3862  AliDebug(1,Form("Input %s cluster entries %d", GetCalorimeterString().Data(), nCaloClusters));
3863 
3864  //----------------------------------------------------
3865  // Fill AOD with PHOS/EMCAL AliAODPWG4Particle objects
3866  //----------------------------------------------------
3867  // Loop on clusters
3868  for(Int_t icalo = 0; icalo < nCaloClusters; icalo++)
3869  {
3870  AliVCluster * calo = (AliVCluster*) (pl->At(icalo));
3871  //printf("calo %d, %f\n",icalo,calo->E());
3872 
3873  //Get the index where the cluster comes, to retrieve the corresponding vertex
3874  Int_t evtIndex = 0 ;
3875  if (GetMixedEvent())
3876  {
3877  evtIndex=GetMixedEvent()->EventIndexForCaloCluster(calo->GetID()) ;
3878  //Get the vertex and check it is not too large in z
3879  if(TMath::Abs(GetVertex(evtIndex)[2])> GetZvertexCut()) continue;
3880  }
3881 
3882  //Cluster selection, not charged, with photon id and in fiducial cut
3884  {
3885  calo->GetMomentum(fMomentum,GetVertex(evtIndex)) ;
3886  }//Assume that come from vertex in straight line
3887  else
3888  {
3889  Double_t vertex[]={0,0,0};
3890  calo->GetMomentum(fMomentum,vertex) ;
3891  }
3892 
3893  //-----------------------------
3894  // Cluster selection
3895  //-----------------------------
3896  Int_t nMaxima = GetCaloUtils()->GetNumberOfLocalMaxima(calo, cells); // NLM
3897  if(!ClusterSelected(calo,nMaxima)) continue;
3898 
3899  //----------------------------
3900  // Create AOD for analysis
3901  //----------------------------
3902  AliAODPWG4Particle aodph = AliAODPWG4Particle(fMomentum);
3903 
3904  //...............................................
3905  // Set the indeces of the original caloclusters (MC, ID), and calorimeter
3906  Int_t label = calo->GetLabel();
3907  aodph.SetLabel(label);
3908  aodph.SetCaloLabel(calo->GetID(),-1);
3909  aodph.SetDetectorTag(GetCalorimeter());
3910  //printf("Index %d, Id %d, iaod %d\n",icalo, calo->GetID(),GetOutputAODBranch()->GetEntriesFast());
3911 
3912  //...............................................
3913  // Set bad channel distance bit
3914  Double_t distBad=calo->GetDistanceToBadChannel() ; //Distance to bad channel
3915  if (distBad > fMinDist3) aodph.SetDistToBad(2) ;
3916  else if(distBad > fMinDist2) aodph.SetDistToBad(1) ;
3917  else aodph.SetDistToBad(0) ;
3918  //printf("DistBad %f Bit %d\n",distBad, aodph.DistToBad());
3919 
3920  //-------------------------------------
3921  // Play with the MC stack if available
3922  //-------------------------------------
3923 
3924  // Check origin of the candidates
3925  Int_t tag = -1;
3926 
3927  if(IsDataMC())
3928  {
3929  tag = GetMCAnalysisUtils()->CheckOrigin(calo->GetLabels(),calo->GetNLabels(),GetReader(),GetCalorimeter());
3930  aodph.SetTag(tag);
3931 
3932  AliDebug(1,Form("Origin of candidate, bit map %d",aodph.GetTag()));
3933  }//Work with stack also
3934 
3935  //--------------------------------------------------------
3936  // Fill some shower shape histograms before PID is applied
3937  //--------------------------------------------------------
3938 
3939  Float_t maxCellFraction = 0;
3940  Int_t absIdMax = GetCaloUtils()->GetMaxEnergyCell(cells, calo, maxCellFraction);
3941  if( absIdMax < 0 ) AliFatal("Wrong absID");
3942 
3943  Int_t largeTimeInCellCluster = kFALSE;
3944  FillShowerShapeHistograms(calo,tag,nMaxima,maxCellFraction,largeTimeInCellCluster);
3945  aodph.SetFiducialArea(largeTimeInCellCluster); // Temporary use of this container, FIXME
3946  //if(largeTimeInCellCluster > 1) printf("Set n cells large time %d, pt %2.2f\n",aodph.GetFiducialArea(),aodph.Pt());
3947 
3948  aodph.SetM02(calo->GetM02());
3949  aodph.SetM20(calo->GetM20());
3950  aodph.SetNLM(nMaxima);
3951 
3952  Float_t time = calo->GetTOF()*1e9;
3953  if(time > 400) time-=fConstantTimeShift; // in case of clusterizer running before (calibrate clusters not cells)
3954  aodph.SetTime(time);
3955 
3956  aodph.SetNCells(calo->GetNCells());
3957  Int_t nSM = GetModuleNumber(calo);
3958  aodph.SetSModNumber(nSM);
3959 
3960  Float_t en = fMomentum.E ();
3961  Float_t pt = fMomentum.Pt();
3962  Float_t eta = fMomentum.Eta();
3963  Float_t phi = GetPhi(fMomentum.Phi());
3964  Int_t ebin = -1;
3965  for(Int_t ie = 0; ie < fNEBinCuts; ie++)
3966  {
3967  if( en >= fEBinCuts[ie] && en < fEBinCuts[ie+1] ) ebin = ie;
3968  }
3969 
3970  Int_t icolAbs = -1, irowAbs = -1;
3972  {
3973  Float_t maxCellFraction = 0;
3974  Int_t absIdMax = GetCaloUtils()->GetMaxEnergyCell(cells,calo,maxCellFraction);
3975 
3976  Int_t icol = -1, irow = -1, iRCU = -1;
3977  GetModuleNumberCellIndexesAbsCaloMap(absIdMax,GetCalorimeter(), icol, irow, iRCU, icolAbs, irowAbs);
3978 
3979  if(ebin>=0 && ebin < fNEBinCuts)
3980  {
3981  fhEBinClusterEtaPhi[ebin]->Fill(eta,phi,GetEventWeight()) ;
3982 
3983  fhEBinClusterColRow[ebin]->Fill(icolAbs,irowAbs,GetEventWeight()) ;
3984  }
3985  }
3986 
3987  //-------------------------------------
3988  // PID selection or bit setting
3989  //-------------------------------------
3990 
3991  //...............................................
3992  // Data, PID check on
3993  if(IsCaloPIDOn())
3994  {
3995  // Get most probable PID, 2 options check bayesian PID weights or redo PID
3996  // By default, redo PID
3997 
3998  aodph.SetIdentifiedParticleType(GetCaloPID()->GetIdentifiedParticleType(calo));
3999 
4000  AliDebug(1,Form("PDG of identified particle %d",aodph.GetIdentifiedParticleType()));
4001 
4002  //If cluster does not pass pid, not photon, skip it.
4003  if(aodph.GetIdentifiedParticleType() != AliCaloPID::kPhoton) continue ;
4004  }
4005 
4006  //...............................................
4007  // Data, PID check off
4008  else
4009  {
4010  // Set PID bits for later selection (AliAnaPi0 for example)
4011  // GetIdentifiedParticleType already called in SetPIDBits.
4012 
4013  GetCaloPID()->SetPIDBits(calo,&aodph, GetCaloUtils(),GetReader()->GetInputEvent());
4014 
4015  AliDebug(1,"PID Bits set");
4016  }
4017 
4018  AliDebug(1,Form("Photon selection cuts passed: pT %3.2f, pdg %d",aodph.Pt(),aodph.GetIdentifiedParticleType()));
4019 
4020  fhClusterCutsE [9]->Fill(en, GetEventWeight());
4021  fhClusterCutsPt[9]->Fill(pt, GetEventWeight());
4022 
4023  //
4024  // Check local cluster activity around the current cluster
4025  //
4026  if(fStudyActivityNearCluster && en > 1.5) // 1.5 GeV cut used on Pb-Pb analysis
4027  ActivityNearCluster(icalo,en,eta,phi,pl);
4028 
4029  //
4030  // Check if other generators contributed to the cluster
4031  //
4034 
4035 
4037  {
4038  if(ebin>=0 && ebin < fNEBinCuts)
4039  {
4040  fhEBinClusterEtaPhiPID[ebin]->Fill(eta,phi,GetEventWeight()) ;
4041 
4042  fhEBinClusterColRowPID[ebin]->Fill(icolAbs,irowAbs,GetEventWeight()) ;
4043  }
4044  }
4045 
4046  if(nSM < GetCaloUtils()->GetNumberOfSuperModulesUsed() && nSM >=0)
4047  {
4048  fhEPhotonSM ->Fill(en, nSM, GetEventWeight());
4049  fhPtPhotonSM->Fill(pt, nSM, GetEventWeight());
4050  }
4051 
4052  fhNLocMax->Fill(calo->E(),nMaxima);
4053 
4054  // Few more control histograms for selected clusters
4055  fhMaxCellDiffClusterE->Fill(en, maxCellFraction , GetEventWeight());
4056  fhNCellsE ->Fill(en, calo->GetNCells() , GetEventWeight());
4057  fhTimePt ->Fill(pt, time , GetEventWeight());
4058 
4059  if(cells)
4060  {
4061  for(Int_t icell = 0; icell < calo->GetNCells(); icell++)
4062  fhCellsE->Fill(en, cells->GetCellAmplitude(calo->GetCellsAbsId()[icell]), GetEventWeight());
4063  }
4064 
4065  // Matching after cuts
4067 
4068  // Fill histograms to undertand pile-up before other cuts applied
4069  // Remember to relax time cuts in the reader
4070  if( IsPileUpAnalysisOn() ) FillPileUpHistograms(calo,cells, absIdMax);
4071 
4072  // Add AOD with photon object to aod branch
4073  AddAODParticle(aodph);
4074  }// loop
4075 
4076  AliDebug(1,Form("End fill AODs, with %d entries",GetOutputAODBranch()->GetEntriesFast()));
4077 }
4078 
4079 //______________________________________________
4080 // Fill histograms with selected clusters/output AOD particles.
4081 //______________________________________________
4083 {
4084  // In case of simulated data, fill acceptance histograms
4086 
4087  // Get vertex
4088  Double_t v[3] = {0,0,0}; //vertex ;
4089  GetReader()->GetVertex(v);
4090  //fhVertex->Fill(v[0],v[1],v[2]);
4091  if(TMath::Abs(v[2]) > GetZvertexCut()) return ; // done elsewhere for Single Event analysis, but there for mixed event
4092 
4093  //----------------------------------
4094  // Loop on stored AOD photons
4095  Int_t naod = GetOutputAODBranch()->GetEntriesFast();
4096  AliDebug(1,Form("AOD branch entries %d", naod));
4097 
4098  Float_t cen = GetEventCentrality();
4099  // printf("++++++++++ GetEventCentrality() %f\n",cen);
4100 
4101  Float_t ep = GetEventPlaneAngle();
4102 
4103  for(Int_t iaod = 0; iaod < naod ; iaod++)
4104  {
4105  AliAODPWG4Particle* ph = (AliAODPWG4Particle*) (GetOutputAODBranch()->At(iaod));
4106  Int_t pdg = ph->GetIdentifiedParticleType();
4107 
4108  AliDebug(2,Form("PDG %d, MC TAG %d, Calorimeter <%d>",ph->GetIdentifiedParticleType(),ph->GetTag(), ph->GetDetectorTag())) ;
4109 
4110  // If PID used, fill histos with photons in Calorimeter GetCalorimeter()
4111  if(IsCaloPIDOn() && pdg != AliCaloPID::kPhoton) continue;
4112 
4113  if(((Int_t) ph->GetDetectorTag()) != GetCalorimeter()) continue;
4114 
4115  AliDebug(2,Form("ID Photon: pt %f, phi %f, eta %f", ph->Pt(),ph->Phi(),ph->Eta())) ;
4116 
4117  //................................
4118  //Fill photon histograms
4119  Float_t ptcluster = ph->Pt();
4120  Float_t phicluster = ph->Phi();
4121  Float_t etacluster = ph->Eta();
4122  Float_t ecluster = ph->E();
4123 
4124  fhEPhoton ->Fill(ecluster , GetEventWeight());
4125  fhPtPhoton ->Fill(ptcluster, GetEventWeight());
4126 
4127  fhPhiPhoton ->Fill(ptcluster, phicluster, GetEventWeight());
4128  fhEtaPhoton ->Fill(ptcluster, etacluster, GetEventWeight());
4129 
4130  if (ecluster > 0.5) fhEtaPhiPhoton ->Fill(etacluster, phicluster, GetEventWeight());
4131  else if(GetMinPt() < 0.5) fhEtaPhi05Photon->Fill(etacluster, phicluster, GetEventWeight());
4132 
4134  {
4135  fhPtCentralityPhoton ->Fill(ptcluster,cen, GetEventWeight()) ;
4136  fhPtEventPlanePhoton ->Fill(ptcluster,ep , GetEventWeight()) ;
4137  }
4138 
4139 // Comment this part, not needed but in case to know how to do it in the future
4140 // // Get original cluster, to recover some information
4141 // AliVCaloCells* cells = 0;
4142 // TObjArray * clusters = 0;
4143 // if(GetCalorimeter() == kEMCAL)
4144 // {
4145 // cells = GetEMCALCells();
4146 // clusters = GetEMCALClusters();
4147 // }
4148 // else
4149 // {
4150 // cells = GetPHOSCells();
4151 // clusters = GetPHOSClusters();
4152 // }
4153 //
4154 // Int_t iclus = -1;
4155 // AliVCluster *cluster = FindCluster(clusters,ph->GetCaloLabel(0),iclus);
4156 // if(cluster)
4157 
4158  //.......................................
4159  // Play with the MC data if available
4160  if(IsDataMC())
4161  {
4162  //....................................................................
4163  // Access MC information in stack if requested, check that it exists.
4164  Int_t label = ph->GetLabel();
4165 
4166  if(label < 0)
4167  {
4168  AliDebug(1,Form("*** bad label ***: label %d", label));
4169  continue;
4170  }
4171 
4172  Float_t eprim = 0;
4173  Float_t ptprim = 0;
4174  Bool_t ok = kFALSE;
4175  Int_t pdg = 0, status = 0, momLabel = -1;
4176 
4177  //fPrimaryMom = GetMCAnalysisUtils()->GetMother(label,GetReader(),ok);
4178  fPrimaryMom = GetMCAnalysisUtils()->GetMother(label,GetReader(), pdg, status, ok, momLabel);
4179 
4180  if(ok)
4181  {
4182  eprim = fPrimaryMom.Energy();
4183  ptprim = fPrimaryMom.Pt();
4184  }
4185 
4186  Int_t tag =ph->GetTag();
4187  Int_t mcParticleTag = -1;
4189  {
4190  fhMCE [kmcPhoton] ->Fill(ecluster , GetEventWeight());
4191  fhMCPt [kmcPhoton] ->Fill(ptcluster, GetEventWeight());
4192 
4193  fhMCPhi[kmcPhoton] ->Fill(ecluster,phicluster, GetEventWeight());
4194  fhMCEta[kmcPhoton] ->Fill(ecluster,etacluster, GetEventWeight());
4195 
4196  fhMC2E [kmcPhoton] ->Fill(ecluster , eprim , GetEventWeight());
4197  fhMC2Pt [kmcPhoton] ->Fill(ptcluster, ptprim, GetEventWeight());
4198 
4199  fhMCDeltaE [kmcPhoton] ->Fill(ecluster , eprim-ecluster , GetEventWeight());
4200  fhMCDeltaPt[kmcPhoton] ->Fill(ptcluster, ptprim-ptcluster, GetEventWeight());
4201 
4202  if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCConversion) &&
4204  {
4205  fhMCE [kmcConversion] ->Fill(ecluster , GetEventWeight());
4206  fhMCPt [kmcConversion] ->Fill(ptcluster, GetEventWeight());
4207 
4208  fhMCPhi[kmcConversion] ->Fill(ecluster, phicluster, GetEventWeight());
4209  fhMCEta[kmcConversion] ->Fill(ecluster, etacluster, GetEventWeight());
4210 
4211  fhMC2E [kmcConversion] ->Fill(ecluster , eprim , GetEventWeight());
4212  fhMC2Pt [kmcConversion] ->Fill(ptcluster, ptprim, GetEventWeight());
4213 
4214  fhMCDeltaE [kmcConversion] ->Fill(ecluster , eprim-ecluster , GetEventWeight());
4215  fhMCDeltaPt[kmcConversion] ->Fill(ptcluster, ptprim-ptcluster, GetEventWeight());
4216 
4217  Int_t pdgD = 0, statusD = 0, daugLabel = -1;
4218  Bool_t okD = kFALSE;
4219 
4220  //fMomentum =
4221  GetMCAnalysisUtils()->GetDaughter(0,momLabel,GetReader(),pdgD, statusD, okD, daugLabel, fProdVertex);
4222 
4223  if(okD)
4224  {
4225  Float_t prodR = TMath::Sqrt(fProdVertex.X()*fProdVertex.X()+fProdVertex.Y()*fProdVertex.Y());
4226 
4227  //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,
4228  // momLabel, label,daugLabel,prodR);
4229 
4231  {
4232 
4233  fhMCConversionVertex->Fill(ptcluster,prodR,GetEventWeight());
4234 
4235  if(GetCalorimeter() == kEMCAL && GetFirstSMCoveredByTRD() >= 0 && ph->GetSModNumber() >= GetFirstSMCoveredByTRD() )
4236  fhMCConversionVertexTRD->Fill(ptcluster,prodR,GetEventWeight());
4237  }
4238 
4239  if ( fFillSSHistograms )
4240  {
4241  Float_t m02 = ph->GetM02();
4242  Float_t m20 = ph->GetM20();
4243 
4244  // conversion vertex vs shower shape
4246  {
4247  Int_t convR = -1;
4248  if ( prodR < 75. ) convR = 0;
4249  else if ( prodR < 275. ) convR = 1;
4250  else if ( prodR < 375. ) convR = 2;
4251  else if ( prodR < 400. ) convR = 3;
4252  else if ( prodR < 430. ) convR = 4;
4253  else convR = 5;
4254 
4255  if ( convR >= 0 )
4256  {
4257  fhMCConversionLambda0Rcut[convR]->Fill(ptcluster,m02,GetEventWeight());
4258  fhMCConversionLambda1Rcut[convR]->Fill(ptcluster,m20,GetEventWeight());
4259 
4260  if ( GetCalorimeter() == kEMCAL && GetFirstSMCoveredByTRD() >= 0 && ph->GetSModNumber() >= GetFirstSMCoveredByTRD() )
4261  {
4262  fhMCConversionLambda0RcutTRD[convR]->Fill(ptcluster,m02,GetEventWeight());
4263  fhMCConversionLambda1RcutTRD[convR]->Fill(ptcluster,m20,GetEventWeight());
4264  }
4265  // //
4266  // // EMCAL SM regions
4267  // //
4268  // if ( GetCalorimeter() == kEMCAL && fFillEMCALRegionSSHistograms )
4269  // {
4270  // // Get original cluster, needed to feed the subregion selection method
4271  //
4272  // Int_t iclus = -1;
4273  // AliVCluster *cluster = FindCluster(GetEMCALClusters(),ph->GetCaloLabel(0),iclus);
4274  //
4275  // Int_t etaRegion = -1, phiRegion = -1;
4276  //
4277  // if ( cluster ) GetCaloUtils()->GetEMCALSubregion(cluster,GetReader()->GetEMCALCells(),etaRegion,phiRegion);
4278  //
4279  // if( etaRegion >= 0 && etaRegion < 4 && phiRegion >=0 && phiRegion < 3 )
4280  // {
4281  // fhLam0EMCALRegionMCConvRcut[etaRegion][phiRegion][convR]->Fill(ptcluster,m02, GetEventWeight());
4282  //
4283  // if ( GetFirstSMCoveredByTRD() >= 0 && ph->GetSModNumber() >= GetFirstSMCoveredByTRD() )
4284  // fhLam0EMCALRegionTRDMCConvRcut[etaRegion][phiRegion][convR]->Fill(ptcluster, m02, GetEventWeight());
4285  //
4286  // } // region found
4287  // } // check region
4288  } // conv region
4289  } // check conv region
4290 
4291  } // fill Sh Sh histograms
4292  } // okD
4293  } // conversion
4294 
4295  if ( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPrompt) )
4296  {
4297  mcParticleTag = kmcPrompt;
4298  }
4299  else if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCFragmentation) )
4300  {
4301  mcParticleTag = kmcFragmentation;
4302  }
4303  else if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCISR) )
4304  {
4305  mcParticleTag = kmcISR;
4306  }
4307  else if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0) )
4308  {
4309  mcParticleTag = kmcPi0;
4310  }
4311  else if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCEta) )
4312  {
4313  mcParticleTag = kmcEta;
4314  }
4315  else if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0Decay) )
4316  {
4317  mcParticleTag = kmcPi0Decay;
4318  }
4319  else if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCEtaDecay) )
4320  {
4321  mcParticleTag = kmcEtaDecay;
4322  }
4323  else
4324  {
4325  mcParticleTag = kmcOtherDecay;
4326  }
4327  }
4328  else if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCAntiNeutron) )
4329  {
4330  mcParticleTag = kmcAntiNeutron;
4331  }
4332  else if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCAntiProton) )
4333  {
4334  mcParticleTag = kmcAntiProton;
4335  }
4336  else if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCNeutron) )
4337  {
4338  mcParticleTag = kmcNeutron;
4339  }
4340  else if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCProton) )
4341  {
4342  mcParticleTag = kmcProton;
4343  }
4344  else if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPion) )
4345  {
4346  mcParticleTag = kmcChPion;
4347  }
4348  else if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCElectron) )
4349  {
4350  mcParticleTag = kmcElectron;
4351  }
4352  else if( fhMCE[kmcOther] )
4353  {
4354  mcParticleTag = kmcOther;
4355 
4356  // printf(" AliAnaPhoton::MakeAnalysisFillHistograms() - Label %d, pT %2.3f Unknown, bits set: ",
4357  // ph->GetLabel(),ph->Pt());
4358  // for(Int_t i = 0; i < 20; i++) {
4359  // if(GetMCAnalysisUtils()->CheckTagBit(tag,i)) printf(" %d, ",i);
4360  // }
4361  // printf("\n");
4362  }
4363 
4364  if(mcParticleTag >= 0 && fhMCE[mcParticleTag])
4365  {
4366  fhMCE [mcParticleTag]->Fill(ecluster , GetEventWeight());
4367  fhMCPt [mcParticleTag]->Fill(ptcluster, GetEventWeight());
4368 
4369  fhMCPhi [mcParticleTag]->Fill(ecluster, phicluster, GetEventWeight());
4370  fhMCEta [mcParticleTag]->Fill(ecluster, etacluster, GetEventWeight());
4371 
4372  fhMC2E [mcParticleTag]->Fill(ecluster , eprim , GetEventWeight());
4373  fhMC2Pt [mcParticleTag]->Fill(ptcluster, ptprim, GetEventWeight());
4374 
4375  fhMCDeltaE [mcParticleTag]->Fill(ecluster , eprim-ecluster , GetEventWeight());
4376  fhMCDeltaPt[mcParticleTag]->Fill(ptcluster, ptprim-ptcluster, GetEventWeight());
4377  }
4378  }// Histograms with MC
4379  }// aod loop
4380 }
4381 
4382 //__________________________________________________
4384 //__________________________________________________
4385 void AliAnaPhoton::Print(const Option_t * opt) const
4386 {
4387  if(! opt)
4388  return;
4389 
4390  printf("**** Print %s %s ****\n", GetName(), GetTitle() ) ;
4392 
4393  printf("Calorimeter = %s\n", GetCalorimeterString().Data()) ;
4394  printf("Min Distance to Bad Channel = %2.1f\n",fMinDist);
4395  printf("Min Distance to Bad Channel 2 = %2.1f\n",fMinDist2);
4396  printf("Min Distance to Bad Channel 3 = %2.1f\n",fMinDist3);
4397  printf("Reject clusters with a track matched = %d\n",fRejectTrackMatch);
4398  printf("Time Cut: %3.1f < TOF < %3.1f\n", fTimeCutMin, fTimeCutMax);
4399  printf("Number of cells in cluster is > %d \n", fNCellsCut);
4400  printf(" \n") ;
4401 }
Bool_t IsPileUpFromSPD() const
Bool_t fFillOnlySimpleSSHisto
Fill selected cluster histograms, selected SS histograms.
Definition: AliAnaPhoton.h:192
Float_t GetHistoPtMax() const
TH2F * fhEBinClusterColRowPID[14]
! Column and row location of cluster max E cell in different energy bins, after PID cut ...
Definition: AliAnaPhoton.h:476
TH2F * fhLam1ETMTRD
! Cluster lambda1 vs E, SM covered by TRD, cut on Track Matching residual
Definition: AliAnaPhoton.h:262
TH1F * fhClusterCutsE[10]
! control histogram on the different photon selection cuts, E
Definition: AliAnaPhoton.h:216
TH2F * fhPhiLam0LowE
! Cluster phi vs lambda0, E<2
Definition: AliAnaPhoton.h:272
TH2F * fhNCellsLam1HighE
! number of cells in cluster vs lambda1, E>2
Definition: AliAnaPhoton.h:268
TH2F * fhMergeGeneratorClusterHijingBkgEPrimRecoDiff[10][fgkNGenTypes]
! Cluster energy, at least 2 generators contributions, one is HIJING, for different generator origins...
Definition: AliAnaPhoton.h:510
TH2F * fhLam1Lam0HighE
! Cluster lambda1 vs lambda0, E>2
Definition: AliAnaPhoton.h:278
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:375
TH1F * fhMCPt[fgkNmcTypes]
! Number of identified photon vs cluster pT coming from MC particle
Definition: AliAnaPhoton.h:303
TH2F * fhTrackMatchedDEtaDPhiNeg[2]
! Eta vs Phi distance between track and cluster, E cluster > 0.5 GeV, after and before photon cuts ...
Definition: AliAnaPhoton.h:373
Int_t fNCellsCut
Accept for the analysis clusters with more than fNCellsCut cells.
Definition: AliAnaPhoton.h:181
TH2F * fhEmbedPi0ELambda0FullBkg
! Lambda0 vs E for embedded photons with less than 10% of the cluster energy
Definition: AliAnaPhoton.h:359
TH2F * fhDispE
! Cluster dispersion vs E
Definition: AliAnaPhoton.h:237
TH2F * fhMCEta[fgkNmcTypes]
! eta of identified photon coming from MC particle
Definition: AliAnaPhoton.h:305
virtual AliMCEvent * GetMC() const
Int_t fNLMCutMax
Remove clusters/cells with number of local maxima larger than this value.
Definition: AliAnaPhoton.h:184
TH2F * fhMCNCellsE[fgkNssTypes]
! NCells per cluster vs energy
Definition: AliAnaPhoton.h:336
TH1F * fhMergeGeneratorClusterHijingBkg[10][fgkNGenTypes]
! Cluster energy, at least 2 generators contributions, one is HIJING, for different generator origins...
Definition: AliAnaPhoton.h:497
TH2F * fhDTimeLam0BinPerSMWeighted[2][20]
! t_max-t_cell, not maximum cluster cell, in a l0 bin per SM, log weight Cell E / Cluster E ...
Definition: AliAnaPhoton.h:433
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:291
TH2F * fhCellClusterEFracAndTime[2][7]
! Cell in Cluster E cell/ E cluster vs cell Time, in a l0 bin, different Pt bins
Definition: AliAnaPhoton.h:445
TH2F * fhLam1ETRD
! Cluster lambda1 vs E, SM covered by TRD
Definition: AliAnaPhoton.h:252
TH2F * fhEClusterSM
! Cluster E distribution per SM, before any selection, after reader
Definition: AliAnaPhoton.h:406
TH2F * fhMaxCellDiffClusterE
! Fraction of energy carried by cell with maximum energy
Definition: AliAnaPhoton.h:220
TH2F * fhPtClusterSM
! Cluster E distribution per SM, before any selection, after reader
Definition: AliAnaPhoton.h:408
double Double_t
Definition: External.C:58
TH2F * fhDispETRD
! Cluster dispersion vs E, SM covered by TRD
Definition: AliAnaPhoton.h:249
Int_t GetHistoShowerShapeBins() const
TH2F * fhLocalRegionClusterMultiplicityPerCentralityAdded[6]
! Cluster multiplicity near cluster, R<0.2, vs centrality percentile, not hijing (added signal) tagge...
Definition: AliAnaPhoton.h:492
Float_t GetHistodEdxMax() const
virtual void AddToHistogramsName(TString add)
virtual AliVCaloCells * GetEMCALCells() const
TH2F * fhLam1PerNLargeTimeInClusterCell[5]
! Cluster lambda1 vs Pt, when any secondary cell has t > 50 ns, per number of large time secondary ce...
Definition: AliAnaPhoton.h:469
TH2F * fhTrackMatchedDEtaDPhiPos[2]
! Eta vs Phi distance between track and cluster, E cluster > 0.5 GeV, after and before ...
Definition: AliAnaPhoton.h:369
TH2F * fhLam0PerNLargeTimeInClusterCell[5]
! Cluster lambda0 vs Pt, when any secondary cell has t > 50 ns, per number of large time secondary ce...
Definition: AliAnaPhoton.h:468
Definition: External.C:236
TH2F * fhMCPtLambda0[fgkNssTypes]
! pT vs Lambda0 from MC particle
Definition: AliAnaPhoton.h:322
TH2F * fhNCellsLam1LowE
! number of cells in cluster vs lambda1
Definition: AliAnaPhoton.h:265
Bool_t fFillTMHisto
Fill track matching plots.
Definition: AliAnaPhoton.h:176
TH2F * fhPtPhotonNPileUpSPDVtx
! photon pt vs number of spd pile-up vertices
Definition: AliAnaPhoton.h:399
Bool_t fFillEBinAcceptanceHisto
Fill histograms with cluster eta-phi distribution and column-row cell, for different energy bins...
Definition: AliAnaPhoton.h:206
TH2F * fhPtPhotonNPileUpTrkVtx
! photon pt vs number of track pile-up vertices
Definition: AliAnaPhoton.h:400
TH2F * fhNCellsDispLowE
! number of cells in cluster vs dispersion
Definition: AliAnaPhoton.h:266
TH2F * fhTrackMatchedDEtaPos[2]
! Eta distance between track and cluster vs cluster E, after and before photon cuts ...
Definition: AliAnaPhoton.h:367
TH2F * fhSphericityE
! shower sphericity in eta vs phi
Definition: AliAnaPhoton.h:288
Int_t GetHistoEDiffBins() const
TH2F * fhDispEtaE
! shower dispersion in eta direction
Definition: AliAnaPhoton.h:282
TH2F * fhPtPhotonNPileUpTrkVtxTimeCut
! photon pt vs number of track pile-up vertices, time cut +- 25 ns
Definition: AliAnaPhoton.h:402
virtual Int_t GetModuleNumberCellIndexesAbsCaloMap(Int_t absId, Int_t calo, Int_t &icol, Int_t &irow, Int_t &iRCU, Int_t &icolAbs, Int_t &irowAbs) const
TH2F * fhEtaPrimMC[fgkNmcPrimTypes]
! Eta of generated photon
Definition: AliAnaPhoton.h:311
AliEMCALRecoUtils * GetEMCALRecoUtils() const
TH2F * fhPtPhotonSM
! photon-like cluster E distribution per SM
Definition: AliAnaPhoton.h:409
TH2F * fhEBinClusterEtaPhiPID[14]
! Eta-Phi location of cluster in different energy bins, after PID cut
Definition: AliAnaPhoton.h:475
Bool_t ReadAODMCParticles() const
Bool_t fFillEMCALRegionSSHistograms
Fill shower shape histograms in EMCal slices.
Definition: AliAnaPhoton.h:188
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:376
Float_t fMinDist2
Cuts on Minimal distance to study acceptance evaluation.
Definition: AliAnaPhoton.h:171
TH1F * fhPtPhotonPileUp[7]
! pT distribution of selected photons
Definition: AliAnaPhoton.h:392
TH1F * fhPtPhoton
! Number of identified photon vs transerse momentum
Definition: AliAnaPhoton.h:225
Bool_t fRejectTrackMatch
If PID on, reject clusters which have an associated TPC track.
Definition: AliAnaPhoton.h:174
TH2F * fhMergeGeneratorClusterHijingAndOtherBkgEPrimRecoDiff[10][fgkNGenTypes]
! Cluster energy, at least 3 generators contributions, one is HIJING, for different generator origins...
Definition: AliAnaPhoton.h:509
void MakeAnalysisFillAOD()
Float_t GetPhi(Float_t phi) const
Shift phi angle in case of negative value 360 degrees. Example TLorenzVector::Phi defined in -pi to p...
TLorentzVector fMomentum2
! Cluster momentum, temporary container
Definition: AliAnaPhoton.h:200
TH2F * fhTrackMatchedDPhiMCNoOverlap[2]
! Phi distance between track and cluster vs cluster E, not other particle overlap, after and before photon cuts
Definition: AliAnaPhoton.h:381
TH2F * fhPtEventPlanePhoton
! event plane vs photon pT
Definition: AliAnaPhoton.h:232
TH2F * fhLam0PtNLM1
! Cluster lambda0 vs pT, for clusters with NLM=1
Definition: AliAnaPhoton.h:244
TH2F * fhLam1PtNLM1
! Cluster lambda0 vs pT, for clusters with NLM=1
Definition: AliAnaPhoton.h:246
virtual AliVEvent * GetInputEvent() const
TH2F * fhLam0PtTRD
! Cluster lambda0 vs pT, SM covered by TRD
Definition: AliAnaPhoton.h:251
TH2F * fhMCConversionVertex
! Conversion distance for photon clusters that have at least a contributor from the conversion...
Definition: AliAnaPhoton.h:411
TH2F * fhEmbedPi0ELambda0FullSignal
! Lambda0 vs E for embedded photons with more than 90% of the cluster energy
Definition: AliAnaPhoton.h:356
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:358
TH2F * fhMCESphericity[fgkNssTypes]
! shower sphericity, eta vs phi
Definition: AliAnaPhoton.h:343
TH2F * fhEtaPhiLargeTimeInClusterCell[2][7]
! Cluster eta/phi, with at least one significant cell with large time
Definition: AliAnaPhoton.h:456
TH2F * fhLocalRegionClusterMultiplicityPerCentralityHijing[6]
! Cluster multiplicity near cluster, R<0.2, vs centrality percentile, hijing tagged mc clusters...
Definition: AliAnaPhoton.h:487
TH2F * fhMC2Pt[fgkNmcTypes]
! pT distribution, Reco vs MC coming from MC particle
Definition: AliAnaPhoton.h:300
TH2F * fhEtaPhiPhoton
! Pseudorapidity vs Phi of identified photon for E > 0.5
Definition: AliAnaPhoton.h:228
virtual Float_t GetZvertexCut() const
Maximal number of events for mixin.
TH2F * fhPhiPhoton
! Azimuthal angle of identified photon vs transerse momentum
Definition: AliAnaPhoton.h:226
virtual Double_t GetEventPlaneAngle() const
TH2F * fhYPrimMC[fgkNmcPrimTypes]
! Rapidity of generated photon
Definition: AliAnaPhoton.h:310
TH2F * fhTrackMatchedDEtaDPhi[2]
! Eta vs Phi distance between track and cluster, E cluster > 0.5 GeV, after and before ...
Definition: AliAnaPhoton.h:365
void FillAcceptanceHistograms()
TH2F * fhSumEtaE
! shower dispersion in eta direction
Definition: AliAnaPhoton.h:284
Float_t GetHistoRatioMax() const
virtual Bool_t IsTrackMatched(AliVCluster *cluster, AliVEvent *event)
TH2F * fhLocalRegionClusterEnergySum[6]
! Sum of energy near the cluster, R<0.2, vs cluster E, for different cocktail merging cases ...
Definition: AliAnaPhoton.h:479
TH2F * fhMCEDispEtaPhiDiff[fgkNssTypes]
! shower dispersion in eta -phi direction
Definition: AliAnaPhoton.h:342
TH2F * fhColRowLam0BinPtBin[2][7]
! Cell hits, not maximum cluster cell, in a l0 bin, different Pt bins
Definition: AliAnaPhoton.h:441
void Print(const Option_t *opt) const
Print some relevant parameters set for the analysis.
TH2F * fhMergeGeneratorClusterEPrimRecoRatio[10][fgkNGenTypes]
! Cluster energy, at least 2 generators contributions, for different generator origins and different ...
Definition: AliAnaPhoton.h:501
TH2F * fhCellClusterELam0BinPerSM[2][20]
! Cell E vs cluster pT, not maximum cluster cell, in a l0 bin per SM
Definition: AliAnaPhoton.h:436
Int_t GetHistoPhiBins() const
TH2F * fhEtaLam0LowE
! Cluster eta vs lambda0, E<2
Definition: AliAnaPhoton.h:271
Float_t fEBinCuts[15]
Energy bins cut.
Definition: AliAnaPhoton.h:207
TH2F * fhPtPhotonNPileUpSPDVtxTimeCut
! photon pt vs number of spd pile-up vertices, time cut +-25 ns
Definition: AliAnaPhoton.h:401
TH2F * fhMCConversionVertexTRD
! Conversion distance for photon clusters that have at least a contributor from the conversion...
Definition: AliAnaPhoton.h:412
Bool_t fFillSSHistograms
Fill shower shape histograms.
Definition: AliAnaPhoton.h:186
TH2F * fhdEdx[2]
! Matched track dEdx vs cluster E, after and before photon cuts
Definition: AliAnaPhoton.h:386
TH2F * fhLam0ETRD
! Cluster lambda0 vs E, SM covered by TRD
Definition: AliAnaPhoton.h:250
Int_t GetCocktailGeneratorBackgroundTag(AliVCluster *clus, TString &genName, TString &genNameBkg)
TH2F * fhLocalRegionClusterEnergySumHijing[6]
! Sum of energy near the cluster, R<0.2, vs cluster E, hijing tagged mc clusters, for different cockt...
Definition: AliAnaPhoton.h:484
TH2F * fhDispETM
! Cluster dispersion vs E, cut on Track Matching residual
Definition: AliAnaPhoton.h:254
TString GetCocktailGenNameToCheck(Int_t i) const
TH2F * fhEmbedPhotonELambda0FullSignal
! Lambda0 vs E for embedded photons with more than 90% of the cluster energy
Definition: AliAnaPhoton.h:351
Float_t GetHistoTrackResidualPhiMin() const
TH2F * fhEtaPhoton
! Pseudorapidity of identified photon vs transerse momentum
Definition: AliAnaPhoton.h:227
TH2F * fhDispSumEtaDiffE
! difference of 2 eta dispersions
Definition: AliAnaPhoton.h:289
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:395