comment out unused variable
[u/mrichter/AliRoot.git] / PWGGA / CaloTrackCorrelations / AliAnaPhoton.cxx
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 //_________________________________________________________________________
17 //
18 // Class for the photon identification.
19 // Clusters from calorimeters are identified as photons
20 // and kept in the AOD. Few histograms produced.
21 // Produces input for other analysis classes like AliAnaPi0,
22 // AliAnaParticleHadronCorrelation ...
23 //
24 // -- Author: Gustavo Conesa (LNF-INFN)
25 //////////////////////////////////////////////////////////////////////////////
26
27
28 // --- ROOT system ---
29 #include <TH2F.h>
30 #include <TH3D.h>
31 #include <TClonesArray.h>
32 #include <TObjString.h>
33 #include "TParticle.h"
34 #include "TDatabasePDG.h"
35
36 // --- Analysis system ---
37 #include "AliAnaPhoton.h"
38 #include "AliCaloTrackReader.h"
39 #include "AliStack.h"
40 #include "AliCaloPID.h"
41 #include "AliMCAnalysisUtils.h"
42 #include "AliFiducialCut.h"
43 #include "AliVCluster.h"
44 #include "AliAODMCParticle.h"
45 #include "AliMixedEvent.h"
46 #include "AliAODEvent.h"
47 #include "AliESDEvent.h"
48
49 // --- Detectors ---
50 #include "AliPHOSGeoUtils.h"
51 #include "AliEMCALGeometry.h"
52
53 ClassImp(AliAnaPhoton)
54
55 //____________________________
56 AliAnaPhoton::AliAnaPhoton() :
57 AliAnaCaloTrackCorrBaseClass(), fCalorimeter(""),
58 fMinDist(0.),                 fMinDist2(0.),                fMinDist3(0.),
59 fRejectTrackMatch(0),         fFillTMHisto(kFALSE),
60 fTimeCutMin(-10000),          fTimeCutMax(10000),
61 fNCellsCut(0),
62 fNLMCutMin(-1),               fNLMCutMax(10),
63 fFillSSHistograms(kFALSE),    fFillOnlySimpleSSHisto(1),
64 fNOriginHistograms(8),        fNPrimaryHistograms(4),
65 fFillPileUpHistograms(0),     fFillEMCALBCHistograms(0),
66 // Histograms
67 fhNCellsE(0),                 fhCellsE(0),   // Control histograms
68 fhMaxCellDiffClusterE(0),     fhTimePt(0),   // Control histograms
69 fhEtaPhi(0),                  fhEtaPhiEMCALBC0(0),
70 fhEtaPhiEMCALBC1(0),          fhEtaPhiEMCALBCN(0),
71 fhTimeTriggerEMCALBCCluster(0),
72 fhTimeTriggerEMCALBCUMCluster(0),
73 fhEtaPhiTriggerEMCALBCClusterOverTh(0),
74 fhEtaPhiTriggerEMCALBCUMClusterOverTh(0),
75 fhEtaPhiTriggerEMCALBCClusterBelowTh1(0),
76 fhEtaPhiTriggerEMCALBCUMClusterBelowTh1(0),
77 fhEtaPhiTriggerEMCALBCClusterBelowTh2(0),
78 fhEtaPhiTriggerEMCALBCUMClusterBelowTh2(0),
79 fhEtaPhiTriggerEMCALBCExotic(0),             fhTimeTriggerEMCALBCExotic(0),
80 fhEtaPhiTriggerEMCALBCUMExotic(0),           fhTimeTriggerEMCALBCUMExotic(0),
81 fhEtaPhiTriggerEMCALBCBad(0),                fhTimeTriggerEMCALBCBad(0),
82 fhEtaPhiTriggerEMCALBCUMBad(0),              fhTimeTriggerEMCALBCUMBad(0),
83 fhEtaPhiTriggerEMCALBCBadExotic(0),          fhTimeTriggerEMCALBCBadExotic(0),
84 fhEtaPhiTriggerEMCALBCUMBadExotic(0),        fhTimeTriggerEMCALBCUMBadExotic(0),
85 fhEtaPhiTriggerEMCALBCExoticCluster(0),      fhTimeTriggerEMCALBCExoticCluster(0),
86 fhEtaPhiTriggerEMCALBCUMExoticCluster(0),    fhTimeTriggerEMCALBCUMExoticCluster(0),
87 fhEtaPhiTriggerEMCALBCBadCluster(0),         fhTimeTriggerEMCALBCBadCluster(0),
88 fhEtaPhiTriggerEMCALBCUMBadCluster(0),       fhTimeTriggerEMCALBCUMBadCluster(0),
89 fhEtaPhiTriggerEMCALBCBadExoticCluster(0),   fhTimeTriggerEMCALBCBadExoticCluster(0),
90 fhEtaPhiTriggerEMCALBCUMBadExoticCluster(0), fhTimeTriggerEMCALBCUMBadExoticCluster(0),
91 fhTimeTriggerEMCALBCBadMaxCell(0),           fhTimeTriggerEMCALBCUMBadMaxCell(0),
92 fhTimeTriggerEMCALBCBadMaxCellExotic(0),     fhTimeTriggerEMCALBCUMBadMaxCellExotic(0),
93 fhEtaPhiTriggerEMCALBCUMReMatchOpenTimeCluster (0), fhTimeTriggerEMCALBCUMReMatchOpenTimeCluster(0),
94 fhEtaPhiTriggerEMCALBCUMReMatchCheckNeighCluster(0),fhTimeTriggerEMCALBCUMReMatchCheckNeighCluster(0),
95 fhEtaPhiTriggerEMCALBCUMReMatchBothCluster(0),      fhTimeTriggerEMCALBCUMReMatchBothCluster(0),
96 fhTimeTriggerEMCALBC0UMReMatchOpenTime(0),
97 fhTimeTriggerEMCALBC0UMReMatchCheckNeigh(0),
98 fhTimeTriggerEMCALBC0UMReMatchBoth(0),
99
100 fhEtaPhiNoTrigger(0),                        fhTimeNoTrigger(0),
101
102 fhEPhoton(0),                 fhPtPhoton(0),
103 fhPhiPhoton(0),               fhEtaPhoton(0),
104 fhEtaPhiPhoton(0),            fhEtaPhi05Photon(0),
105 fhEtaPhiPhotonEMCALBC0(0),    fhEtaPhiPhotonEMCALBC1(0),   fhEtaPhiPhotonEMCALBCN(0),
106 fhTimePhotonTriggerEMCALBC0UMReMatchOpenTime(0),
107 fhTimePhotonTriggerEMCALBC0UMReMatchCheckNeigh(0),
108 fhTimePhotonTriggerEMCALBC0UMReMatchBoth(0),
109
110 fhPtCentralityPhoton(0),      fhPtEventPlanePhoton(0),
111
112 // Shower shape histograms
113 fhNLocMax(0),
114 fhDispE(0),                   fhLam0E(0),                   fhLam1E(0),
115 fhDispETRD(0),                fhLam0ETRD(0),                fhLam1ETRD(0),
116 fhDispETM(0),                 fhLam0ETM(0),                 fhLam1ETM(0),
117 fhDispETMTRD(0),              fhLam0ETMTRD(0),              fhLam1ETMTRD(0),
118
119 fhNCellsLam0LowE(0),          fhNCellsLam1LowE(0),          fhNCellsDispLowE(0),
120 fhNCellsLam0HighE(0),         fhNCellsLam1HighE(0),         fhNCellsDispHighE(0),
121
122 fhEtaLam0LowE(0),             fhPhiLam0LowE(0),
123 fhEtaLam0HighE(0),            fhPhiLam0HighE(0),
124 fhLam0DispLowE(0),            fhLam0DispHighE(0),
125 fhLam1Lam0LowE(0),            fhLam1Lam0HighE(0),
126 fhDispLam1LowE(0),            fhDispLam1HighE(0),
127 fhDispEtaE(0),                fhDispPhiE(0),
128 fhSumEtaE(0),                 fhSumPhiE(0),                 fhSumEtaPhiE(0),
129 fhDispEtaPhiDiffE(0),         fhSphericityE(0),
130 fhDispSumEtaDiffE(0),         fhDispSumPhiDiffE(0),
131
132 // MC histograms
133 fhMCPhotonELambda0NoOverlap(0),       fhMCPhotonELambda0TwoOverlap(0),      fhMCPhotonELambda0NOverlap(0),
134 // Embedding
135 fhEmbeddedSignalFractionEnergy(0),
136 fhEmbedPhotonELambda0FullSignal(0),   fhEmbedPhotonELambda0MostlySignal(0),
137 fhEmbedPhotonELambda0MostlyBkg(0),    fhEmbedPhotonELambda0FullBkg(0),
138 fhEmbedPi0ELambda0FullSignal(0),      fhEmbedPi0ELambda0MostlySignal(0),
139 fhEmbedPi0ELambda0MostlyBkg(0),       fhEmbedPi0ELambda0FullBkg(0),
140 // PileUp
141 fhTimePtNoCut(0),                     fhTimePtSPD(0),
142 fhTimePtPhotonNoCut(0),               fhTimePtPhotonSPD(0),
143 fhTimeNPileUpVertSPD(0),              fhTimeNPileUpVertTrack(0),
144 fhTimeNPileUpVertContributors(0),
145 fhTimePileUpMainVertexZDistance(0),   fhTimePileUpMainVertexZDiamond(0),
146 fhClusterMultSPDPileUp(),             fhClusterMultNoPileUp(),
147 fhEtaPhiBC0(0),  fhEtaPhiBCPlus(0),   fhEtaPhiBCMinus(0),
148 fhEtaPhiBC0PileUpSPD(0),
149 fhEtaPhiBCPlusPileUpSPD(0),           fhEtaPhiBCMinusPileUpSPD(0),
150 fhPtNPileUpSPDVtx(0),                 fhPtNPileUpTrkVtx(0),
151 fhPtNPileUpSPDVtxTimeCut(0),          fhPtNPileUpTrkVtxTimeCut(0),
152 fhPtNPileUpSPDVtxTimeCut2(0),         fhPtNPileUpTrkVtxTimeCut2(0),
153 fhPtPhotonNPileUpSPDVtx(0),           fhPtPhotonNPileUpTrkVtx(0),
154 fhPtPhotonNPileUpSPDVtxTimeCut(0),    fhPtPhotonNPileUpTrkVtxTimeCut(0),
155 fhPtPhotonNPileUpSPDVtxTimeCut2(0),   fhPtPhotonNPileUpTrkVtxTimeCut2(0),
156 fhEClusterSM(0),                      fhEPhotonSM(0),
157 fhPtClusterSM(0),                     fhPtPhotonSM(0)
158 {
159   //default ctor
160   
161   for(Int_t i = 0; i < 14; i++)
162   {
163     fhMCPt     [i] = 0;
164     fhMCE      [i] = 0;
165     fhMCPhi    [i] = 0;
166     fhMCEta    [i] = 0;
167     fhMCDeltaE [i] = 0;
168     fhMCDeltaPt[i] = 0;
169     fhMC2E     [i] = 0;
170     fhMC2Pt    [i] = 0;
171   }
172   
173   for(Int_t i = 0; i < 7; i++)
174   {
175     fhPtPrimMC [i] = 0;
176     fhEPrimMC  [i] = 0;
177     fhPhiPrimMC[i] = 0;
178     fhEtaPrimMC[i] = 0;
179     fhYPrimMC  [i] = 0;
180     
181     fhPtPrimMCAcc [i] = 0;
182     fhEPrimMCAcc  [i] = 0;
183     fhPhiPrimMCAcc[i] = 0;
184     fhEtaPrimMCAcc[i] = 0;
185     fhYPrimMCAcc  [i] = 0;
186     
187     fhDispEtaDispPhi[i] = 0;
188     fhLambda0DispPhi[i] = 0;
189     fhLambda0DispEta[i] = 0;
190     
191     fhPtPileUp       [i] = 0;
192     fhPtChargedPileUp[i] = 0;
193     fhPtPhotonPileUp [i] = 0;
194     
195     fhLambda0PileUp       [i] = 0;
196     fhLambda0ChargedPileUp[i] = 0;
197     
198     fhClusterEFracLongTimePileUp  [i] = 0;
199     
200     fhClusterCellTimePileUp       [i] = 0;
201     fhClusterTimeDiffPileUp       [i] = 0;
202     fhClusterTimeDiffChargedPileUp[i] = 0;
203     fhClusterTimeDiffPhotonPileUp [i] = 0;
204     
205     for(Int_t j = 0; j < 6; j++)
206     {
207       fhMCDispEtaDispPhi[i][j] = 0;
208       fhMCLambda0DispEta[i][j] = 0;
209       fhMCLambda0DispPhi[i][j] = 0;
210     }
211   }
212   
213   for(Int_t i = 0; i < 6; i++)
214   {
215     fhMCELambda0    [i]                  = 0;
216     fhMCELambda1    [i]                  = 0;
217     fhMCEDispersion [i]                  = 0;
218     fhMCNCellsE     [i]                  = 0;
219     fhMCMaxCellDiffClusterE[i]           = 0;
220     fhLambda0DispEta[i]                  = 0;
221     fhLambda0DispPhi[i]                  = 0;
222     
223     fhMCLambda0vsClusterMaxCellDiffE0[i] = 0;
224     fhMCLambda0vsClusterMaxCellDiffE2[i] = 0;
225     fhMCLambda0vsClusterMaxCellDiffE6[i] = 0;
226     fhMCNCellsvsClusterMaxCellDiffE0 [i] = 0;
227     fhMCNCellsvsClusterMaxCellDiffE2 [i] = 0;
228     fhMCNCellsvsClusterMaxCellDiffE6 [i] = 0;
229     
230     fhMCEDispEta       [i]               = 0;
231     fhMCEDispPhi       [i]               = 0;
232     fhMCESumEtaPhi     [i]               = 0;
233     fhMCEDispEtaPhiDiff[i]               = 0;
234     fhMCESphericity    [i]               = 0;
235   }
236   
237   for(Int_t i = 0; i < 5; i++)
238   {
239     fhClusterCuts[i] = 0;
240   }
241   
242   // Track matching residuals
243   for(Int_t i = 0; i < 2; i++)
244   {
245     fhTrackMatchedDEta   [i] = 0;             fhTrackMatchedDPhi   [i] = 0;         fhTrackMatchedDEtaDPhi   [i] = 0;
246     fhTrackMatchedDEtaNeg[i] = 0;             fhTrackMatchedDPhiNeg[i] = 0;         fhTrackMatchedDEtaDPhiNeg[i] = 0;
247     fhTrackMatchedDEtaPos[i] = 0;             fhTrackMatchedDPhiPos[i] = 0;         fhTrackMatchedDEtaDPhiPos[i] = 0;
248     fhTrackMatchedDEtaTRD[i] = 0;             fhTrackMatchedDPhiTRD[i] = 0;
249     fhTrackMatchedDEtaMCOverlap[i] = 0;       fhTrackMatchedDPhiMCOverlap[i] = 0;
250     fhTrackMatchedDEtaMCNoOverlap[i] = 0;     fhTrackMatchedDPhiMCNoOverlap[i] = 0;
251     fhTrackMatchedDEtaMCConversion[i] = 0;    fhTrackMatchedDPhiMCConversion[i] = 0;
252     fhTrackMatchedMCParticle[i] = 0;          fhTrackMatchedMCParticle[i] = 0;
253     fhdEdx[i] = 0;                            fhEOverP[i] = 0;
254     fhEOverPTRD[i] = 0;
255   }
256   
257   for(Int_t i = 0; i < 4; i++)
258   {
259     fhClusterMultSPDPileUp[i] = 0;
260     fhClusterMultNoPileUp [i] = 0;
261   }
262   
263   for(Int_t i = 0; i < 11; i++)
264   {
265     fhEtaPhiTriggerEMCALBC             [i] = 0 ;
266     fhTimeTriggerEMCALBC               [i] = 0 ;
267     fhEtaPhiTriggerEMCALBCUM           [i] = 0 ;
268     fhTimeTriggerEMCALBCUM             [i] = 0 ;
269     
270     fhEtaPhiPhotonTriggerEMCALBC       [i] = 0 ;
271     fhTimePhotonTriggerEMCALBC         [i] = 0 ;
272     fhEtaPhiPhotonTriggerEMCALBCUM     [i] = 0 ;
273     fhTimePhotonTriggerEMCALBCUM       [i] = 0 ;
274     
275     fhTimePhotonTriggerEMCALBCPileUpSPD[i] = 0 ;
276     fhTimeTriggerEMCALBCPileUpSPD      [i] = 0 ;
277     
278     fhEtaPhiTriggerEMCALBCCluster      [i] = 0 ;
279     fhEtaPhiTriggerEMCALBCUMCluster    [i] = 0 ;    
280   }
281   
282   //Initialize parameters
283   InitParameters();
284   
285 }
286
287 //_________________________________________________________________________________________
288 Bool_t  AliAnaPhoton::ClusterSelected(AliVCluster* calo, TLorentzVector mom, Int_t nMaxima)
289 {
290   //Select clusters if they pass different cuts
291   
292   Float_t ptcluster  = mom.Pt();
293   Float_t ecluster   = mom.E();
294   Float_t l0cluster  = calo->GetM02();
295   Float_t etacluster = mom.Eta();
296   Float_t phicluster = mom.Phi();
297   if(phicluster < 0) phicluster+=TMath::TwoPi();
298   Float_t tofcluster   = calo->GetTOF()*1.e9;
299   
300   Bool_t matched = IsTrackMatched(calo,GetReader()->GetInputEvent());
301   
302   if(GetDebug() > 2)
303     printf("AliAnaPhoton::ClusterSelected() - Current Event %d; Before selection : E %2.2f, pT %2.2f, phi %2.2f, eta %2.2f\n",
304            GetReader()->GetEventNumber(),
305            ecluster,ptcluster, phicluster*TMath::RadToDeg(),etacluster);
306   
307   fhClusterCuts[1]->Fill(ecluster);
308   
309   if(ecluster > 0.5) fhEtaPhi->Fill(etacluster, phicluster);
310   
311   Int_t   nSM  = GetModuleNumber(calo);
312   if(nSM < GetCaloUtils()->GetNumberOfSuperModulesUsed() && nSM >=0)
313   {
314     fhEClusterSM ->Fill(ecluster ,nSM);
315     fhPtClusterSM->Fill(ptcluster,nSM);
316   }
317   
318   FillEMCALTriggerClusterBCHistograms(calo->GetID(),ecluster,tofcluster,etacluster,phicluster);
319   
320   //.......................................
321   //If too small or big energy, skip it
322   if(ecluster < GetMinEnergy() || ecluster > GetMaxEnergy() ) return kFALSE ;
323   
324   if(GetDebug() > 2) printf("\t Cluster %d Pass E Cut \n",calo->GetID());
325   
326   fhClusterCuts[2]->Fill(ecluster);
327   
328   FillClusterPileUpHistograms(calo,matched,ptcluster,etacluster,phicluster,l0cluster);
329   
330   //.......................................
331   // TOF cut, BE CAREFUL WITH THIS CUT
332   Double_t tof = calo->GetTOF()*1e9;
333   if(tof < fTimeCutMin || tof > fTimeCutMax) return kFALSE;
334   
335   if(GetDebug() > 2)  printf("\t Cluster %d Pass Time Cut \n",calo->GetID());
336   
337   fhClusterCuts[3]->Fill(ecluster);
338   
339   //.......................................
340   if(calo->GetNCells() <= fNCellsCut && GetReader()->GetDataType() != AliCaloTrackReader::kMC) return kFALSE;
341   
342   if(GetDebug() > 2) printf("\t Cluster %d Pass NCell Cut \n",calo->GetID());
343   
344   fhClusterCuts[4]->Fill(ecluster);
345   
346   if(nMaxima < fNLMCutMin || nMaxima > fNLMCutMax) return kFALSE ;
347   if(GetDebug() > 2) printf(" \t Cluster %d pass NLM %d of out of range \n",calo->GetID(), nMaxima);
348   
349   fhClusterCuts[5]->Fill(ecluster);
350   
351   //.......................................
352   //Check acceptance selection
353   if(IsFiducialCutOn())
354   {
355     Bool_t in = GetFiducialCut()->IsInFiducialCut(mom,fCalorimeter) ;
356     if(! in ) return kFALSE ;
357   }
358   
359   if(GetDebug() > 2) printf("\t Fiducial cut passed \n");
360   
361   fhClusterCuts[6]->Fill(ecluster);
362   
363   //.......................................
364   //Skip matched clusters with tracks
365   
366   // Fill matching residual histograms before PID cuts
367   if(fFillTMHisto) FillTrackMatchingResidualHistograms(calo,0);
368   
369   if(fRejectTrackMatch)
370   {
371     if(matched)
372     {
373       if(GetDebug() > 2) printf("\t Reject track-matched clusters\n");
374       return kFALSE ;
375     }
376     else
377       if(GetDebug() > 2)  printf(" Track-matching cut passed \n");
378   }// reject matched clusters
379   
380   fhClusterCuts[7]->Fill(ecluster);
381   
382   if(fFillPileUpHistograms)
383   {
384     if(GetReader()->IsPileUpFromSPD())               {fhPtChargedPileUp[0]->Fill(ptcluster); fhLambda0ChargedPileUp[0]->Fill(ecluster,l0cluster); }
385     if(GetReader()->IsPileUpFromEMCal())             {fhPtChargedPileUp[1]->Fill(ptcluster); fhLambda0ChargedPileUp[1]->Fill(ecluster,l0cluster); }
386     if(GetReader()->IsPileUpFromSPDOrEMCal())        {fhPtChargedPileUp[2]->Fill(ptcluster); fhLambda0ChargedPileUp[2]->Fill(ecluster,l0cluster); }
387     if(GetReader()->IsPileUpFromSPDAndEMCal())       {fhPtChargedPileUp[3]->Fill(ptcluster); fhLambda0ChargedPileUp[3]->Fill(ecluster,l0cluster); }
388     if(GetReader()->IsPileUpFromSPDAndNotEMCal())    {fhPtChargedPileUp[4]->Fill(ptcluster); fhLambda0ChargedPileUp[4]->Fill(ecluster,l0cluster); }
389     if(GetReader()->IsPileUpFromEMCalAndNotSPD())    {fhPtChargedPileUp[5]->Fill(ptcluster); fhLambda0ChargedPileUp[5]->Fill(ecluster,l0cluster); }
390     if(GetReader()->IsPileUpFromNotSPDAndNotEMCal()) {fhPtChargedPileUp[6]->Fill(ptcluster); fhLambda0ChargedPileUp[6]->Fill(ecluster,l0cluster); }
391   }
392   
393   //.......................................
394   //Check Distance to Bad channel, set bit.
395   Double_t distBad=calo->GetDistanceToBadChannel() ; //Distance to bad channel
396   if(distBad < 0.) distBad=9999. ; //workout strange convension dist = -1. ;
397   if(distBad < fMinDist)
398   {//In bad channel (PHOS cristal size 2.2x2.2 cm), EMCAL ( cell units )
399     return kFALSE ;
400   }
401   else if(GetDebug() > 2) printf("\t Bad channel cut passed %4.2f > %2.2f \n",distBad, fMinDist);
402   
403   fhClusterCuts[8]->Fill(ecluster);
404   
405   if(GetDebug() > 0)
406     printf("AliAnaPhoton::ClusterSelected() Current Event %d; After  selection : E %2.2f, pT %2.2f, phi %2.2f, eta %2.2f\n",
407            GetReader()->GetEventNumber(),
408            ecluster, ptcluster,mom.Phi()*TMath::RadToDeg(),mom.Eta());
409   
410   //All checks passed, cluster selected
411   return kTRUE;
412   
413 }
414
415 //___________________________________________
416 void AliAnaPhoton::FillAcceptanceHistograms()
417 {
418   //Fill acceptance histograms if MC data is available
419   
420   Double_t photonY   = -100 ;
421   Double_t photonE   = -1 ;
422   Double_t photonPt  = -1 ;
423   Double_t photonPhi =  100 ;
424   Double_t photonEta = -1 ;
425   
426   Int_t    pdg       =  0 ;
427   Int_t    tag       =  0 ;
428   Int_t    mcIndex   =  0 ;
429   Bool_t   inacceptance = kFALSE;
430   
431   if(GetReader()->ReadStack())
432   {
433     AliStack * stack = GetMCStack();
434     if(stack)
435     {
436       for(Int_t i=0 ; i<stack->GetNtrack(); i++)
437       {
438         if(GetReader()->AcceptOnlyHIJINGLabels() && !GetReader()->IsHIJINGLabel(i)) continue ;
439         
440         TParticle * prim = stack->Particle(i) ;
441         pdg = prim->GetPdgCode();
442         //printf("i %d, %s %d  %s %d \n",i, stack->Particle(i)->GetName(), stack->Particle(i)->GetPdgCode(),
443         //                             prim->GetName(), prim->GetPdgCode());
444         
445         if(pdg == 22)
446         {
447           // Get tag of this particle photon from fragmentation, decay, prompt ...
448           tag = GetMCAnalysisUtils()->CheckOrigin(i,GetReader());
449           if(!GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPhoton))
450           {
451             //A conversion photon from a hadron, skip this kind of photon
452             // printf("AliAnaPhoton::FillAcceptanceHistograms() - not a photon, weird!\n ");
453             // GetMCAnalysisUtils()->PrintMCTag(tag);
454             
455             return;
456           }
457           
458           //Get photon kinematics
459           if(prim->Energy() == TMath::Abs(prim->Pz()))  continue ; //Protection against floating point exception
460           
461           photonY   = 0.5*TMath::Log((prim->Energy()-prim->Pz())/(prim->Energy()+prim->Pz())) ;
462           photonE   = prim->Energy() ;
463           photonPt  = prim->Pt() ;
464           photonPhi = TMath::RadToDeg()*prim->Phi() ;
465           if(photonPhi < 0) photonPhi+=TMath::TwoPi();
466           photonEta = prim->Eta() ;
467           
468           //Check if photons hit the Calorimeter
469           TLorentzVector lv;
470           prim->Momentum(lv);
471           inacceptance = kFALSE;
472           if     (fCalorimeter == "PHOS")
473           {
474             if(GetPHOSGeometry() && GetCaloUtils()->IsPHOSGeoMatrixSet())
475             {
476               Int_t mod ;
477               Double_t x,z ;
478               if(GetPHOSGeometry()->ImpactOnEmc(prim,mod,z,x))
479                 inacceptance = kTRUE;
480               if(GetDebug() > 2) printf("In %s Real acceptance? %d\n",fCalorimeter.Data(),inacceptance);
481             }
482             else
483             {
484               if(GetFiducialCut()->IsInFiducialCut(lv,fCalorimeter))
485                 inacceptance = kTRUE ;
486               if(GetDebug() > 2) printf("In %s fiducial cut acceptance? %d\n",fCalorimeter.Data(),inacceptance);
487             }
488           }
489           else if(fCalorimeter == "EMCAL" && GetCaloUtils()->IsEMCALGeoMatrixSet())
490           {
491             if(GetEMCALGeometry())
492             {
493               Int_t absID=0;
494               
495               GetEMCALGeometry()->GetAbsCellIdFromEtaPhi(prim->Eta(),prim->Phi(),absID);
496               
497               if( absID >= 0)
498                 inacceptance = kTRUE;
499               
500               //                  if(GetEMCALGeometry()->Impact(phot1) && GetEMCALGeometry()->Impact(phot2))
501               //                    inacceptance = kTRUE;
502               if(GetDebug() > 2) printf("In %s Real acceptance? %d\n",fCalorimeter.Data(),inacceptance);
503             }
504             else
505             {
506               if(GetFiducialCut()->IsInFiducialCut(lv,fCalorimeter))
507                 inacceptance = kTRUE ;
508               if(GetDebug() > 2) printf("In %s fiducial cut acceptance? %d\n",fCalorimeter.Data(),inacceptance);
509             }
510           }       //In EMCAL
511           
512           //Fill histograms
513           fhYPrimMC[kmcPPhoton]->Fill(photonPt, photonY) ;
514           if(TMath::Abs(photonY) < 1.0)
515           {
516             fhEPrimMC  [kmcPPhoton]->Fill(photonE ) ;
517             fhPtPrimMC [kmcPPhoton]->Fill(photonPt) ;
518             fhPhiPrimMC[kmcPPhoton]->Fill(photonE , photonPhi) ;
519             fhEtaPrimMC[kmcPPhoton]->Fill(photonE , photonEta) ;
520           }
521           if(inacceptance)
522           {
523             fhEPrimMCAcc  [kmcPPhoton]->Fill(photonE ) ;
524             fhPtPrimMCAcc [kmcPPhoton]->Fill(photonPt) ;
525             fhPhiPrimMCAcc[kmcPPhoton]->Fill(photonE , photonPhi) ;
526             fhEtaPrimMCAcc[kmcPPhoton]->Fill(photonE , photonEta) ;
527             fhYPrimMCAcc  [kmcPPhoton]->Fill(photonE , photonY) ;
528           }//Accepted
529           
530           //Origin of photon
531           if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPrompt) && fhEPrimMC[kmcPPrompt])
532           {
533             mcIndex = kmcPPrompt;
534           }
535           else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCFragmentation) && fhEPrimMC[kmcPFragmentation])
536           {
537             mcIndex = kmcPFragmentation ;
538           }
539           else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCISR) && fhEPrimMC[kmcPISR])
540           {
541             mcIndex = kmcPISR;
542           }
543           else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0Decay)&& fhEPrimMC[kmcPPi0Decay])
544           {
545             mcIndex = kmcPPi0Decay;
546           }
547           else if( (GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCEtaDecay) ||
548                     GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCOtherDecay)) && fhEPrimMC[kmcPOtherDecay])
549           {
550             mcIndex = kmcPOtherDecay;
551           }
552           else if(fhEPrimMC[kmcPOther])
553           {
554             mcIndex = kmcPOther;
555           }//Other origin
556           
557           fhYPrimMC[mcIndex]->Fill(photonPt, photonY) ;
558           if(TMath::Abs(photonY) < 1.0)
559           {
560             fhEPrimMC  [mcIndex]->Fill(photonE ) ;
561             fhPtPrimMC [mcIndex]->Fill(photonPt) ;
562             fhPhiPrimMC[mcIndex]->Fill(photonE , photonPhi) ;
563             fhEtaPrimMC[mcIndex]->Fill(photonE , photonEta) ;
564           }
565           
566           if(inacceptance)
567           {
568             fhEPrimMCAcc  [mcIndex]->Fill(photonE ) ;
569             fhPtPrimMCAcc [mcIndex]->Fill(photonPt) ;
570             fhPhiPrimMCAcc[mcIndex]->Fill(photonE , photonPhi) ;
571             fhEtaPrimMCAcc[mcIndex]->Fill(photonE , photonEta) ;
572             fhYPrimMCAcc  [mcIndex]->Fill(photonE , photonY) ;
573           }//Accepted
574           
575         }// Primary photon
576       }//loop on primaries
577     }//stack exists and data is MC
578   }//read stack
579   else if(GetReader()->ReadAODMCParticles())
580   {
581     TClonesArray * mcparticles = GetReader()->GetAODMCParticles();
582     if(mcparticles)
583     {
584       Int_t nprim = mcparticles->GetEntriesFast();
585       
586       for(Int_t i=0; i < nprim; i++)
587       {
588         if(GetReader()->AcceptOnlyHIJINGLabels() && !GetReader()->IsHIJINGLabel(i)) continue ;
589         
590         AliAODMCParticle * prim = (AliAODMCParticle *) mcparticles->At(i);
591         
592         pdg = prim->GetPdgCode();
593         
594         if(pdg == 22)
595         {
596           // Get tag of this particle photon from fragmentation, decay, prompt ...
597           tag = GetMCAnalysisUtils()->CheckOrigin(i,GetReader());
598           if(!GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPhoton))
599           {
600             //A conversion photon from a hadron, skip this kind of photon
601             //            printf("AliAnaPhoton::FillAcceptanceHistograms() - not a photon, weird!\n ");
602             //            GetMCAnalysisUtils()->PrintMCTag(tag);
603             
604             return;
605           }
606           
607           //Get photon kinematics
608           if(prim->E() == TMath::Abs(prim->Pz()))  continue ; //Protection against floating point exception
609           
610           photonY   = 0.5*TMath::Log((prim->E()-prim->Pz())/(prim->E()+prim->Pz())) ;
611           photonE   = prim->E() ;
612           photonPt  = prim->Pt() ;
613           photonPhi = prim->Phi() ;
614           if(photonPhi < 0) photonPhi+=TMath::TwoPi();
615           photonEta = prim->Eta() ;
616           
617           //Check if photons hit the Calorimeter
618           TLorentzVector lv;
619           lv.SetPxPyPzE(prim->Px(),prim->Py(),prim->Pz(),prim->E());
620           inacceptance = kFALSE;
621           if     (fCalorimeter == "PHOS")
622           {
623             if(GetPHOSGeometry() && GetCaloUtils()->IsPHOSGeoMatrixSet())
624             {
625               Int_t mod ;
626               Double_t x,z ;
627               Double_t vtx[]={prim->Xv(),prim->Yv(),prim->Zv()};
628               if(GetPHOSGeometry()->ImpactOnEmc(vtx, prim->Theta(),prim->Phi(),mod,z,x))
629                 inacceptance = kTRUE;
630               if(GetDebug() > 2) printf("In %s Real acceptance? %d\n",fCalorimeter.Data(),inacceptance);
631             }
632             else
633             {
634               if(GetFiducialCut()->IsInFiducialCut(lv,fCalorimeter))
635                 inacceptance = kTRUE ;
636               if(GetDebug() > 2) printf("In %s fiducial cut acceptance? %d\n",fCalorimeter.Data(),inacceptance);
637             }
638           }
639           else if(fCalorimeter == "EMCAL" && GetCaloUtils()->IsEMCALGeoMatrixSet())
640           {
641             if(GetEMCALGeometry())
642             {
643               Int_t absID=0;
644               
645               GetEMCALGeometry()->GetAbsCellIdFromEtaPhi(prim->Eta(),prim->Phi(),absID);
646               
647               if( absID >= 0)
648                 inacceptance = kTRUE;
649               
650               if(GetDebug() > 2) printf("In %s Real acceptance? %d\n",fCalorimeter.Data(),inacceptance);
651             }
652             else
653             {
654               if(GetFiducialCut()->IsInFiducialCut(lv,fCalorimeter))
655                 inacceptance = kTRUE ;
656               if(GetDebug() > 2) printf("In %s fiducial cut acceptance? %d\n",fCalorimeter.Data(),inacceptance);
657             }
658           }       //In EMCAL
659           
660           //Fill histograms
661           
662           fhYPrimMC[kmcPPhoton]->Fill(photonPt, photonY) ;
663           if(TMath::Abs(photonY) < 1.0)
664           {
665             fhEPrimMC  [kmcPPhoton]->Fill(photonE ) ;
666             fhPtPrimMC [kmcPPhoton]->Fill(photonPt) ;
667             fhPhiPrimMC[kmcPPhoton]->Fill(photonE , photonPhi) ;
668             fhEtaPrimMC[kmcPPhoton]->Fill(photonE , photonEta) ;
669           }
670           
671           if(inacceptance)
672           {
673             fhEPrimMCAcc[kmcPPhoton]  ->Fill(photonE ) ;
674             fhPtPrimMCAcc[kmcPPhoton] ->Fill(photonPt) ;
675             fhPhiPrimMCAcc[kmcPPhoton]->Fill(photonE , photonPhi) ;
676             fhEtaPrimMCAcc[kmcPPhoton]->Fill(photonE , photonEta) ;
677             fhYPrimMCAcc[kmcPPhoton]  ->Fill(photonE , photonY) ;
678           }//Accepted
679           
680           //Origin of photon
681           if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPrompt) && fhEPrimMC[kmcPPrompt])
682           {
683             mcIndex = kmcPPrompt;
684           }
685           else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCFragmentation) && fhEPrimMC[kmcPFragmentation])
686           {
687             mcIndex = kmcPFragmentation ;
688           }
689           else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCISR) && fhEPrimMC[kmcPISR])
690           {
691             mcIndex = kmcPISR;
692           }
693           else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0Decay)&& fhEPrimMC[kmcPPi0Decay])
694           {
695             mcIndex = kmcPPi0Decay;
696           }
697           else if( (GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCEtaDecay) ||
698                     GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCOtherDecay)) && fhEPrimMC[kmcPOtherDecay])
699           {
700             mcIndex = kmcPOtherDecay;
701           }
702           else if(fhEPrimMC[kmcPOther])
703           {
704             mcIndex = kmcPOther;
705           }//Other origin
706           
707           fhYPrimMC[mcIndex]->Fill(photonPt, photonY) ;
708           if(TMath::Abs(photonY) < 1.0)
709           {
710             fhEPrimMC  [mcIndex]->Fill(photonE ) ;
711             fhPtPrimMC [mcIndex]->Fill(photonPt) ;
712             fhPhiPrimMC[mcIndex]->Fill(photonE , photonPhi) ;
713             fhEtaPrimMC[mcIndex]->Fill(photonE , photonEta) ;
714           }
715           if(inacceptance)
716           {
717             fhEPrimMCAcc  [mcIndex]->Fill(photonE ) ;
718             fhPtPrimMCAcc [mcIndex]->Fill(photonPt) ;
719             fhPhiPrimMCAcc[mcIndex]->Fill(photonE , photonPhi) ;
720             fhEtaPrimMCAcc[mcIndex]->Fill(photonE , photonEta) ;
721             fhYPrimMCAcc  [mcIndex]->Fill(photonE , photonY) ;
722           }//Accepted
723           
724         }// Primary photon
725       }//loop on primaries
726       
727     }//kmc array exists and data is MC
728   }     // read AOD MC
729 }
730
731 //________________________________________________________________________________________________________________
732 void  AliAnaPhoton::FillEMCALTriggerClusterBCHistograms(Int_t idcalo,       Float_t ecluster,  Float_t tofcluster,
733                                                         Float_t etacluster, Float_t phicluster)
734
735 {
736   // Fill trigger related histograms
737   
738   if(!fFillEMCALBCHistograms || fCalorimeter!="EMCAL") return ;
739   
740   Float_t tofclusterUS = TMath::Abs(tofcluster);
741   
742   if(ecluster > 2)
743   {
744     if      (tofclusterUS < 25) fhEtaPhiEMCALBC0->Fill(etacluster, phicluster);
745     else if (tofclusterUS < 75) fhEtaPhiEMCALBC1->Fill(etacluster, phicluster);
746     else                        fhEtaPhiEMCALBCN->Fill(etacluster, phicluster);
747   }
748   
749   Int_t  bc     = GetReader()->GetTriggerClusterBC();
750   Int_t  id     = GetReader()->GetTriggerClusterId();
751   Bool_t badMax = GetReader()->IsBadMaxCellTriggerEvent();
752   
753   Int_t histoBC = bc+5;
754   if(GetReader()->AreBadTriggerEventsRemoved()) histoBC=0; // histograms created only for one BC since the others where rejected
755
756   if(id==-2)
757   {
758     //printf("AliAnaPhoton::ClusterSelected() - No trigger found bc=%d\n",bc);
759     fhEtaPhiNoTrigger->Fill(etacluster, phicluster);
760     fhTimeNoTrigger  ->Fill(ecluster, tofcluster);
761   }
762   else if(TMath::Abs(bc) < 6)
763   {
764     if(!GetReader()->IsBadCellTriggerEvent() && !GetReader()->IsExoticEvent() )
765     {
766       if(GetReader()->IsTriggerMatched())
767       {
768         if(ecluster > 2) fhEtaPhiTriggerEMCALBC[histoBC]->Fill(etacluster, phicluster);
769         fhTimeTriggerEMCALBC[histoBC]->Fill(ecluster, tofcluster);
770         if(GetReader()->IsPileUpFromSPD()) fhTimeTriggerEMCALBCPileUpSPD[histoBC]->Fill(ecluster, tofcluster);
771         
772         if(idcalo ==  GetReader()->GetTriggerClusterId())
773         {
774           fhEtaPhiTriggerEMCALBCCluster[histoBC]->Fill(etacluster, phicluster);
775           fhTimeTriggerEMCALBCCluster        ->Fill(ecluster, tofcluster);
776           
777           if(bc==0)
778           {
779             Float_t threshold = GetReader()->GetEventTriggerL1Threshold() ;
780             if(GetReader()->IsEventEMCALL0()) threshold = GetReader()->GetEventTriggerL0Threshold() ;
781             
782             if(ecluster > threshold)
783               fhEtaPhiTriggerEMCALBCClusterOverTh->Fill(etacluster, phicluster);
784             else if(ecluster > threshold-1)
785               fhEtaPhiTriggerEMCALBCClusterBelowTh1->Fill(etacluster, phicluster);
786             else
787               fhEtaPhiTriggerEMCALBCClusterBelowTh2->Fill(etacluster, phicluster);
788           }
789         }
790       }
791       else
792       {
793         if(ecluster > 2) fhEtaPhiTriggerEMCALBCUM[histoBC]->Fill(etacluster, phicluster);
794         fhTimeTriggerEMCALBCUM[histoBC]->Fill(ecluster, tofcluster);
795         
796         if(bc==0)
797         {
798           if(GetReader()->IsTriggerMatchedOpenCuts(0)) fhTimeTriggerEMCALBC0UMReMatchOpenTime   ->Fill(ecluster, tofcluster);
799           if(GetReader()->IsTriggerMatchedOpenCuts(1)) fhTimeTriggerEMCALBC0UMReMatchCheckNeigh ->Fill(ecluster, tofcluster);
800           if(GetReader()->IsTriggerMatchedOpenCuts(2)) fhTimeTriggerEMCALBC0UMReMatchBoth       ->Fill(ecluster, tofcluster);
801         }
802         
803         if(idcalo ==  GetReader()->GetTriggerClusterId())
804         {
805           fhEtaPhiTriggerEMCALBCUMCluster[histoBC]->Fill(etacluster, phicluster);
806           fhTimeTriggerEMCALBCUMCluster->Fill(ecluster, tofcluster);
807           if(bc==0)
808           {
809             Float_t threshold = GetReader()->GetEventTriggerL1Threshold() ;
810             if(GetReader()->IsEventEMCALL0()) threshold = GetReader()->GetEventTriggerL0Threshold() ;
811             
812             if(ecluster > threshold)
813               fhEtaPhiTriggerEMCALBCUMClusterOverTh->Fill(etacluster, phicluster);
814             else if(ecluster > threshold-1)
815               fhEtaPhiTriggerEMCALBCUMClusterBelowTh1->Fill(etacluster, phicluster);
816             else
817               fhEtaPhiTriggerEMCALBCUMClusterBelowTh2->Fill(etacluster, phicluster);
818             
819             if(GetReader()->IsTriggerMatchedOpenCuts(0))
820             {
821               fhEtaPhiTriggerEMCALBCUMReMatchOpenTimeCluster->Fill(etacluster, phicluster);
822               fhTimeTriggerEMCALBCUMReMatchOpenTimeCluster  ->Fill(ecluster, tofcluster);
823             }
824             if(GetReader()->IsTriggerMatchedOpenCuts(1))
825             {
826               fhEtaPhiTriggerEMCALBCUMReMatchCheckNeighCluster->Fill(etacluster, phicluster);
827               fhTimeTriggerEMCALBCUMReMatchCheckNeighCluster  ->Fill(ecluster, tofcluster);
828             }
829             if(GetReader()->IsTriggerMatchedOpenCuts(2))
830             {
831               fhEtaPhiTriggerEMCALBCUMReMatchBothCluster->Fill(etacluster, phicluster);
832               fhTimeTriggerEMCALBCUMReMatchBothCluster  ->Fill(ecluster, tofcluster);
833             }
834             
835           }
836         }
837       }
838     }// neither bad nor exotic
839     else if(GetReader()->IsBadCellTriggerEvent() && GetReader()->IsExoticEvent())
840     {
841       if(GetReader()->IsTriggerMatched())
842       {
843         if(ecluster > 2) fhEtaPhiTriggerEMCALBCBadExotic->Fill(etacluster, phicluster);
844         fhTimeTriggerEMCALBCBadExotic->Fill(ecluster, tofcluster);
845         if(badMax)  fhTimeTriggerEMCALBCBadMaxCellExotic->Fill(ecluster, tofcluster);
846       }
847       else
848       {
849         if(ecluster > 2) fhEtaPhiTriggerEMCALBCUMBadExotic->Fill(etacluster, phicluster);
850         fhTimeTriggerEMCALBCUMBadExotic->Fill(ecluster, tofcluster);
851         if(badMax)  fhTimeTriggerEMCALBCUMBadMaxCellExotic->Fill(ecluster, tofcluster);
852         
853       }
854     }// Bad and exotic cluster trigger
855     else if(GetReader()->IsBadCellTriggerEvent() )
856     {
857       if(GetReader()->IsTriggerMatched())
858       {
859         if(ecluster > 2) fhEtaPhiTriggerEMCALBCBad->Fill(etacluster, phicluster);
860         fhTimeTriggerEMCALBCBad->Fill(ecluster, tofcluster);
861         if(badMax)  fhTimeTriggerEMCALBCBadMaxCell->Fill(ecluster, tofcluster);
862       }
863       else
864       {
865         if(ecluster > 2) fhEtaPhiTriggerEMCALBCUMBad->Fill(etacluster, phicluster);
866         fhTimeTriggerEMCALBCUMBad->Fill(ecluster, tofcluster);
867         if(badMax)  fhTimeTriggerEMCALBCUMBadMaxCell->Fill(ecluster, tofcluster);
868       }
869     }// Bad cluster trigger
870     else if(GetReader()->IsExoticEvent() )
871     {
872       if(GetReader()->IsTriggerMatched())
873       {
874         if(ecluster > 2) fhEtaPhiTriggerEMCALBCExotic->Fill(etacluster, phicluster);
875         fhTimeTriggerEMCALBCExotic->Fill(ecluster, tofcluster);
876       }
877       else
878       {
879         if(ecluster > 2) fhEtaPhiTriggerEMCALBCUMExotic->Fill(etacluster, phicluster);
880         fhTimeTriggerEMCALBCUMExotic->Fill(ecluster, tofcluster);
881       }
882     }
883   }
884   else if(TMath::Abs(bc) >= 6)
885     printf("AliAnaPhoton::ClusterSelected() - Trigger BC not expected = %d\n",bc);
886   
887 }
888
889 //_________________________________________________________________________________________________________
890 void  AliAnaPhoton::FillClusterPileUpHistograms(AliVCluster * calo, Bool_t matched,     Float_t ptcluster,
891                                                 Float_t etacluster, Float_t phicluster, Float_t l0cluster)
892 {
893   // Fill some histograms related to pile up before any cluster cut is applied
894   
895   if(!fFillPileUpHistograms) return ;
896   
897   // Get the fraction of the cluster energy that carries the cell with highest energy and its absId
898   AliVCaloCells* cells = 0;
899   if(fCalorimeter == "EMCAL") cells = GetEMCALCells();
900   else                        cells = GetPHOSCells();
901   
902   Float_t maxCellFraction = 0.;
903   Int_t absIdMax = GetCaloUtils()->GetMaxEnergyCell(cells, calo,maxCellFraction);
904   
905   Double_t tmax  = cells->GetCellTime(absIdMax);
906   GetCaloUtils()->RecalibrateCellTime(tmax, fCalorimeter, absIdMax,GetReader()->GetInputEvent()->GetBunchCrossNumber());
907   tmax*=1.e9;
908   
909   Bool_t okPhoton = kFALSE;
910   if( GetCaloPID()->GetIdentifiedParticleType(calo)== AliCaloPID::kPhoton) okPhoton = kTRUE;
911   
912   Float_t clusterLongTimePt = 0;
913   Float_t clusterOKTimePt   = 0;
914     
915   //Loop on cells inside cluster, max cell must be over 100 MeV and time in BC=0
916   if(cells->GetCellAmplitude(absIdMax) > 0.1 && TMath::Abs(tmax) < 30)
917   {
918     for (Int_t ipos = 0; ipos < calo->GetNCells(); ipos++)
919     {
920       Int_t absId  = calo->GetCellsAbsId()[ipos];
921       
922       if( absId == absIdMax ) continue ;
923       
924       Double_t time  = cells->GetCellTime(absId);
925       Float_t  amp   = cells->GetCellAmplitude(absId);
926       Int_t    bc    = GetReader()->GetInputEvent()->GetBunchCrossNumber();
927       GetCaloUtils()->GetEMCALRecoUtils()->AcceptCalibrateCell(absId,bc,amp,time,cells);
928       time*=1e9;
929       
930       Float_t diff = (tmax-time);
931       
932       if(GetReader()->IsInTimeWindow(time,amp)) clusterOKTimePt   += amp;
933       else                                      clusterLongTimePt += amp;
934       
935       if( cells->GetCellAmplitude(absIdMax) < 0.1 ) continue ;
936       
937       if(GetReader()->IsPileUpFromSPD())
938       {
939         fhClusterCellTimePileUp[0]->Fill(ptcluster, time);
940         fhClusterTimeDiffPileUp[0]->Fill(ptcluster, diff);
941         if(!matched)
942         {
943           fhClusterTimeDiffChargedPileUp[0]->Fill(ptcluster, diff);
944           if(okPhoton)  fhClusterTimeDiffPhotonPileUp[0]->Fill(ptcluster, diff);
945         }
946       }
947       
948       if(GetReader()->IsPileUpFromEMCal())
949       {
950         fhClusterCellTimePileUp[1]->Fill(ptcluster, time);
951         fhClusterTimeDiffPileUp[1]->Fill(ptcluster, diff);
952         if(!matched)
953         {
954           fhClusterTimeDiffChargedPileUp[1]->Fill(ptcluster, diff);
955           if(okPhoton)  fhClusterTimeDiffPhotonPileUp[1]->Fill(ptcluster, diff);
956         }
957       }
958       
959       if(GetReader()->IsPileUpFromSPDOrEMCal())
960       {
961         fhClusterCellTimePileUp[2]->Fill(ptcluster, time);
962         fhClusterTimeDiffPileUp[2]->Fill(ptcluster, diff);
963         if(!matched)
964         {
965           fhClusterTimeDiffChargedPileUp[2]->Fill(ptcluster, diff);
966           if(okPhoton)  fhClusterTimeDiffPhotonPileUp[2]->Fill(ptcluster, diff);
967         }
968       }
969       
970       if(GetReader()->IsPileUpFromSPDAndEMCal())
971       {
972         fhClusterCellTimePileUp[3]->Fill(ptcluster, time);
973         fhClusterTimeDiffPileUp[3]->Fill(ptcluster, diff);
974         if(!matched)
975         {
976           fhClusterTimeDiffChargedPileUp[3]->Fill(ptcluster, diff);
977           if(okPhoton)  fhClusterTimeDiffPhotonPileUp[3]->Fill(ptcluster, diff);
978         }
979       }
980       
981       if(GetReader()->IsPileUpFromSPDAndNotEMCal())
982       {
983         fhClusterCellTimePileUp[4]->Fill(ptcluster, time);
984         fhClusterTimeDiffPileUp[4]->Fill(ptcluster, diff);
985         if(!matched)
986         {
987           fhClusterTimeDiffChargedPileUp[4]->Fill(ptcluster, diff);
988           if(okPhoton)  fhClusterTimeDiffPhotonPileUp[4]->Fill(ptcluster, diff);
989         }
990       }
991       
992       if(GetReader()->IsPileUpFromEMCalAndNotSPD())
993       {
994         fhClusterCellTimePileUp[5]->Fill(ptcluster, time);
995         fhClusterTimeDiffPileUp[5]->Fill(ptcluster, diff);
996         if(!matched)
997         {
998           fhClusterTimeDiffChargedPileUp[5]->Fill(ptcluster, diff);
999           if(okPhoton)  fhClusterTimeDiffPhotonPileUp[5]->Fill(ptcluster, diff);
1000         }
1001       }
1002       
1003       if(GetReader()->IsPileUpFromNotSPDAndNotEMCal())
1004       {
1005         fhClusterCellTimePileUp[6]->Fill(ptcluster, time);
1006         fhClusterTimeDiffPileUp[6]->Fill(ptcluster, diff);
1007         if(!matched)
1008         {
1009           fhClusterTimeDiffChargedPileUp[6]->Fill(ptcluster, diff);
1010           if(okPhoton)  fhClusterTimeDiffPhotonPileUp[6]->Fill(ptcluster, diff);
1011         }
1012       }
1013     }//loop
1014     
1015     
1016     Float_t frac = 0;
1017     if(clusterLongTimePt+clusterOKTimePt > 0.001)
1018       frac = clusterLongTimePt/(clusterLongTimePt+clusterOKTimePt);
1019     //printf("E long %f, E OK %f, Fraction large time %f, E %f\n",clusterLongTimePt,clusterOKTimePt,frac,ptcluster);
1020     
1021     if(GetReader()->IsPileUpFromSPD())               {fhPtPileUp[0]->Fill(ptcluster); fhLambda0PileUp[0]->Fill(ptcluster,l0cluster); fhClusterEFracLongTimePileUp[0]->Fill(ptcluster,frac);}
1022     if(GetReader()->IsPileUpFromEMCal())             {fhPtPileUp[1]->Fill(ptcluster); fhLambda0PileUp[1]->Fill(ptcluster,l0cluster); fhClusterEFracLongTimePileUp[1]->Fill(ptcluster,frac);}
1023     if(GetReader()->IsPileUpFromSPDOrEMCal())        {fhPtPileUp[2]->Fill(ptcluster); fhLambda0PileUp[2]->Fill(ptcluster,l0cluster); fhClusterEFracLongTimePileUp[2]->Fill(ptcluster,frac);}
1024     if(GetReader()->IsPileUpFromSPDAndEMCal())       {fhPtPileUp[3]->Fill(ptcluster); fhLambda0PileUp[3]->Fill(ptcluster,l0cluster); fhClusterEFracLongTimePileUp[3]->Fill(ptcluster,frac);}
1025     if(GetReader()->IsPileUpFromSPDAndNotEMCal())    {fhPtPileUp[4]->Fill(ptcluster); fhLambda0PileUp[4]->Fill(ptcluster,l0cluster); fhClusterEFracLongTimePileUp[4]->Fill(ptcluster,frac);}
1026     if(GetReader()->IsPileUpFromEMCalAndNotSPD())    {fhPtPileUp[5]->Fill(ptcluster); fhLambda0PileUp[5]->Fill(ptcluster,l0cluster); fhClusterEFracLongTimePileUp[5]->Fill(ptcluster,frac);}
1027     if(GetReader()->IsPileUpFromNotSPDAndNotEMCal()) {fhPtPileUp[6]->Fill(ptcluster); fhLambda0PileUp[6]->Fill(ptcluster,l0cluster); fhClusterEFracLongTimePileUp[6]->Fill(ptcluster,frac);}
1028     
1029     fhEtaPhiBC0->Fill(etacluster,phicluster);
1030     if(GetReader()->IsPileUpFromSPD()) fhEtaPhiBC0PileUpSPD    ->Fill(etacluster,phicluster);
1031   }
1032   
1033   else if (tmax > 25)         {fhEtaPhiBCPlus ->Fill(etacluster,phicluster); if(GetReader()->IsPileUpFromSPD()) fhEtaPhiBCPlusPileUpSPD ->Fill(etacluster,phicluster); }
1034   else if (tmax <-25)         {fhEtaPhiBCMinus->Fill(etacluster,phicluster); if(GetReader()->IsPileUpFromSPD()) fhEtaPhiBCMinusPileUpSPD->Fill(etacluster,phicluster); }
1035 }
1036
1037 //_______________________________________________
1038 void AliAnaPhoton::FillPileUpHistogramsPerEvent()
1039 {
1040   // Fill some histograms per event to understand pile-up
1041   // Open the time cut in the reader to be more meaningful
1042   
1043   if(!fFillPileUpHistograms) return;
1044   
1045   AliVEvent * event = GetReader()->GetInputEvent();
1046         
1047   AliESDEvent* esdEv = dynamic_cast<AliESDEvent*> (event);
1048   AliAODEvent* aodEv = dynamic_cast<AliAODEvent*> (event);
1049         
1050   // N pile up vertices
1051   Int_t nVtxSPD = -1;
1052   Int_t nVtxTrk = -1;
1053   TLorentzVector mom;
1054         
1055   if      (esdEv)
1056   {
1057                 nVtxSPD = esdEv->GetNumberOfPileupVerticesSPD();
1058                 nVtxTrk = esdEv->GetNumberOfPileupVerticesTracks();
1059   }//ESD
1060   else if (aodEv)
1061   {
1062                 nVtxSPD = aodEv->GetNumberOfPileupVerticesSPD();
1063                 nVtxTrk = aodEv->GetNumberOfPileupVerticesTracks();
1064   }//AOD
1065         
1066         
1067         // Get the appropriate list of clusters
1068         TClonesArray * clusterList = 0;
1069         TString  clusterListName   = GetReader()->GetEMCALClusterListName();
1070         if     (event->FindListObject(clusterListName))
1071                 clusterList = dynamic_cast<TClonesArray*> (event->FindListObject(clusterListName));
1072         else if(GetReader()->GetOutputEvent())
1073                 clusterList = dynamic_cast<TClonesArray*> (GetReader()->GetOutputEvent()->FindListObject(clusterListName));
1074         
1075         // Loop on clusters, get the maximum energy cluster as reference
1076   Int_t nclusters = 0;
1077         if(clusterList) nclusters = clusterList->GetEntriesFast();
1078         else            nclusters = event->GetNumberOfCaloClusters();
1079         
1080   Int_t   idMax = 0;
1081   Float_t  eMax = 0;
1082   Float_t  tMax = 0;
1083   for(Int_t iclus = 0; iclus < nclusters ; iclus++)
1084   {
1085                 AliVCluster * clus = 0;
1086                 if(clusterList) clus = (AliVCluster*) (clusterList->At(iclus));
1087     else            clus = GetReader()->GetInputEvent()->GetCaloCluster(iclus);
1088                 
1089                 if(!clus)            continue;
1090                 
1091                 if(!clus->IsEMCAL()) continue;
1092                 
1093                 Float_t tof = clus->GetTOF()*1e9;
1094                 if(clus->E() > eMax && TMath::Abs(tof) < 30)
1095     {
1096       eMax  = clus->E();
1097                         tMax  = tof;
1098       idMax = iclus;
1099     }
1100           
1101                 clus->GetMomentum(mom,GetVertex(0));
1102                 Float_t pt = mom.Pt();
1103           
1104                 fhPtNPileUpSPDVtx->Fill(pt,nVtxSPD);
1105                 fhPtNPileUpTrkVtx->Fill(pt,nVtxTrk);
1106     
1107                 if(TMath::Abs(tof) < 30)
1108                 {
1109                         fhPtNPileUpSPDVtxTimeCut->Fill(pt,nVtxSPD);
1110                         fhPtNPileUpTrkVtxTimeCut->Fill(pt,nVtxTrk);
1111                 }
1112     
1113     if(tof < 75 && tof > -30)
1114     {
1115       fhPtNPileUpSPDVtxTimeCut2->Fill(pt,nVtxSPD);
1116       fhPtNPileUpTrkVtxTimeCut2->Fill(pt,nVtxTrk);
1117     }
1118     
1119     fhTimePtNoCut->Fill(pt,tof);
1120     if(GetReader()->IsPileUpFromSPD()) fhTimePtSPD->Fill(pt,tof);
1121
1122   }
1123         
1124   if(eMax < 5) return;
1125   
1126   // Loop again on clusters to compare this max cluster t and the rest of the clusters, if E > 0.3
1127   Int_t n20  = 0;
1128   Int_t n40  = 0;
1129   Int_t n    = 0;
1130   Int_t nOK  = 0;
1131   
1132   for(Int_t iclus = 0; iclus < nclusters ; iclus++)
1133   {
1134                 AliVCluster * clus = 0;
1135                 if(clusterList) clus = (AliVCluster*) (clusterList->At(iclus));
1136     else            clus = GetReader()->GetInputEvent()->GetCaloCluster(iclus);
1137                 
1138                 if(!clus)            continue;
1139                 
1140                 if(!clus->IsEMCAL()) continue;
1141     
1142     if(clus->E() < 0.3 || iclus==idMax) continue;
1143     
1144     Float_t tdiff = TMath::Abs(tMax-clus->GetTOF()*1e9);
1145     n++;
1146     if(tdiff < 25) nOK++;
1147     else
1148     {
1149       n20++;
1150       if(tdiff > 40 ) n40++;
1151     }
1152   }
1153   
1154   // Check pile-up and fill histograms depending on the different cluster multiplicities
1155   if(GetReader()->IsPileUpFromSPD())
1156   {
1157     fhClusterMultSPDPileUp[0]->Fill(eMax,n  );
1158     fhClusterMultSPDPileUp[1]->Fill(eMax,nOK);
1159     fhClusterMultSPDPileUp[2]->Fill(eMax,n20);
1160     fhClusterMultSPDPileUp[3]->Fill(eMax,n40);
1161   }
1162   else
1163   {
1164     fhClusterMultNoPileUp[0]->Fill(eMax,n  );
1165     fhClusterMultNoPileUp[1]->Fill(eMax,nOK);
1166     fhClusterMultNoPileUp[2]->Fill(eMax,n20);
1167     fhClusterMultNoPileUp[3]->Fill(eMax,n40);
1168   }
1169   
1170 }
1171
1172
1173 //_________________________________________________________________________________________________
1174 void AliAnaPhoton::FillPileUpHistograms(Float_t energy, Float_t pt, Float_t time)
1175 {
1176   // Fill some histograms to understand pile-up
1177   if(!fFillPileUpHistograms) return;
1178   
1179   //printf("E %f, time %f\n",energy,time);
1180   AliVEvent * event = GetReader()->GetInputEvent();
1181   
1182   if(GetReader()->IsPileUpFromSPD())               fhPtPhotonPileUp[0]->Fill(pt);
1183   if(GetReader()->IsPileUpFromEMCal())             fhPtPhotonPileUp[1]->Fill(pt);
1184   if(GetReader()->IsPileUpFromSPDOrEMCal())        fhPtPhotonPileUp[2]->Fill(pt);
1185   if(GetReader()->IsPileUpFromSPDAndEMCal())       fhPtPhotonPileUp[3]->Fill(pt);
1186   if(GetReader()->IsPileUpFromSPDAndNotEMCal())    fhPtPhotonPileUp[4]->Fill(pt);
1187   if(GetReader()->IsPileUpFromEMCalAndNotSPD())    fhPtPhotonPileUp[5]->Fill(pt);
1188   if(GetReader()->IsPileUpFromNotSPDAndNotEMCal()) fhPtPhotonPileUp[6]->Fill(pt);
1189   
1190   fhTimePtPhotonNoCut->Fill(pt,time);
1191   if(GetReader()->IsPileUpFromSPD()) fhTimePtPhotonSPD->Fill(pt,time);
1192   
1193   if(energy < 8) return; // Fill time figures for high energy clusters not too close to trigger threshold
1194   
1195   AliESDEvent* esdEv = dynamic_cast<AliESDEvent*> (event);
1196   AliAODEvent* aodEv = dynamic_cast<AliAODEvent*> (event);
1197   
1198   // N pile up vertices
1199   Int_t nVtxSPD = -1;
1200   Int_t nVtxTrk = -1;
1201   
1202   if      (esdEv)
1203   {
1204     nVtxSPD = esdEv->GetNumberOfPileupVerticesSPD();
1205     nVtxTrk = esdEv->GetNumberOfPileupVerticesTracks();
1206     
1207   }//ESD
1208   else if (aodEv)
1209   {
1210     nVtxSPD = aodEv->GetNumberOfPileupVerticesSPD();
1211     nVtxTrk = aodEv->GetNumberOfPileupVerticesTracks();
1212   }//AOD
1213   
1214   fhTimeNPileUpVertSPD  ->Fill(time,nVtxSPD);
1215   fhTimeNPileUpVertTrack->Fill(time,nVtxTrk);
1216   
1217   fhPtPhotonNPileUpSPDVtx->Fill(pt,nVtxSPD);
1218   fhPtPhotonNPileUpTrkVtx->Fill(pt,nVtxTrk);
1219         
1220   if(TMath::Abs(time) < 25)
1221   {
1222     fhPtPhotonNPileUpSPDVtxTimeCut->Fill(pt,nVtxSPD);
1223     fhPtPhotonNPileUpTrkVtxTimeCut->Fill(pt,nVtxTrk);
1224   }
1225         
1226   if(time < 75 && time > -25)
1227   {
1228     fhPtPhotonNPileUpSPDVtxTimeCut2->Fill(pt,nVtxSPD);
1229     fhPtPhotonNPileUpTrkVtxTimeCut2->Fill(pt,nVtxTrk);
1230   }
1231   
1232   //printf("Is SPD %d, Is SPD Multi %d, n spd %d, n track %d\n",
1233   //       GetReader()->IsPileUpFromSPD(),event->IsPileupFromSPDInMultBins(),nVtxSPD,nVtxTrk);
1234   
1235   Int_t ncont = -1;
1236   Float_t z1 = -1, z2 = -1;
1237   Float_t diamZ = -1;
1238   for(Int_t iVert=0; iVert<nVtxSPD;iVert++)
1239   {
1240     if      (esdEv)
1241     {
1242       const AliESDVertex* pv=esdEv->GetPileupVertexSPD(iVert);
1243       ncont=pv->GetNContributors();
1244       z1 = esdEv->GetPrimaryVertexSPD()->GetZ();
1245       z2 = pv->GetZ();
1246       diamZ = esdEv->GetDiamondZ();
1247     }//ESD
1248     else if (aodEv)
1249     {
1250       AliAODVertex *pv=aodEv->GetVertex(iVert);
1251       if(pv->GetType()!=AliAODVertex::kPileupSPD) continue;
1252       ncont=pv->GetNContributors();
1253       z1=aodEv->GetPrimaryVertexSPD()->GetZ();
1254       z2=pv->GetZ();
1255       diamZ = aodEv->GetDiamondZ();
1256     }// AOD
1257     
1258     Double_t distZ  = TMath::Abs(z2-z1);
1259     diamZ  = TMath::Abs(z2-diamZ);
1260     
1261     fhTimeNPileUpVertContributors  ->Fill(time,ncont);
1262     fhTimePileUpMainVertexZDistance->Fill(time,distZ);
1263     fhTimePileUpMainVertexZDiamond ->Fill(time,diamZ);
1264     
1265   }// loop
1266 }
1267
1268 //____________________________________________________________________________________
1269 void  AliAnaPhoton::FillShowerShapeHistograms(AliVCluster* cluster, Int_t mcTag)
1270 {
1271   //Fill cluster Shower Shape histograms
1272   
1273   if(!fFillSSHistograms || GetMixedEvent()) return;
1274   
1275   Float_t energy  = cluster->E();
1276   Int_t   ncells  = cluster->GetNCells();
1277   Float_t lambda0 = cluster->GetM02();
1278   Float_t lambda1 = cluster->GetM20();
1279   Float_t disp    = cluster->GetDispersion()*cluster->GetDispersion();
1280   
1281   TLorentzVector mom;
1282   if(GetReader()->GetDataType() != AliCaloTrackReader::kMC)
1283   {
1284     cluster->GetMomentum(mom,GetVertex(0)) ;
1285   }//Assume that come from vertex in straight line
1286   else
1287   {
1288     Double_t vertex[]={0,0,0};
1289     cluster->GetMomentum(mom,vertex) ;
1290   }
1291   
1292   Float_t eta = mom.Eta();
1293   Float_t phi = mom.Phi();
1294   if(phi < 0) phi+=TMath::TwoPi();
1295   
1296   fhLam0E ->Fill(energy,lambda0);
1297   fhLam1E ->Fill(energy,lambda1);
1298   fhDispE ->Fill(energy,disp);
1299   
1300   if(fCalorimeter == "EMCAL" && GetModuleNumber(cluster) > 5)
1301   {
1302     fhLam0ETRD->Fill(energy,lambda0);
1303     fhLam1ETRD->Fill(energy,lambda1);
1304     fhDispETRD->Fill(energy,disp);
1305   }
1306   
1307   Float_t l0   = 0., l1   = 0.;
1308   Float_t dispp= 0., dEta = 0., dPhi    = 0.;
1309   Float_t sEta = 0., sPhi = 0., sEtaPhi = 0.;
1310   if(fCalorimeter == "EMCAL" && !fFillOnlySimpleSSHisto)
1311   {
1312     GetCaloUtils()->GetEMCALRecoUtils()->RecalculateClusterShowerShapeParameters(GetEMCALGeometry(), GetReader()->GetInputEvent()->GetEMCALCells(), cluster,
1313                                                                                  l0, l1, dispp, dEta, dPhi, sEta, sPhi, sEtaPhi);
1314     //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",
1315     //       l0, l1, dispp, dEta, dPhi, sEta, sPhi, sEtaPhi );
1316     //printf("AliAnaPhoton::FillShowerShapeHistogram - dispersion %f, dispersion eta+phi %f \n",
1317     //       disp, dPhi+dEta );
1318     fhDispEtaE        -> Fill(energy,dEta);
1319     fhDispPhiE        -> Fill(energy,dPhi);
1320     fhSumEtaE         -> Fill(energy,sEta);
1321     fhSumPhiE         -> Fill(energy,sPhi);
1322     fhSumEtaPhiE      -> Fill(energy,sEtaPhi);
1323     fhDispEtaPhiDiffE -> Fill(energy,dPhi-dEta);
1324     if(dEta+dPhi>0)fhSphericityE     -> Fill(energy,(dPhi-dEta)/(dEta+dPhi));
1325     if(dEta+sEta>0)fhDispSumEtaDiffE -> Fill(energy,(dEta-sEta)/((dEta+sEta)/2.));
1326     if(dPhi+sPhi>0)fhDispSumPhiDiffE -> Fill(energy,(dPhi-sPhi)/((dPhi+sPhi)/2.));
1327     
1328     Int_t ebin = -1;
1329     if      (energy < 2 ) ebin = 0;
1330     else if (energy < 4 ) ebin = 1;
1331     else if (energy < 6 ) ebin = 2;
1332     else if (energy < 10) ebin = 3;
1333     else if (energy < 15) ebin = 4;
1334     else if (energy < 20) ebin = 5;
1335     else                  ebin = 6;
1336     
1337     fhDispEtaDispPhi[ebin]->Fill(dEta   ,dPhi);
1338     fhLambda0DispEta[ebin]->Fill(lambda0,dEta);
1339     fhLambda0DispPhi[ebin]->Fill(lambda0,dPhi);
1340     
1341   }
1342   
1343   // if track-matching was of, check effect of track-matching residual cut
1344   
1345   if(!fRejectTrackMatch)
1346   {
1347     Float_t dZ  = cluster->GetTrackDz();
1348     Float_t dR  = cluster->GetTrackDx();
1349     if(cluster->IsEMCAL() && GetCaloUtils()->IsRecalculationOfClusterTrackMatchingOn())
1350     {
1351       dR = 2000., dZ = 2000.;
1352       GetCaloUtils()->GetEMCALRecoUtils()->GetMatchedResiduals(cluster->GetID(),dZ,dR);
1353     }
1354     
1355     if(TMath::Abs(dZ) < 0.05 && TMath::Abs(dR) < 0.05)
1356     {
1357       fhLam0ETM ->Fill(energy,lambda0);
1358       fhLam1ETM ->Fill(energy,lambda1);
1359       fhDispETM ->Fill(energy,disp);
1360       
1361       if(fCalorimeter == "EMCAL" && GetModuleNumber(cluster) > 5)
1362       {
1363         fhLam0ETMTRD->Fill(energy,lambda0);
1364         fhLam1ETMTRD->Fill(energy,lambda1);
1365         fhDispETMTRD->Fill(energy,disp);
1366       }
1367     }
1368   }// if track-matching was of, check effect of matching residual cut
1369   
1370   
1371   if(!fFillOnlySimpleSSHisto)
1372   {
1373     if(energy < 2)
1374     {
1375       fhNCellsLam0LowE ->Fill(ncells,lambda0);
1376       fhNCellsLam1LowE ->Fill(ncells,lambda1);
1377       fhNCellsDispLowE ->Fill(ncells,disp);
1378       
1379       fhLam1Lam0LowE  ->Fill(lambda1,lambda0);
1380       fhLam0DispLowE  ->Fill(lambda0,disp);
1381       fhDispLam1LowE  ->Fill(disp,lambda1);
1382       fhEtaLam0LowE   ->Fill(eta,lambda0);
1383       fhPhiLam0LowE   ->Fill(phi,lambda0);
1384     }
1385     else
1386     {
1387       fhNCellsLam0HighE ->Fill(ncells,lambda0);
1388       fhNCellsLam1HighE ->Fill(ncells,lambda1);
1389       fhNCellsDispHighE ->Fill(ncells,disp);
1390       
1391       fhLam1Lam0HighE  ->Fill(lambda1,lambda0);
1392       fhLam0DispHighE  ->Fill(lambda0,disp);
1393       fhDispLam1HighE  ->Fill(disp,lambda1);
1394       fhEtaLam0HighE   ->Fill(eta, lambda0);
1395       fhPhiLam0HighE   ->Fill(phi, lambda0);
1396     }
1397   }
1398   
1399   if(IsDataMC())
1400   {
1401     AliVCaloCells* cells = 0;
1402     if(fCalorimeter == "EMCAL") cells = GetEMCALCells();
1403     else                        cells = GetPHOSCells();
1404     
1405     //Fill histograms to check shape of embedded clusters
1406     Float_t fraction = 0;
1407     // printf("check embedding %i\n",GetReader()->IsEmbeddedClusterSelectionOn());
1408     
1409     if(GetReader()->IsEmbeddedClusterSelectionOn())
1410     {//Only working for EMCAL
1411       //        printf("embedded\n");
1412       Float_t clusterE = 0; // recalculate in case corrections applied.
1413       Float_t cellE    = 0;
1414       for(Int_t icell  = 0; icell < cluster->GetNCells(); icell++)
1415       {
1416         cellE    = cells->GetCellAmplitude(cluster->GetCellAbsId(icell));
1417         clusterE+=cellE;
1418         fraction+=cellE*cluster->GetCellAmplitudeFraction(icell);
1419       }
1420       
1421       //Fraction of total energy due to the embedded signal
1422       fraction/=clusterE;
1423       
1424       if(GetDebug() > 1 )
1425         printf("AliAnaPhoton::FillShowerShapeHistogram() - Energy fraction of embedded signal %2.3f, Energy %2.3f\n",fraction, clusterE);
1426       
1427       fhEmbeddedSignalFractionEnergy->Fill(clusterE,fraction);
1428       
1429     }  // embedded fraction
1430     
1431     // Get the fraction of the cluster energy that carries the cell with highest energy
1432     Float_t maxCellFraction = 0.;
1433     Int_t absID = GetCaloUtils()->GetMaxEnergyCell(cells, cluster,maxCellFraction);
1434     
1435     if( absID < 0 ) AliFatal("Wrong absID");
1436       
1437     // Check the origin and fill histograms
1438     
1439     Int_t mcIndex = -1;
1440     
1441     if( GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPhoton)     &&
1442        !GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCConversion) &&
1443        !GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPi0)        &&
1444        !GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCEta))
1445     {
1446       mcIndex = kmcssPhoton ;
1447       
1448       if(!GetReader()->IsEmbeddedClusterSelectionOn())
1449       {
1450         //Check particle overlaps in cluster
1451         
1452         // Compare the primary depositing more energy with the rest,
1453         // if no photon/electron as comon ancestor (conversions), count as other particle
1454         const UInt_t nlabels = cluster->GetNLabels();
1455         Int_t overpdg[nlabels];
1456         Int_t noverlaps = GetMCAnalysisUtils()->GetNOverlaps(cluster->GetLabels(), nlabels,mcTag,-1,GetReader(),overpdg);
1457
1458         //printf("N overlaps %d \n",noverlaps);
1459         
1460         if(noverlaps == 0)
1461         {
1462           fhMCPhotonELambda0NoOverlap  ->Fill(energy, lambda0);
1463         }
1464         else if(noverlaps == 1)
1465         {
1466           fhMCPhotonELambda0TwoOverlap ->Fill(energy, lambda0);
1467         }
1468         else if(noverlaps > 1)
1469         {
1470           fhMCPhotonELambda0NOverlap   ->Fill(energy, lambda0);
1471         }
1472         else
1473         {
1474           printf("AliAnaPhoton::FillShowerShapeHistogram() - n overlaps = %d!!\n", noverlaps);
1475         }
1476       }//No embedding
1477       
1478       //Fill histograms to check shape of embedded clusters
1479       if(GetReader()->IsEmbeddedClusterSelectionOn())
1480       {
1481         if     (fraction > 0.9)
1482         {
1483           fhEmbedPhotonELambda0FullSignal   ->Fill(energy, lambda0);
1484         }
1485         else if(fraction > 0.5)
1486         {
1487           fhEmbedPhotonELambda0MostlySignal ->Fill(energy, lambda0);
1488         }
1489         else if(fraction > 0.1)
1490         {
1491           fhEmbedPhotonELambda0MostlyBkg    ->Fill(energy, lambda0);
1492         }
1493         else
1494         {
1495           fhEmbedPhotonELambda0FullBkg      ->Fill(energy, lambda0);
1496         }
1497       } // embedded
1498       
1499     }//photon   no conversion
1500     else if  ( GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPhoton)     &&
1501                GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCConversion) &&
1502               !GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPi0)        &&
1503               !GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCEta))
1504     {
1505       mcIndex = kmcssConversion ;
1506     }//conversion photon
1507
1508     else if  ( GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCElectron))
1509     {
1510       mcIndex = kmcssElectron ;
1511     }//electron
1512     else if  ( GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPi0)  )
1513     {
1514       mcIndex = kmcssPi0 ;
1515       
1516       //Fill histograms to check shape of embedded clusters
1517       if(GetReader()->IsEmbeddedClusterSelectionOn())
1518       {
1519         if     (fraction > 0.9)
1520         {
1521           fhEmbedPi0ELambda0FullSignal   ->Fill(energy, lambda0);
1522         }
1523         else if(fraction > 0.5)
1524         {
1525           fhEmbedPi0ELambda0MostlySignal ->Fill(energy, lambda0);
1526         }
1527         else if(fraction > 0.1)
1528         {
1529           fhEmbedPi0ELambda0MostlyBkg    ->Fill(energy, lambda0);
1530         }
1531         else
1532         {
1533           fhEmbedPi0ELambda0FullBkg      ->Fill(energy, lambda0);
1534         }
1535       } // embedded
1536       
1537     }//pi0
1538     else if  ( GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCEta)  )
1539     {
1540       mcIndex = kmcssEta ;
1541     }//eta
1542     else
1543     {
1544       mcIndex = kmcssOther ;
1545     }//other particles
1546     
1547     fhMCELambda0           [mcIndex]->Fill(energy, lambda0);
1548     fhMCELambda1           [mcIndex]->Fill(energy, lambda1);
1549     fhMCEDispersion        [mcIndex]->Fill(energy, disp);
1550     fhMCNCellsE            [mcIndex]->Fill(energy, ncells);
1551     fhMCMaxCellDiffClusterE[mcIndex]->Fill(energy, maxCellFraction);
1552     
1553     if(!fFillOnlySimpleSSHisto)
1554     {
1555       if     (energy < 2.)
1556       {
1557         fhMCLambda0vsClusterMaxCellDiffE0[mcIndex]->Fill(lambda0, maxCellFraction);
1558         fhMCNCellsvsClusterMaxCellDiffE0 [mcIndex]->Fill(ncells,  maxCellFraction);
1559       }
1560       else if(energy < 6.)
1561       {
1562         fhMCLambda0vsClusterMaxCellDiffE2[mcIndex]->Fill(lambda0, maxCellFraction);
1563         fhMCNCellsvsClusterMaxCellDiffE2 [mcIndex]->Fill(ncells,  maxCellFraction);
1564       }
1565       else
1566       {
1567         fhMCLambda0vsClusterMaxCellDiffE6[mcIndex]->Fill(lambda0, maxCellFraction);
1568         fhMCNCellsvsClusterMaxCellDiffE6 [mcIndex]->Fill(ncells,  maxCellFraction);
1569       }
1570       
1571       if(fCalorimeter == "EMCAL")
1572       {
1573         fhMCEDispEta        [mcIndex]-> Fill(energy,dEta);
1574         fhMCEDispPhi        [mcIndex]-> Fill(energy,dPhi);
1575         fhMCESumEtaPhi      [mcIndex]-> Fill(energy,sEtaPhi);
1576         fhMCEDispEtaPhiDiff [mcIndex]-> Fill(energy,dPhi-dEta);
1577         if(dEta+dPhi>0)fhMCESphericity[mcIndex]-> Fill(energy,(dPhi-dEta)/(dEta+dPhi));
1578         
1579         Int_t ebin = -1;
1580         if      (energy < 2 ) ebin = 0;
1581         else if (energy < 4 ) ebin = 1;
1582         else if (energy < 6 ) ebin = 2;
1583         else if (energy < 10) ebin = 3;
1584         else if (energy < 15) ebin = 4;
1585         else if (energy < 20) ebin = 5;
1586         else                  ebin = 6;
1587         
1588         fhMCDispEtaDispPhi[ebin][mcIndex]->Fill(dEta   ,dPhi);
1589         fhMCLambda0DispEta[ebin][mcIndex]->Fill(lambda0,dEta);
1590         fhMCLambda0DispPhi[ebin][mcIndex]->Fill(lambda0,dPhi);
1591       }
1592     }
1593   }//MC data
1594   
1595 }
1596
1597 //__________________________________________________________________________
1598 void AliAnaPhoton::FillTrackMatchingResidualHistograms(AliVCluster* cluster,
1599                                                        Int_t cut)
1600 {
1601   // If selected, fill histograms with residuals of matched clusters, help to define track matching cut
1602   // Residual filled for different cuts 0 (No cut), after 1 PID cut
1603   
1604   Float_t dZ = cluster->GetTrackDz();
1605   Float_t dR = cluster->GetTrackDx();
1606   
1607   if(cluster->IsEMCAL() && GetCaloUtils()->IsRecalculationOfClusterTrackMatchingOn())
1608   {
1609     dR = 2000., dZ = 2000.;
1610     GetCaloUtils()->GetEMCALRecoUtils()->GetMatchedResiduals(cluster->GetID(),dZ,dR);
1611   }
1612   
1613   AliVTrack *track = GetCaloUtils()->GetMatchedTrack(cluster, GetReader()->GetInputEvent());
1614   
1615   Bool_t positive = kFALSE;
1616   if(track) positive = (track->Charge()>0);
1617   
1618   if(fhTrackMatchedDEta[cut] && TMath::Abs(dR) < 999)
1619   {
1620     fhTrackMatchedDEta[cut]->Fill(cluster->E(),dZ);
1621     fhTrackMatchedDPhi[cut]->Fill(cluster->E(),dR);
1622     if(cluster->E() > 0.5) fhTrackMatchedDEtaDPhi[cut]->Fill(dZ,dR);
1623
1624     if(track)
1625     {
1626       if(positive)
1627       {
1628         fhTrackMatchedDEtaPos[cut]->Fill(cluster->E(),dZ);
1629         fhTrackMatchedDPhiPos[cut]->Fill(cluster->E(),dR);
1630         if(cluster->E() > 0.5) fhTrackMatchedDEtaDPhiPos[cut]->Fill(dZ,dR);
1631       }
1632       else
1633       {
1634         fhTrackMatchedDEtaNeg[cut]->Fill(cluster->E(),dZ);
1635         fhTrackMatchedDPhiNeg[cut]->Fill(cluster->E(),dR);
1636         if(cluster->E() > 0.5) fhTrackMatchedDEtaDPhiNeg[cut]->Fill(dZ,dR);
1637       }
1638     }
1639     
1640     Int_t nSMod = GetModuleNumber(cluster);
1641     
1642     if(fCalorimeter=="EMCAL" &&  nSMod > 5)
1643     {
1644       fhTrackMatchedDEtaTRD[cut]->Fill(cluster->E(),dZ);
1645       fhTrackMatchedDPhiTRD[cut]->Fill(cluster->E(),dR);
1646     }
1647     
1648     // Check dEdx and E/p of matched clusters
1649
1650     if(TMath::Abs(dZ) < 0.05 && TMath::Abs(dR) < 0.05)
1651     {
1652       if(track)
1653       {
1654         Float_t dEdx   = track->GetTPCsignal();
1655         Float_t eOverp = cluster->E()/track->P();
1656         
1657         fhdEdx[cut]  ->Fill(cluster->E(), dEdx);
1658         fhEOverP[cut]->Fill(cluster->E(), eOverp);
1659         
1660         if(fCalorimeter=="EMCAL" && nSMod > 5)
1661           fhEOverPTRD[cut]->Fill(cluster->E(), eOverp);
1662         
1663         
1664       }
1665       else
1666         printf("AliAnaPhoton::FillTrackMatchingResidualHistograms() - Residual OK but (dR, dZ)= (%2.4f,%2.4f) no track associated WHAT? \n", dR,dZ);
1667       
1668       
1669       
1670       if(IsDataMC())
1671       {
1672         
1673         Int_t tag = GetMCAnalysisUtils()->CheckOrigin(cluster->GetLabels(),cluster->GetNLabels(),GetReader());
1674         
1675         if  ( !GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCConversion)  )
1676         {
1677           if       ( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0)      ||
1678                     GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCEta)       ) fhTrackMatchedMCParticle[cut]->Fill(cluster->E(), 2.5 );
1679           else if  ( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPhoton)    ) fhTrackMatchedMCParticle[cut]->Fill(cluster->E(), 0.5 );
1680           else if  ( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCElectron)  ) fhTrackMatchedMCParticle[cut]->Fill(cluster->E(), 1.5 );
1681           else                                                                                 fhTrackMatchedMCParticle[cut]->Fill(cluster->E(), 3.5 );
1682           
1683           // Check if several particles contributed to cluster and discard overlapped mesons
1684           if(!GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0) ||
1685              !GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCEta))
1686           {
1687             if(cluster->GetNLabels()==1)
1688             {
1689               fhTrackMatchedDEtaMCNoOverlap[cut]->Fill(cluster->E(),dZ);
1690               fhTrackMatchedDPhiMCNoOverlap[cut]->Fill(cluster->E(),dR);
1691             }
1692             else
1693             {
1694               fhTrackMatchedDEtaMCOverlap[cut]->Fill(cluster->E(),dZ);
1695               fhTrackMatchedDPhiMCOverlap[cut]->Fill(cluster->E(),dR);
1696             }
1697             
1698           }// Check overlaps
1699           
1700         }
1701         else
1702         {
1703           if       ( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0)      ||
1704                     GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCEta)       ) fhTrackMatchedMCParticle[cut]->Fill(cluster->E(), 6.5 );
1705           else if  ( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPhoton)    ) fhTrackMatchedMCParticle[cut]->Fill(cluster->E(), 4.5 );
1706           else if  ( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCElectron)  ) fhTrackMatchedMCParticle[cut]->Fill(cluster->E(), 5.5 );
1707           else                                                                                 fhTrackMatchedMCParticle[cut]->Fill(cluster->E(), 7.5 );
1708           
1709           // Check if several particles contributed to cluster
1710           if(!GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0) ||
1711              !GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCEta))
1712           {
1713             fhTrackMatchedDEtaMCConversion[cut]->Fill(cluster->E(),dZ);
1714             fhTrackMatchedDPhiMCConversion[cut]->Fill(cluster->E(),dR);
1715             
1716           }// Check overlaps
1717           
1718         }
1719         
1720       } // MC
1721       
1722     } // residuals window
1723     
1724   } // Small residual
1725   
1726 }
1727
1728 //___________________________________________
1729 TObjString *  AliAnaPhoton::GetAnalysisCuts()
1730 {
1731   //Save parameters used for analysis
1732   TString parList ; //this will be list of parameters used for this analysis.
1733   const Int_t buffersize = 255;
1734   char onePar[buffersize] ;
1735   
1736   snprintf(onePar,buffersize,"--- AliAnaPhoton ---\n") ;
1737   parList+=onePar ;
1738   snprintf(onePar,buffersize,"Calorimeter: %s\n",fCalorimeter.Data()) ;
1739   parList+=onePar ;
1740   snprintf(onePar,buffersize,"fMinDist =%2.2f (Minimal distance to bad channel to accept cluster) \n",fMinDist) ;
1741   parList+=onePar ;
1742   snprintf(onePar,buffersize,"fMinDist2=%2.2f (Cuts on Minimal distance to study acceptance evaluation) \n",fMinDist2) ;
1743   parList+=onePar ;
1744   snprintf(onePar,buffersize,"fMinDist3=%2.2f (One more cut on distance used for acceptance-efficiency study) \n",fMinDist3) ;
1745   parList+=onePar ;
1746   snprintf(onePar,buffersize,"fRejectTrackMatch: %d\n",fRejectTrackMatch) ;
1747   parList+=onePar ;
1748   
1749   //Get parameters set in base class.
1750   parList += GetBaseParametersList() ;
1751   
1752   //Get parameters set in PID class.
1753   parList += GetCaloPID()->GetPIDParametersList() ;
1754   
1755   //Get parameters set in FiducialCut class (not available yet)
1756   //parlist += GetFidCut()->GetFidCutParametersList()
1757   
1758   return new TObjString(parList) ;
1759 }
1760
1761 //________________________________________________________________________
1762 TList *  AliAnaPhoton::GetCreateOutputObjects()
1763 {
1764   // Create histograms to be saved in output file and
1765   // store them in outputContainer
1766   TList * outputContainer = new TList() ;
1767   outputContainer->SetName("PhotonHistos") ;
1768         
1769   Int_t nptbins  = GetHistogramRanges()->GetHistoPtBins();  Float_t ptmax  = GetHistogramRanges()->GetHistoPtMax();  Float_t ptmin  = GetHistogramRanges()->GetHistoPtMin();
1770   Int_t nphibins = GetHistogramRanges()->GetHistoPhiBins(); Float_t phimax = GetHistogramRanges()->GetHistoPhiMax(); Float_t phimin = GetHistogramRanges()->GetHistoPhiMin();
1771   Int_t netabins = GetHistogramRanges()->GetHistoEtaBins(); Float_t etamax = GetHistogramRanges()->GetHistoEtaMax(); Float_t etamin = GetHistogramRanges()->GetHistoEtaMin();
1772   Int_t ssbins   = GetHistogramRanges()->GetHistoShowerShapeBins();  Float_t ssmax   = GetHistogramRanges()->GetHistoShowerShapeMax();  Float_t ssmin   = GetHistogramRanges()->GetHistoShowerShapeMin();
1773   Int_t nbins    = GetHistogramRanges()->GetHistoNClusterCellBins(); Int_t   nmax    = GetHistogramRanges()->GetHistoNClusterCellMax(); Int_t   nmin    = GetHistogramRanges()->GetHistoNClusterCellMin();
1774   Int_t ntimebins= GetHistogramRanges()->GetHistoTimeBins();         Float_t timemax = GetHistogramRanges()->GetHistoTimeMax();         Float_t timemin = GetHistogramRanges()->GetHistoTimeMin();
1775   
1776   Int_t   nresetabins = GetHistogramRanges()->GetHistoTrackResidualEtaBins();
1777   Float_t resetamax   = GetHistogramRanges()->GetHistoTrackResidualEtaMax();
1778   Float_t resetamin   = GetHistogramRanges()->GetHistoTrackResidualEtaMin();
1779   Int_t   nresphibins = GetHistogramRanges()->GetHistoTrackResidualPhiBins();
1780   Float_t resphimax   = GetHistogramRanges()->GetHistoTrackResidualPhiMax();
1781   Float_t resphimin   = GetHistogramRanges()->GetHistoTrackResidualPhiMin();
1782   
1783   Int_t   ndedxbins   = GetHistogramRanges()->GetHistodEdxBins();
1784   Float_t dedxmax     = GetHistogramRanges()->GetHistodEdxMax();
1785   Float_t dedxmin     = GetHistogramRanges()->GetHistodEdxMin();
1786   Int_t   nPoverEbins = GetHistogramRanges()->GetHistoPOverEBins();
1787   Float_t pOverEmax   = GetHistogramRanges()->GetHistoPOverEMax();
1788   Float_t pOverEmin   = GetHistogramRanges()->GetHistoPOverEMin();
1789   
1790   Int_t bin[] = {0,2,4,6,10,15,20,100}; // energy bins for SS studies
1791   
1792   Int_t nTrigBC  = 1;
1793   Int_t iBCShift = 0;
1794   if(!GetReader()->AreBadTriggerEventsRemoved())
1795   {
1796     nTrigBC = 11;
1797     iBCShift = 5;
1798   }
1799   
1800   TString cut[] = {"Open","Reader","E","Time","NCells","NLM","Fidutial","Matching","Bad","PID"};
1801   for (Int_t i = 0; i < 10 ;  i++)
1802   {
1803     fhClusterCuts[i] = new TH1F(Form("hCut_%d_%s", i, cut[i].Data()),
1804                                 Form("Number of clusters that pass cuts <= %d, %s", i, cut[i].Data()),
1805                                 nptbins,ptmin,ptmax);
1806     fhClusterCuts[i]->SetYTitle("dN/dE ");
1807     fhClusterCuts[i]->SetXTitle("E (GeV)");
1808     outputContainer->Add(fhClusterCuts[i]) ;
1809   }
1810   
1811   fhEClusterSM = new TH2F("hEClusterSM","Raw clusters E and super-module number",
1812                           nptbins,ptmin,ptmax,
1813                           GetCaloUtils()->GetNumberOfSuperModulesUsed(),0,GetCaloUtils()->GetNumberOfSuperModulesUsed());
1814   fhEClusterSM->SetYTitle("SuperModule ");
1815   fhEClusterSM->SetXTitle("E (GeV)");
1816   outputContainer->Add(fhEClusterSM) ;
1817
1818   fhPtClusterSM = new TH2F("hPtClusterSM","Raw clusters p_{T} and super-module number",
1819                           nptbins,ptmin,ptmax,
1820                           GetCaloUtils()->GetNumberOfSuperModulesUsed(),0,GetCaloUtils()->GetNumberOfSuperModulesUsed());
1821   fhPtClusterSM->SetYTitle("SuperModule ");
1822   fhPtClusterSM->SetXTitle("E (GeV)");
1823   outputContainer->Add(fhPtClusterSM) ;
1824   
1825   fhEPhotonSM = new TH2F("hEPhotonSM","Selected clusters E and super-module number",
1826                           nptbins,ptmin,ptmax,
1827                           GetCaloUtils()->GetNumberOfSuperModulesUsed(),0,GetCaloUtils()->GetNumberOfSuperModulesUsed());
1828   fhEPhotonSM->SetYTitle("SuperModule ");
1829   fhEPhotonSM->SetXTitle("E (GeV)");
1830   outputContainer->Add(fhEPhotonSM) ;
1831   
1832   fhPtPhotonSM = new TH2F("hPtPhotonSM","Selected clusters p_{T} and super-module number",
1833                            nptbins,ptmin,ptmax,
1834                            GetCaloUtils()->GetNumberOfSuperModulesUsed(),0,GetCaloUtils()->GetNumberOfSuperModulesUsed());
1835   fhPtPhotonSM->SetYTitle("SuperModule ");
1836   fhPtPhotonSM->SetXTitle("E (GeV)");
1837   outputContainer->Add(fhPtPhotonSM) ;
1838   
1839   fhNCellsE  = new TH2F ("hNCellsE","# of cells in cluster vs E of clusters", nptbins,ptmin,ptmax, nbins,nmin,nmax);
1840   fhNCellsE->SetXTitle("E (GeV)");
1841   fhNCellsE->SetYTitle("# of cells in cluster");
1842   outputContainer->Add(fhNCellsE);
1843   
1844   fhCellsE  = new TH2F ("hCellsE","energy of cells in cluster vs E of clusters", nptbins,ptmin,ptmax, nptbins*2,ptmin,ptmax);
1845   fhCellsE->SetXTitle("E_{cluster} (GeV)");
1846   fhCellsE->SetYTitle("E_{cell} (GeV)");
1847   outputContainer->Add(fhCellsE);
1848   
1849   fhTimePt  = new TH2F ("hTimePt","time of cluster vs pT of clusters", nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
1850   fhTimePt->SetXTitle("p_{T} (GeV/c)");
1851   fhTimePt->SetYTitle("time (ns)");
1852   outputContainer->Add(fhTimePt);
1853   
1854   fhMaxCellDiffClusterE  = new TH2F ("hMaxCellDiffClusterE","energy vs difference of cluster energy - max cell energy / cluster energy, good clusters",
1855                                      nptbins,ptmin,ptmax, 500,0,1.);
1856   fhMaxCellDiffClusterE->SetXTitle("E_{cluster} (GeV) ");
1857   fhMaxCellDiffClusterE->SetYTitle("(E_{cluster} - E_{cell max})/ E_{cluster}");
1858   outputContainer->Add(fhMaxCellDiffClusterE);
1859   
1860   fhEPhoton  = new TH1F("hEPhoton","Number of #gamma over calorimeter vs energy",nptbins,ptmin,ptmax);
1861   fhEPhoton->SetYTitle("N");
1862   fhEPhoton->SetXTitle("E_{#gamma}(GeV)");
1863   outputContainer->Add(fhEPhoton) ;
1864   
1865   fhPtPhoton  = new TH1F("hPtPhoton","Number of #gamma over calorimeter vs p_{T}",nptbins,ptmin,ptmax);
1866   fhPtPhoton->SetYTitle("N");
1867   fhPtPhoton->SetXTitle("p_{T #gamma}(GeV/c)");
1868   outputContainer->Add(fhPtPhoton) ;
1869   
1870   fhPtCentralityPhoton  = new TH2F("hPtCentralityPhoton","centrality vs p_{T}",nptbins,ptmin,ptmax, 100,0,100);
1871   fhPtCentralityPhoton->SetYTitle("Centrality");
1872   fhPtCentralityPhoton->SetXTitle("p_{T}(GeV/c)");
1873   outputContainer->Add(fhPtCentralityPhoton) ;
1874   
1875   fhPtEventPlanePhoton  = new TH2F("hPtEventPlanePhoton","centrality vs p_{T}",nptbins,ptmin,ptmax, 100,0,TMath::Pi());
1876   fhPtEventPlanePhoton->SetYTitle("Event plane angle (rad)");
1877   fhPtEventPlanePhoton->SetXTitle("p_{T} (GeV/c)");
1878   outputContainer->Add(fhPtEventPlanePhoton) ;
1879   
1880   fhEtaPhi  = new TH2F
1881   ("hEtaPhi","cluster,E > 0.5 GeV, #eta vs #phi",netabins,etamin,etamax,nphibins,phimin,phimax);
1882   fhEtaPhi->SetYTitle("#phi (rad)");
1883   fhEtaPhi->SetXTitle("#eta");
1884   outputContainer->Add(fhEtaPhi) ;
1885   
1886   if(fCalorimeter=="EMCAL" && fFillEMCALBCHistograms)
1887   {
1888     fhEtaPhiEMCALBC0  = new TH2F
1889     ("hEtaPhiEMCALBC0","cluster,E > 2 GeV, #eta vs #phi, for clusters with |time| < 25 ns, EMCAL-BC=0",netabins,etamin,etamax,nphibins,phimin,phimax);
1890     fhEtaPhiEMCALBC0->SetYTitle("#phi (rad)");
1891     fhEtaPhiEMCALBC0->SetXTitle("#eta");
1892     outputContainer->Add(fhEtaPhiEMCALBC0) ;
1893     
1894     fhEtaPhiEMCALBC1  = new TH2F
1895     ("hEtaPhiEMCALBC1","cluster,E > 2 GeV, #eta vs #phi, for clusters with 25 < |time| < 75 ns, EMCAL-BC=1",netabins,etamin,etamax,nphibins,phimin,phimax);
1896     fhEtaPhiEMCALBC1->SetYTitle("#phi (rad)");
1897     fhEtaPhiEMCALBC1->SetXTitle("#eta");
1898     outputContainer->Add(fhEtaPhiEMCALBC1) ;
1899     
1900     fhEtaPhiEMCALBCN  = new TH2F
1901     ("hEtaPhiEMCALBCN","cluster,E > 2 GeV, #eta vs #phi, for clusters with |time| > 75 ns, EMCAL-BC>1",netabins,etamin,etamax,nphibins,phimin,phimax);
1902     fhEtaPhiEMCALBCN->SetYTitle("#phi (rad)");
1903     fhEtaPhiEMCALBCN->SetXTitle("#eta");
1904     outputContainer->Add(fhEtaPhiEMCALBCN) ;
1905     
1906     for(Int_t i = 0; i < nTrigBC; i++)
1907     {
1908       fhEtaPhiTriggerEMCALBC[i] = new TH2F
1909       (Form("hEtaPhiTriggerEMCALBC%d",i-iBCShift),
1910        Form("cluster E > 2 GeV, #eta vs #phi, Trigger EMCAL-BC=%d",i-iBCShift),
1911        netabins,etamin,etamax,nphibins,phimin,phimax);
1912       fhEtaPhiTriggerEMCALBC[i]->SetYTitle("#phi (rad)");
1913       fhEtaPhiTriggerEMCALBC[i]->SetXTitle("#eta");
1914       outputContainer->Add(fhEtaPhiTriggerEMCALBC[i]) ;
1915       
1916       fhTimeTriggerEMCALBC[i] = new TH2F
1917       (Form("hTimeTriggerEMCALBC%d",i-iBCShift),
1918        Form("cluster time vs E of clusters, Trigger EMCAL-BC=%d",i-iBCShift),
1919        nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
1920       fhTimeTriggerEMCALBC[i]->SetXTitle("E (GeV)");
1921       fhTimeTriggerEMCALBC[i]->SetYTitle("time (ns)");
1922       outputContainer->Add(fhTimeTriggerEMCALBC[i]);
1923       
1924       fhTimeTriggerEMCALBCPileUpSPD[i] = new TH2F
1925       (Form("hTimeTriggerEMCALBC%dPileUpSPD",i-iBCShift),
1926        Form("cluster time vs E of clusters, Trigger EMCAL-BC=%d",i-iBCShift),
1927        nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
1928       fhTimeTriggerEMCALBCPileUpSPD[i]->SetXTitle("E (GeV)");
1929       fhTimeTriggerEMCALBCPileUpSPD[i]->SetYTitle("time (ns)");
1930       outputContainer->Add(fhTimeTriggerEMCALBCPileUpSPD[i]);
1931       
1932       fhEtaPhiTriggerEMCALBCUM[i] = new TH2F
1933       (Form("hEtaPhiTriggerEMCALBC%d_UnMatch",i-iBCShift),
1934        Form("cluster E > 2 GeV, #eta vs #phi, unmatched trigger EMCAL-BC=%d",i-iBCShift),
1935        netabins,etamin,etamax,nphibins,phimin,phimax);
1936       fhEtaPhiTriggerEMCALBCUM[i]->SetYTitle("#phi (rad)");
1937       fhEtaPhiTriggerEMCALBCUM[i]->SetXTitle("#eta");
1938       outputContainer->Add(fhEtaPhiTriggerEMCALBCUM[i]) ;
1939       
1940       fhTimeTriggerEMCALBCUM[i] = new TH2F
1941       (Form("hTimeTriggerEMCALBC%d_UnMatch",i-iBCShift),
1942        Form("cluster time vs E of clusters, unmatched trigger EMCAL-BC=%d",i-iBCShift),
1943        nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
1944       fhTimeTriggerEMCALBCUM[i]->SetXTitle("E (GeV)");
1945       fhTimeTriggerEMCALBCUM[i]->SetYTitle("time (ns)");
1946       outputContainer->Add(fhTimeTriggerEMCALBCUM[i]);
1947       
1948       fhEtaPhiTriggerEMCALBCCluster[i] = new TH2F
1949       (Form("hEtaPhiTriggerEMCALBC%d_OnlyTrigger",i-iBCShift),
1950        Form("trigger cluster, #eta vs #phi, Trigger EMCAL-BC=%d",i-iBCShift),
1951        netabins,etamin,etamax,nphibins,phimin,phimax);
1952       fhEtaPhiTriggerEMCALBCCluster[i]->SetYTitle("#phi (rad)");
1953       fhEtaPhiTriggerEMCALBCCluster[i]->SetXTitle("#eta");
1954       outputContainer->Add(fhEtaPhiTriggerEMCALBCCluster[i]) ;
1955             
1956       fhEtaPhiTriggerEMCALBCUMCluster[i] = new TH2F
1957       (Form("hEtaPhiTriggerEMCALBC%d_OnlyTrigger_UnMatch",i-iBCShift),
1958        Form("trigger cluster, #eta vs #phi, unmatched trigger EMCAL-BC=%d",i-iBCShift),
1959        netabins,etamin,etamax,nphibins,phimin,phimax);
1960       fhEtaPhiTriggerEMCALBCUMCluster[i]->SetYTitle("#phi (rad)");
1961       fhEtaPhiTriggerEMCALBCUMCluster[i]->SetXTitle("#eta");
1962       outputContainer->Add(fhEtaPhiTriggerEMCALBCUMCluster[i]) ;
1963     }
1964     
1965     fhTimeTriggerEMCALBCCluster = new TH2F("hTimeTriggerEMCALBC_OnlyTrigger",
1966                                            "trigger cluster time vs E of clusters",
1967                                            nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
1968     fhTimeTriggerEMCALBCCluster->SetXTitle("E (GeV)");
1969     fhTimeTriggerEMCALBCCluster->SetYTitle("time (ns)");
1970     outputContainer->Add(fhTimeTriggerEMCALBCCluster);
1971     
1972     fhTimeTriggerEMCALBCUMCluster = new TH2F("hTimeTriggerEMCALBC_OnlyTrigger_UnMatch",
1973                                              "trigger cluster time vs E of clusters, unmatched trigger",
1974                                              nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
1975     fhTimeTriggerEMCALBCUMCluster->SetXTitle("E (GeV)");
1976     fhTimeTriggerEMCALBCUMCluster->SetYTitle("time (ns)");
1977     outputContainer->Add(fhTimeTriggerEMCALBCUMCluster);
1978     
1979     fhEtaPhiTriggerEMCALBCClusterOverTh = new TH2F
1980     ("hEtaPhiTriggerEMCALBC0_OnlyTrigger_OverThreshold",
1981      "trigger cluster E > trigger threshold, #eta vs #phi, Trigger EMCAL-BC=0",
1982      netabins,etamin,etamax,nphibins,phimin,phimax);
1983     fhEtaPhiTriggerEMCALBCClusterOverTh->SetYTitle("#phi (rad)");
1984     fhEtaPhiTriggerEMCALBCClusterOverTh->SetXTitle("#eta");
1985     outputContainer->Add(fhEtaPhiTriggerEMCALBCClusterOverTh) ;
1986     
1987     fhEtaPhiTriggerEMCALBCUMClusterOverTh = new TH2F
1988     ("hEtaPhiTriggerEMCALBC0_OnlyTrigger_OverThreshold_UnMatch",
1989      "trigger cluster E > trigger threshold, #eta vs #phi, unmatched trigger EMCAL-BC=0",
1990      netabins,etamin,etamax,nphibins,phimin,phimax);
1991     fhEtaPhiTriggerEMCALBCUMClusterOverTh->SetYTitle("#phi (rad)");
1992     fhEtaPhiTriggerEMCALBCUMClusterOverTh->SetXTitle("#eta");
1993     outputContainer->Add(fhEtaPhiTriggerEMCALBCUMClusterOverTh) ;
1994     
1995     fhEtaPhiTriggerEMCALBCClusterBelowTh1 = new TH2F
1996     ("hEtaPhiTriggerEMCALBC0_OnlyTrigger_BelowThreshold1",
1997      "trigger cluster thresh-1 < E < thres, #eta vs #phi, Trigger EMCAL-BC=0",
1998      netabins,etamin,etamax,nphibins,phimin,phimax);
1999     fhEtaPhiTriggerEMCALBCClusterBelowTh1->SetYTitle("#phi (rad)");
2000     fhEtaPhiTriggerEMCALBCClusterBelowTh1->SetXTitle("#eta");
2001     outputContainer->Add(fhEtaPhiTriggerEMCALBCClusterBelowTh1) ;
2002     
2003     fhEtaPhiTriggerEMCALBCUMClusterBelowTh1 = new TH2F
2004     ("hEtaPhiTriggerEMCALBC0_OnlyTrigger_BelowThreshold1_UnMatch",
2005      "trigger cluster thresh-1 < E < thres, #eta vs #phi, unmatched trigger EMCAL-BC=0",
2006      netabins,etamin,etamax,nphibins,phimin,phimax);
2007     fhEtaPhiTriggerEMCALBCUMClusterBelowTh1->SetYTitle("#phi (rad)");
2008     fhEtaPhiTriggerEMCALBCUMClusterBelowTh1->SetXTitle("#eta");
2009     outputContainer->Add(fhEtaPhiTriggerEMCALBCUMClusterBelowTh1) ;
2010     
2011     fhEtaPhiTriggerEMCALBCClusterBelowTh2 = new TH2F
2012     ("hEtaPhiTriggerEMCALBC0_OnlyTrigger_BelowThreshold2",
2013      "trigger cluster thresh-2 < E < thres, #eta vs #phi, Trigger EMCAL-BC=0",
2014      netabins,etamin,etamax,nphibins,phimin,phimax);
2015     fhEtaPhiTriggerEMCALBCClusterBelowTh2->SetYTitle("#phi (rad)");
2016     fhEtaPhiTriggerEMCALBCClusterBelowTh2->SetXTitle("#eta");
2017     outputContainer->Add(fhEtaPhiTriggerEMCALBCClusterBelowTh2) ;
2018     
2019     fhEtaPhiTriggerEMCALBCUMClusterBelowTh2 = new TH2F
2020     ("hEtaPhiTriggerEMCALBC0_OnlyTrigger_BelowThreshold2_UnMatch",
2021      "trigger cluster thresh-2 < E < thres, #eta vs #phi, unmatched trigger EMCAL-BC=0",
2022      netabins,etamin,etamax,nphibins,phimin,phimax);
2023     fhEtaPhiTriggerEMCALBCUMClusterBelowTh2->SetYTitle("#phi (rad)");
2024     fhEtaPhiTriggerEMCALBCUMClusterBelowTh2->SetXTitle("#eta");
2025     outputContainer->Add(fhEtaPhiTriggerEMCALBCUMClusterBelowTh2) ;
2026     
2027     if(!GetReader()->AreBadTriggerEventsRemoved())
2028     {
2029       fhEtaPhiTriggerEMCALBCExotic = new TH2F
2030       ("hEtaPhiTriggerExotic",
2031        "cluster E > 2 GeV, #eta vs #phi, Trigger Exotic",
2032        netabins,etamin,etamax,nphibins,phimin,phimax);
2033       fhEtaPhiTriggerEMCALBCExotic->SetYTitle("#phi (rad)");
2034       fhEtaPhiTriggerEMCALBCExotic->SetXTitle("#eta");
2035       outputContainer->Add(fhEtaPhiTriggerEMCALBCExotic) ;
2036       
2037       fhTimeTriggerEMCALBCExotic = new TH2F
2038       ("hTimeTriggerExotic",
2039        "cluster time vs E of clusters, Trigger Exotic ",
2040        nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
2041       fhTimeTriggerEMCALBCExotic->SetXTitle("E (GeV)");
2042       fhTimeTriggerEMCALBCExotic->SetYTitle("time (ns)");
2043       outputContainer->Add(fhTimeTriggerEMCALBCExotic);
2044       
2045       fhEtaPhiTriggerEMCALBCUMExotic = new TH2F
2046       ("hEtaPhiTriggerExotic_UnMatch",
2047        "cluster E > 2 GeV, #eta vs #phi, unmatched trigger Exotic",
2048        netabins,etamin,etamax,nphibins,phimin,phimax);
2049       fhEtaPhiTriggerEMCALBCUMExotic->SetYTitle("#phi (rad)");
2050       fhEtaPhiTriggerEMCALBCUMExotic->SetXTitle("#eta");
2051       outputContainer->Add(fhEtaPhiTriggerEMCALBCUMExotic) ;
2052       
2053       fhTimeTriggerEMCALBCUMExotic = new TH2F
2054       ("hTimeTriggerExotic_UnMatch",
2055        "cluster time vs E of clusters, unmatched trigger Exotic",
2056        nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
2057       fhTimeTriggerEMCALBCUMExotic->SetXTitle("E (GeV)");
2058       fhTimeTriggerEMCALBCUMExotic->SetYTitle("time (ns)");
2059       outputContainer->Add(fhTimeTriggerEMCALBCUMExotic);
2060       
2061       fhEtaPhiTriggerEMCALBCExoticCluster = new TH2F
2062       ("hEtaPhiTriggerExotic_OnlyTrigger",
2063        "trigger cluster E > 2 GeV, #eta vs #phi, Trigger Exotic",
2064        netabins,etamin,etamax,nphibins,phimin,phimax);
2065       fhEtaPhiTriggerEMCALBCExoticCluster->SetYTitle("#phi (rad)");
2066       fhEtaPhiTriggerEMCALBCExoticCluster->SetXTitle("#eta");
2067       outputContainer->Add(fhEtaPhiTriggerEMCALBCExoticCluster) ;
2068       
2069       fhTimeTriggerEMCALBCExoticCluster = new TH2F
2070       ("hTimeTriggerExotic_OnlyTrigger",
2071        "trigger cluster time vs E of clusters, Trigger Exotic",
2072        nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
2073       fhTimeTriggerEMCALBCExoticCluster->SetXTitle("E (GeV)");
2074       fhTimeTriggerEMCALBCExoticCluster->SetYTitle("time (ns)");
2075       outputContainer->Add(fhTimeTriggerEMCALBCExoticCluster);
2076       
2077       fhEtaPhiTriggerEMCALBCUMExoticCluster = new TH2F
2078       ("hEtaPhiTriggerExotic_OnlyTrigger_UnMatch",
2079        "trigger cluster E > 2 GeV, #eta vs #phi, unmatched trigger Exotic",
2080        netabins,etamin,etamax,nphibins,phimin,phimax);
2081       fhEtaPhiTriggerEMCALBCUMExoticCluster->SetYTitle("#phi (rad)");
2082       fhEtaPhiTriggerEMCALBCUMExoticCluster->SetXTitle("#eta");
2083       outputContainer->Add(fhEtaPhiTriggerEMCALBCUMExoticCluster) ;
2084       
2085       fhTimeTriggerEMCALBCUMExoticCluster = new TH2F
2086       ("hTimeTriggerExotic_OnlyTrigger_UnMatch",
2087        "trigger cluster time vs E of clusters, unmatched trigger Exotic",
2088        nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
2089       fhTimeTriggerEMCALBCUMExoticCluster->SetXTitle("E (GeV)");
2090       fhTimeTriggerEMCALBCUMExoticCluster->SetYTitle("time (ns)");
2091       outputContainer->Add(fhTimeTriggerEMCALBCUMExoticCluster);
2092       
2093       fhEtaPhiTriggerEMCALBCBad = new TH2F
2094       ("hEtaPhiTriggerBad",
2095        "cluster E > 2 GeV, #eta vs #phi, Trigger Bad",
2096        netabins,etamin,etamax,nphibins,phimin,phimax);
2097       fhEtaPhiTriggerEMCALBCBad->SetYTitle("#phi (rad)");
2098       fhEtaPhiTriggerEMCALBCBad->SetXTitle("#eta");
2099       outputContainer->Add(fhEtaPhiTriggerEMCALBCBad) ;
2100       
2101       fhTimeTriggerEMCALBCBad = new TH2F
2102       ("hTimeTriggerBad",
2103        "cluster time vs E of clusters, Trigger Bad ",
2104        nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
2105       fhTimeTriggerEMCALBCBad->SetXTitle("E (GeV)");
2106       fhTimeTriggerEMCALBCBad->SetYTitle("time (ns)");
2107       outputContainer->Add(fhTimeTriggerEMCALBCBad);
2108       
2109       fhEtaPhiTriggerEMCALBCUMBad = new TH2F
2110       ("hEtaPhiTriggerBad_UnMatch",
2111        "cluster E > 2 GeV, #eta vs #phi, unmatched trigger Bad",
2112        netabins,etamin,etamax,nphibins,phimin,phimax);
2113       fhEtaPhiTriggerEMCALBCUMBad->SetYTitle("#phi (rad)");
2114       fhEtaPhiTriggerEMCALBCUMBad->SetXTitle("#eta");
2115       outputContainer->Add(fhEtaPhiTriggerEMCALBCUMBad) ;
2116       
2117       fhTimeTriggerEMCALBCUMBad = new TH2F
2118       ("hTimeTriggerBad_UnMatch",
2119        "cluster time vs E of clusters, unmatched trigger Bad",
2120        nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
2121       fhTimeTriggerEMCALBCUMBad->SetXTitle("E (GeV)");
2122       fhTimeTriggerEMCALBCUMBad->SetYTitle("time (ns)");
2123       outputContainer->Add(fhTimeTriggerEMCALBCUMBad);
2124       
2125       fhEtaPhiTriggerEMCALBCBadCluster = new TH2F
2126       ("hEtaPhiTriggerBad_OnlyTrigger",
2127        "trigger cluster E > 2 GeV, #eta vs #phi, Trigger Bad",
2128        netabins,etamin,etamax,nphibins,phimin,phimax);
2129       fhEtaPhiTriggerEMCALBCBadCluster->SetYTitle("#phi (rad)");
2130       fhEtaPhiTriggerEMCALBCBadCluster->SetXTitle("#eta");
2131       outputContainer->Add(fhEtaPhiTriggerEMCALBCBadCluster) ;
2132       
2133       fhTimeTriggerEMCALBCBadCluster = new TH2F
2134       ("hTimeTriggerBad_OnlyTrigger",
2135        "trigger cluster time vs E of clusters, Trigger Bad",
2136        nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
2137       fhTimeTriggerEMCALBCBadCluster->SetXTitle("E (GeV)");
2138       fhTimeTriggerEMCALBCBadCluster->SetYTitle("time (ns)");
2139       outputContainer->Add(fhTimeTriggerEMCALBCBadCluster);
2140       
2141       fhEtaPhiTriggerEMCALBCUMBadCluster = new TH2F
2142       ("hEtaPhiTriggerBad_OnlyTrigger_UnMatch",
2143        "trigger cluster E > 2 GeV, #eta vs #phi, unmatched trigger Bad",
2144        netabins,etamin,etamax,nphibins,phimin,phimax);
2145       fhEtaPhiTriggerEMCALBCUMBadCluster->SetYTitle("#phi (rad)");
2146       fhEtaPhiTriggerEMCALBCUMBadCluster->SetXTitle("#eta");
2147       outputContainer->Add(fhEtaPhiTriggerEMCALBCUMBadCluster) ;
2148       
2149       fhTimeTriggerEMCALBCUMBadCluster = new TH2F
2150       ("hTimeTriggerBad_OnlyTrigger_UnMatch",
2151        "trigger cluster time vs E of clusters, unmatched trigger Bad",
2152        nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
2153       fhTimeTriggerEMCALBCUMBadCluster->SetXTitle("E (GeV)");
2154       fhTimeTriggerEMCALBCUMBadCluster->SetYTitle("time (ns)");
2155       outputContainer->Add(fhTimeTriggerEMCALBCUMBadCluster);
2156       
2157       fhEtaPhiTriggerEMCALBCBadExotic = new TH2F
2158       ("hEtaPhiTriggerBadExotic",
2159        "cluster E > 2 GeV, #eta vs #phi, Trigger Bad&Exotic",
2160        netabins,etamin,etamax,nphibins,phimin,phimax);
2161       fhEtaPhiTriggerEMCALBCBadExotic->SetYTitle("#phi (rad)");
2162       fhEtaPhiTriggerEMCALBCBadExotic->SetXTitle("#eta");
2163       outputContainer->Add(fhEtaPhiTriggerEMCALBCBadExotic) ;
2164       
2165       fhTimeTriggerEMCALBCBadExotic = new TH2F
2166       ("hTimeTriggerBadExotic",
2167        "cluster time vs E of clusters, Trigger Bad&Exotic ",
2168        nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
2169       fhTimeTriggerEMCALBCBadExotic->SetXTitle("E (GeV)");
2170       fhTimeTriggerEMCALBCBadExotic->SetYTitle("time (ns)");
2171       outputContainer->Add(fhTimeTriggerEMCALBCBadExotic);
2172       
2173       fhEtaPhiTriggerEMCALBCUMBadExotic = new TH2F
2174       ("hEtaPhiTriggerBadExotic_UnMatch",
2175        "cluster E > 2 GeV, #eta vs #phi, unmatched trigger Bad&Exotic",
2176        netabins,etamin,etamax,nphibins,phimin,phimax);
2177       fhEtaPhiTriggerEMCALBCUMBadExotic->SetYTitle("#phi (rad)");
2178       fhEtaPhiTriggerEMCALBCUMBadExotic->SetXTitle("#eta");
2179       outputContainer->Add(fhEtaPhiTriggerEMCALBCUMBadExotic) ;
2180       
2181       fhTimeTriggerEMCALBCUMBadExotic = new TH2F
2182       ("hTimeTriggerBadExotic_UnMatch",
2183        "cluster time vs E of clusters, unmatched trigger Bad&Exotic",
2184        nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
2185       fhTimeTriggerEMCALBCUMBadExotic->SetXTitle("E (GeV)");
2186       fhTimeTriggerEMCALBCUMBadExotic->SetYTitle("time (ns)");
2187       outputContainer->Add(fhTimeTriggerEMCALBCUMBadExotic);
2188       
2189       fhEtaPhiTriggerEMCALBCBadExoticCluster = new TH2F
2190       ("hEtaPhiTriggerBadExotic_OnlyTrigger",
2191        "trigger cluster E > 2 GeV, #eta vs #phi, Trigger Bad&Exotic",
2192        netabins,etamin,etamax,nphibins,phimin,phimax);
2193       fhEtaPhiTriggerEMCALBCBadExoticCluster->SetYTitle("#phi (rad)");
2194       fhEtaPhiTriggerEMCALBCBadExoticCluster->SetXTitle("#eta");
2195       outputContainer->Add(fhEtaPhiTriggerEMCALBCBadExoticCluster) ;
2196       
2197       fhTimeTriggerEMCALBCBadExoticCluster = new TH2F
2198       ("hTimeTriggerBadExotic_OnlyTrigger",
2199        "trigger cluster time vs E of clusters, Trigger Bad&Exotic",
2200        nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
2201       fhTimeTriggerEMCALBCBadExoticCluster->SetXTitle("E (GeV)");
2202       fhTimeTriggerEMCALBCBadExoticCluster->SetYTitle("time (ns)");
2203       outputContainer->Add(fhTimeTriggerEMCALBCBadExoticCluster);
2204       
2205       fhEtaPhiTriggerEMCALBCUMBadExoticCluster = new TH2F
2206       ("hEtaPhiTriggerBadExotic_OnlyTrigger_UnMatch",
2207        "trigger cluster E > 2 GeV, #eta vs #phi, unmatched trigger Bad&Exotic",
2208        netabins,etamin,etamax,nphibins,phimin,phimax);
2209       fhEtaPhiTriggerEMCALBCUMBadExoticCluster->SetYTitle("#phi (rad)");
2210       fhEtaPhiTriggerEMCALBCUMBadExoticCluster->SetXTitle("#eta");
2211       outputContainer->Add(fhEtaPhiTriggerEMCALBCUMBadExoticCluster) ;
2212       
2213       fhTimeTriggerEMCALBCUMBadExoticCluster = new TH2F
2214       ("hTimeTriggerBadExotic_OnlyTrigger_UnMatch",
2215        "trigger cluster time vs E of clusters, unmatched trigger Bad&Exotic",
2216        nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
2217       fhTimeTriggerEMCALBCUMBadExoticCluster->SetXTitle("E (GeV)");
2218       fhTimeTriggerEMCALBCUMBadExoticCluster->SetYTitle("time (ns)");
2219       outputContainer->Add(fhTimeTriggerEMCALBCUMBadExoticCluster);
2220       
2221       fhTimeTriggerEMCALBCBadMaxCell = new TH2F
2222       ("hTimeTriggerBadMaxCell",
2223        "cluster time vs E of clusters, Trigger BadMaxCell",
2224        nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
2225       fhTimeTriggerEMCALBCBadMaxCell->SetXTitle("E (GeV)");
2226       fhTimeTriggerEMCALBCBadMaxCell->SetYTitle("time (ns)");
2227       outputContainer->Add(fhTimeTriggerEMCALBCBadMaxCell);
2228       
2229       fhTimeTriggerEMCALBCUMBadMaxCell = new TH2F
2230       ("hTimeTriggerBadMaxCell_UnMatch",
2231        "cluster time vs E of clusters, unmatched trigger BadMaxCell",
2232        nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
2233       fhTimeTriggerEMCALBCUMBadMaxCell->SetXTitle("E (GeV)");
2234       fhTimeTriggerEMCALBCUMBadMaxCell->SetYTitle("time (ns)");
2235       outputContainer->Add(fhTimeTriggerEMCALBCUMBadMaxCell);
2236       
2237       
2238       fhTimeTriggerEMCALBCBadMaxCellExotic = new TH2F
2239       ("hTimeTriggerBadMaxCellExotic",
2240        "cluster time vs E of clusters, Trigger BadMaxCell&Exotic",
2241        nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
2242       fhTimeTriggerEMCALBCBadMaxCellExotic->SetXTitle("E (GeV)");
2243       fhTimeTriggerEMCALBCBadMaxCellExotic->SetYTitle("time (ns)");
2244       outputContainer->Add(fhTimeTriggerEMCALBCBadMaxCellExotic);
2245       
2246       fhTimeTriggerEMCALBCUMBadMaxCellExotic = new TH2F
2247       ("hTimeTriggerBadMaxCellExotic_UnMatch",
2248        "cluster time vs E of clusters, unmatched trigger BadMaxCell&Exotic",
2249        nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
2250       fhTimeTriggerEMCALBCUMBadMaxCellExotic->SetXTitle("E (GeV)");
2251       fhTimeTriggerEMCALBCUMBadMaxCellExotic->SetYTitle("time (ns)");
2252       outputContainer->Add(fhTimeTriggerEMCALBCUMBadMaxCellExotic);
2253       
2254       fhTimeNoTrigger = new TH2F
2255       ("hTimeNoTrigger",
2256        "events with no foundable trigger, time vs e of clusters",
2257        nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
2258       fhTimeNoTrigger->SetXTitle("E (GeV)");
2259       fhTimeNoTrigger->SetYTitle("time (ns)");
2260       outputContainer->Add(fhTimeNoTrigger);
2261       
2262       fhEtaPhiNoTrigger = new TH2F
2263       ("hEtaPhiNoTrigger",
2264        "events with no foundable trigger, eta vs phi of clusters",
2265        netabins,etamin,etamax,nphibins,phimin,phimax);
2266       fhEtaPhiNoTrigger->SetYTitle("#phi (rad)");
2267       fhEtaPhiNoTrigger->SetXTitle("#eta");
2268       outputContainer->Add(fhEtaPhiNoTrigger) ;
2269     }
2270     
2271     fhEtaPhiTriggerEMCALBCUMReMatchOpenTimeCluster = new TH2F("hEtaPhiTriggerEMCALBC0_OnlyTrigger_UnMatch_ReMatch_OpenTime",
2272                                                        "cluster E > 2 GeV, #eta vs #phi, Trigger EMCAL-BC=0, un match, rematch open time",
2273                                                        netabins,etamin,etamax,nphibins,phimin,phimax);
2274     fhEtaPhiTriggerEMCALBCUMReMatchOpenTimeCluster->SetYTitle("#phi (rad)");
2275     fhEtaPhiTriggerEMCALBCUMReMatchOpenTimeCluster->SetXTitle("#eta");
2276     outputContainer->Add(fhEtaPhiTriggerEMCALBCUMReMatchOpenTimeCluster) ;
2277     
2278     fhTimeTriggerEMCALBCUMReMatchOpenTimeCluster = new TH2F("hTimeTrigger_OnlyTrigger_UnMatch_ReMatch_OpenTime",
2279                                                      "cluster time vs E of clusters, no match, rematch open time",
2280                                                      nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
2281     fhTimeTriggerEMCALBCUMReMatchOpenTimeCluster->SetXTitle("E (GeV)");
2282     fhTimeTriggerEMCALBCUMReMatchOpenTimeCluster->SetYTitle("time (ns)");
2283     outputContainer->Add(fhTimeTriggerEMCALBCUMReMatchOpenTimeCluster);
2284     
2285     
2286     fhEtaPhiTriggerEMCALBCUMReMatchCheckNeighCluster = new TH2F("hEtaPhiTriggerEMCALBC0_OnlyTrigger_UnMatch_ReMatch_CheckNeighbours",
2287                                                          "cluster E > 2 GeV, #eta vs #phi, Trigger EMCAL-BC=0, un match, rematch with neighbour patches",
2288                                                          netabins,etamin,etamax,nphibins,phimin,phimax);
2289     fhEtaPhiTriggerEMCALBCUMReMatchCheckNeighCluster->SetYTitle("#phi (rad)");
2290     fhEtaPhiTriggerEMCALBCUMReMatchCheckNeighCluster->SetXTitle("#eta");
2291     outputContainer->Add(fhEtaPhiTriggerEMCALBCUMReMatchCheckNeighCluster) ;
2292     
2293     fhTimeTriggerEMCALBCUMReMatchCheckNeighCluster = new TH2F("hTimeTrigger_OnlyTrigger_UnMatch_ReMatch_CheckNeighbours",
2294                                                        "cluster time vs E of clusters, no match, rematch with neigbour parches",
2295                                                        nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
2296     fhTimeTriggerEMCALBCUMReMatchCheckNeighCluster->SetXTitle("E (GeV)");
2297     fhTimeTriggerEMCALBCUMReMatchCheckNeighCluster->SetYTitle("time (ns)");
2298     outputContainer->Add(fhTimeTriggerEMCALBCUMReMatchCheckNeighCluster);
2299     
2300     fhEtaPhiTriggerEMCALBCUMReMatchBothCluster = new TH2F("hEtaPhiTriggerEMCALBC0_OnlyTrigger_UnMatch_ReMatch_Both",
2301                                                    "cluster E > 2 GeV, #eta vs #phi, Trigger EMCAL-BC=0, un match, rematch open time and neighbour",
2302                                                    netabins,etamin,etamax,nphibins,phimin,phimax);
2303     fhEtaPhiTriggerEMCALBCUMReMatchBothCluster->SetYTitle("#phi (rad)");
2304     fhEtaPhiTriggerEMCALBCUMReMatchBothCluster->SetXTitle("#eta");
2305     outputContainer->Add(fhEtaPhiTriggerEMCALBCUMReMatchBothCluster) ;
2306     
2307     fhTimeTriggerEMCALBCUMReMatchBothCluster = new TH2F("hTimeTrigger_OnlyTrigger_UnMatch_ReMatch_Both",
2308                                                  "cluster time vs E of clusters, no match, rematch open time and neigbour",
2309                                                  nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
2310     fhTimeTriggerEMCALBCUMReMatchBothCluster->SetXTitle("E (GeV)");
2311     fhTimeTriggerEMCALBCUMReMatchBothCluster->SetYTitle("time (ns)");
2312     outputContainer->Add(fhTimeTriggerEMCALBCUMReMatchBothCluster);
2313     
2314     fhTimeTriggerEMCALBC0UMReMatchOpenTime = new TH2F("hTimeTriggerBC0_UnMatch_ReMatch_OpenTime",
2315                                                             "cluster time vs E of clusters, no match, rematch open time",
2316                                                             nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
2317     fhTimeTriggerEMCALBC0UMReMatchOpenTime->SetXTitle("E (GeV)");
2318     fhTimeTriggerEMCALBC0UMReMatchOpenTime->SetYTitle("time (ns)");
2319     outputContainer->Add(fhTimeTriggerEMCALBC0UMReMatchOpenTime);
2320     
2321         
2322     fhTimeTriggerEMCALBC0UMReMatchCheckNeigh = new TH2F("hTimeTriggerBC0_UnMatch_ReMatch_CheckNeighbours",
2323                                                               "cluster time vs E of clusters, no match, rematch with neigbour parches",
2324                                                               nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
2325     fhTimeTriggerEMCALBC0UMReMatchCheckNeigh->SetXTitle("E (GeV)");
2326     fhTimeTriggerEMCALBC0UMReMatchCheckNeigh->SetYTitle("time (ns)");
2327     outputContainer->Add(fhTimeTriggerEMCALBC0UMReMatchCheckNeigh);
2328         
2329     fhTimeTriggerEMCALBC0UMReMatchBoth = new TH2F("hTimeTriggerBC0_UnMatch_ReMatch_Both",
2330                                                         "cluster time vs E of clusters, no match, rematch open time and neigbour",
2331                                                         nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
2332     fhTimeTriggerEMCALBC0UMReMatchBoth->SetXTitle("E (GeV)");
2333     fhTimeTriggerEMCALBC0UMReMatchBoth->SetYTitle("time (ns)");
2334     outputContainer->Add(fhTimeTriggerEMCALBC0UMReMatchBoth);
2335     
2336   }
2337   
2338   fhPhiPhoton  = new TH2F
2339   ("hPhiPhoton","#phi_{#gamma} vs p_{T}",nptbins,ptmin,ptmax,nphibins,phimin,phimax);
2340   fhPhiPhoton->SetYTitle("#phi (rad)");
2341   fhPhiPhoton->SetXTitle("p_{T #gamma} (GeV/c)");
2342   outputContainer->Add(fhPhiPhoton) ;
2343   
2344   fhEtaPhoton  = new TH2F
2345   ("hEtaPhoton","#eta_{#gamma} vs p_{T}",nptbins,ptmin,ptmax,netabins,etamin,etamax);
2346   fhEtaPhoton->SetYTitle("#eta");
2347   fhEtaPhoton->SetXTitle("p_{T #gamma} (GeV/c)");
2348   outputContainer->Add(fhEtaPhoton) ;
2349   
2350   fhEtaPhiPhoton  = new TH2F
2351   ("hEtaPhiPhoton","#eta vs #phi",netabins,etamin,etamax,nphibins,phimin,phimax);
2352   fhEtaPhiPhoton->SetYTitle("#phi (rad)");
2353   fhEtaPhiPhoton->SetXTitle("#eta");
2354   outputContainer->Add(fhEtaPhiPhoton) ;
2355   if(GetMinPt() < 0.5)
2356   {
2357     fhEtaPhi05Photon  = new TH2F
2358     ("hEtaPhi05Photon","#eta vs #phi, E < 0.5",netabins,etamin,etamax,nphibins,phimin,phimax);
2359     fhEtaPhi05Photon->SetYTitle("#phi (rad)");
2360     fhEtaPhi05Photon->SetXTitle("#eta");
2361     outputContainer->Add(fhEtaPhi05Photon) ;
2362   }
2363   
2364   if(fCalorimeter=="EMCAL" && fFillEMCALBCHistograms)
2365   {
2366     fhEtaPhiPhotonEMCALBC0  = new TH2F
2367     ("hEtaPhiPhotonEMCALBC0","identified photon, E > 2 GeV, #eta vs #phi, for clusters with |time| < 25 ns, EMCAL-BC=0",netabins,etamin,etamax,nphibins,phimin,phimax);
2368     fhEtaPhiPhotonEMCALBC0->SetYTitle("#phi (rad)");
2369     fhEtaPhiPhotonEMCALBC0->SetXTitle("#eta");
2370     outputContainer->Add(fhEtaPhiPhotonEMCALBC0) ;
2371     
2372     fhEtaPhiPhotonEMCALBC1  = new TH2F
2373     ("hEtaPhiPhotonEMCALBC1","identified photon, E > 2 GeV, #eta vs #phi, for clusters with 25 < |time| < 75 ns, EMCAL-BC=1",netabins,etamin,etamax,nphibins,phimin,phimax);
2374     fhEtaPhiPhotonEMCALBC1->SetYTitle("#phi (rad)");
2375     fhEtaPhiPhotonEMCALBC1->SetXTitle("#eta");
2376     outputContainer->Add(fhEtaPhiPhotonEMCALBC1) ;
2377     
2378     fhEtaPhiPhotonEMCALBCN  = new TH2F
2379     ("hEtaPhiPhotonEMCALBCN","identified photon, E > 2 GeV, #eta vs #phi, for clusters with |time| > 75 ns, EMCAL-BC>1",netabins,etamin,etamax,nphibins,phimin,phimax);
2380     fhEtaPhiPhotonEMCALBCN->SetYTitle("#phi (rad)");
2381     fhEtaPhiPhotonEMCALBCN->SetXTitle("#eta");
2382     outputContainer->Add(fhEtaPhiPhotonEMCALBCN) ;
2383     
2384     for(Int_t i = 0; i < nTrigBC; i++)
2385     {
2386       fhEtaPhiPhotonTriggerEMCALBC[i] = new TH2F
2387       (Form("hEtaPhiPhotonTriggerEMCALBC%d",i-iBCShift),
2388        Form("photon E > 2 GeV, #eta vs #phi, PhotonTrigger EMCAL-BC=%d",i-iBCShift),
2389        netabins,etamin,etamax,nphibins,phimin,phimax);
2390       fhEtaPhiPhotonTriggerEMCALBC[i]->SetYTitle("#phi (rad)");
2391       fhEtaPhiPhotonTriggerEMCALBC[i]->SetXTitle("#eta");
2392       outputContainer->Add(fhEtaPhiPhotonTriggerEMCALBC[i]) ;
2393       
2394       fhTimePhotonTriggerEMCALBC[i] = new TH2F
2395       (Form("hTimePhotonTriggerEMCALBC%d",i-iBCShift),
2396        Form("photon time vs E of clusters, PhotonTrigger EMCAL-BC=%d",i-iBCShift),
2397        nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
2398       fhTimePhotonTriggerEMCALBC[i]->SetXTitle("E (GeV)");
2399       fhTimePhotonTriggerEMCALBC[i]->SetYTitle("time (ns)");
2400       outputContainer->Add(fhTimePhotonTriggerEMCALBC[i]);
2401       
2402       fhTimePhotonTriggerEMCALBCPileUpSPD[i] = new TH2F
2403       (Form("hTimePhotonTriggerEMCALBC%dPileUpSPD",i-iBCShift),
2404        Form("photon time vs E, PhotonTrigger EMCAL-BC=%d",i-iBCShift),
2405        nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
2406       fhTimePhotonTriggerEMCALBCPileUpSPD[i]->SetXTitle("E (GeV)");
2407       fhTimePhotonTriggerEMCALBCPileUpSPD[i]->SetYTitle("time (ns)");
2408       outputContainer->Add(fhTimePhotonTriggerEMCALBCPileUpSPD[i]);
2409       
2410       fhEtaPhiPhotonTriggerEMCALBCUM[i] = new TH2F
2411       (Form("hEtaPhiPhotonTriggerEMCALBC%d_UnMatch",i-iBCShift),
2412        Form("photon E > 2 GeV, #eta vs #phi, unmatched trigger EMCAL-BC=%d",i-iBCShift),
2413        netabins,etamin,etamax,nphibins,phimin,phimax);
2414       fhEtaPhiPhotonTriggerEMCALBCUM[i]->SetYTitle("#phi (rad)");
2415       fhEtaPhiPhotonTriggerEMCALBCUM[i]->SetXTitle("#eta");
2416       outputContainer->Add(fhEtaPhiPhotonTriggerEMCALBCUM[i]) ;
2417       
2418       fhTimePhotonTriggerEMCALBCUM[i] = new TH2F
2419       (Form("hTimePhotonTriggerEMCALBC%d_UnMatch",i-iBCShift),
2420        Form("photon time vs E, unmatched trigger EMCAL-BC=%d",i-iBCShift),
2421        nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
2422       fhTimePhotonTriggerEMCALBCUM[i]->SetXTitle("E (GeV)");
2423       fhTimePhotonTriggerEMCALBCUM[i]->SetYTitle("time (ns)");
2424       outputContainer->Add(fhTimePhotonTriggerEMCALBCUM[i]);
2425       
2426     }
2427     
2428     fhTimePhotonTriggerEMCALBC0UMReMatchOpenTime = new TH2F("hTimePhotonTriggerBC0_UnMatch_ReMatch_OpenTime",
2429                                                       "cluster time vs E of photons, no match, rematch open time",
2430                                                       nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
2431     fhTimePhotonTriggerEMCALBC0UMReMatchOpenTime->SetXTitle("E (GeV)");
2432     fhTimePhotonTriggerEMCALBC0UMReMatchOpenTime->SetYTitle("time (ns)");
2433     outputContainer->Add(fhTimePhotonTriggerEMCALBC0UMReMatchOpenTime);
2434     
2435     
2436     fhTimePhotonTriggerEMCALBC0UMReMatchCheckNeigh = new TH2F("hTimePhotonTriggerBC0_UnMatch_ReMatch_CheckNeighbours",
2437                                                         "cluster time vs E of photons, no match, rematch with neigbour parches",
2438                                                         nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
2439     fhTimePhotonTriggerEMCALBC0UMReMatchCheckNeigh->SetXTitle("E (GeV)");
2440     fhTimePhotonTriggerEMCALBC0UMReMatchCheckNeigh->SetYTitle("time (ns)");
2441     outputContainer->Add(fhTimePhotonTriggerEMCALBC0UMReMatchCheckNeigh);
2442     
2443     fhTimePhotonTriggerEMCALBC0UMReMatchBoth = new TH2F("hTimePhotonTriggerBC0_UnMatch_ReMatch_Both",
2444                                                   "cluster time vs E of photons, no match, rematch open time and neigbour",
2445                                                   nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
2446     fhTimePhotonTriggerEMCALBC0UMReMatchBoth->SetXTitle("E (GeV)");
2447     fhTimePhotonTriggerEMCALBC0UMReMatchBoth->SetYTitle("time (ns)");
2448     outputContainer->Add(fhTimePhotonTriggerEMCALBC0UMReMatchBoth);
2449
2450   }
2451   
2452   fhNLocMax = new TH2F("hNLocMax","Number of local maxima in cluster",
2453                        nptbins,ptmin,ptmax,10,0,10);
2454   fhNLocMax ->SetYTitle("N maxima");
2455   fhNLocMax ->SetXTitle("E (GeV)");
2456   outputContainer->Add(fhNLocMax) ;
2457   
2458   //Shower shape
2459   if(fFillSSHistograms)
2460   {
2461     fhLam0E  = new TH2F ("hLam0E","#lambda_{0}^{2} vs E", nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
2462     fhLam0E->SetYTitle("#lambda_{0}^{2}");
2463     fhLam0E->SetXTitle("E (GeV)");
2464     outputContainer->Add(fhLam0E);
2465     
2466     fhLam1E  = new TH2F ("hLam1E","#lambda_{1}^{2} vs E", nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
2467     fhLam1E->SetYTitle("#lambda_{1}^{2}");
2468     fhLam1E->SetXTitle("E (GeV)");
2469     outputContainer->Add(fhLam1E);
2470     
2471     fhDispE  = new TH2F ("hDispE"," dispersion^{2} vs E", nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
2472     fhDispE->SetYTitle("D^{2}");
2473     fhDispE->SetXTitle("E (GeV) ");
2474     outputContainer->Add(fhDispE);
2475     
2476     if(!fRejectTrackMatch)
2477     {
2478       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);
2479       fhLam0ETM->SetYTitle("#lambda_{0}^{2}");
2480       fhLam0ETM->SetXTitle("E (GeV)");
2481       outputContainer->Add(fhLam0ETM);
2482       
2483       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);
2484       fhLam1ETM->SetYTitle("#lambda_{1}^{2}");
2485       fhLam1ETM->SetXTitle("E (GeV)");
2486       outputContainer->Add(fhLam1ETM);
2487       
2488       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);
2489       fhDispETM->SetYTitle("D^{2}");
2490       fhDispETM->SetXTitle("E (GeV) ");
2491       outputContainer->Add(fhDispETM);
2492     }
2493     
2494     if(fCalorimeter == "EMCAL")
2495     {
2496       fhLam0ETRD  = new TH2F ("hLam0ETRD","#lambda_{0}^{2} vs E, EMCAL SM covered by TRD", nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
2497       fhLam0ETRD->SetYTitle("#lambda_{0}^{2}");
2498       fhLam0ETRD->SetXTitle("E (GeV)");
2499       outputContainer->Add(fhLam0ETRD);
2500       
2501       fhLam1ETRD  = new TH2F ("hLam1ETRD","#lambda_{1}^{2} vs E, EMCAL SM covered by TRD", nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
2502       fhLam1ETRD->SetYTitle("#lambda_{1}^{2}");
2503       fhLam1ETRD->SetXTitle("E (GeV)");
2504       outputContainer->Add(fhLam1ETRD);
2505       
2506       fhDispETRD  = new TH2F ("hDispETRD"," dispersion^{2} vs E, EMCAL SM covered by TRD", nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
2507       fhDispETRD->SetYTitle("Dispersion^{2}");
2508       fhDispETRD->SetXTitle("E (GeV) ");
2509       outputContainer->Add(fhDispETRD);
2510       
2511       if(!fRejectTrackMatch)
2512       {
2513         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);
2514         fhLam0ETMTRD->SetYTitle("#lambda_{0}^{2}");
2515         fhLam0ETMTRD->SetXTitle("E (GeV)");
2516         outputContainer->Add(fhLam0ETMTRD);
2517         
2518         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);
2519         fhLam1ETMTRD->SetYTitle("#lambda_{1}^{2}");
2520         fhLam1ETMTRD->SetXTitle("E (GeV)");
2521         outputContainer->Add(fhLam1ETMTRD);
2522         
2523         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);
2524         fhDispETMTRD->SetYTitle("Dispersion^{2}");
2525         fhDispETMTRD->SetXTitle("E (GeV) ");
2526         outputContainer->Add(fhDispETMTRD);
2527       }
2528     }
2529     
2530     if(!fFillOnlySimpleSSHisto)
2531     {
2532       fhNCellsLam0LowE  = new TH2F ("hNCellsLam0LowE","N_{cells} in cluster vs #lambda_{0}^{2}, E < 2 GeV", nbins,nmin, nmax, ssbins,ssmin,ssmax);
2533       fhNCellsLam0LowE->SetXTitle("N_{cells}");
2534       fhNCellsLam0LowE->SetYTitle("#lambda_{0}^{2}");
2535       outputContainer->Add(fhNCellsLam0LowE);
2536       
2537       fhNCellsLam0HighE  = new TH2F ("hNCellsLam0HighE","N_{cells} in cluster vs #lambda_{0}^{2}, E > 2 GeV", nbins,nmin, nmax, ssbins,ssmin,ssmax);
2538       fhNCellsLam0HighE->SetXTitle("N_{cells}");
2539       fhNCellsLam0HighE->SetYTitle("#lambda_{0}^{2}");
2540       outputContainer->Add(fhNCellsLam0HighE);
2541       
2542       fhNCellsLam1LowE  = new TH2F ("hNCellsLam1LowE","N_{cells} in cluster vs #lambda_{1}^{2}, E < 2 GeV", nbins,nmin, nmax, ssbins,ssmin,ssmax);
2543       fhNCellsLam1LowE->SetXTitle("N_{cells}");
2544       fhNCellsLam1LowE->SetYTitle("#lambda_{0}^{2}");
2545       outputContainer->Add(fhNCellsLam1LowE);
2546       
2547       fhNCellsLam1HighE  = new TH2F ("hNCellsLam1HighE","N_{cells} in cluster vs #lambda_{1}^{2}, E > 2 GeV", nbins,nmin, nmax, ssbins,ssmin,ssmax);
2548       fhNCellsLam1HighE->SetXTitle("N_{cells}");
2549       fhNCellsLam1HighE->SetYTitle("#lambda_{0}^{2}");
2550       outputContainer->Add(fhNCellsLam1HighE);
2551       
2552       fhNCellsDispLowE  = new TH2F ("hNCellsDispLowE","N_{cells} in cluster vs dispersion^{2}, E < 2 GeV", nbins,nmin, nmax, ssbins,ssmin,ssmax);
2553       fhNCellsDispLowE->SetXTitle("N_{cells}");
2554       fhNCellsDispLowE->SetYTitle("D^{2}");
2555       outputContainer->Add(fhNCellsDispLowE);
2556       
2557       fhNCellsDispHighE  = new TH2F ("hNCellsDispHighE","N_{cells} in cluster vs dispersion^{2}, E < 2 GeV", nbins,nmin, nmax, ssbins,ssmin,ssmax);
2558       fhNCellsDispHighE->SetXTitle("N_{cells}");
2559       fhNCellsDispHighE->SetYTitle("D^{2}");
2560       outputContainer->Add(fhNCellsDispHighE);
2561       
2562       fhEtaLam0LowE  = new TH2F ("hEtaLam0LowE","#eta vs #lambda_{0}^{2}, E < 2 GeV", netabins,etamin,etamax, ssbins,ssmin,ssmax);
2563       fhEtaLam0LowE->SetYTitle("#lambda_{0}^{2}");
2564       fhEtaLam0LowE->SetXTitle("#eta");
2565       outputContainer->Add(fhEtaLam0LowE);
2566       
2567       fhPhiLam0LowE  = new TH2F ("hPhiLam0LowE","#phi vs #lambda_{0}^{2}, E < 2 GeV", nphibins,phimin,phimax, ssbins,ssmin,ssmax);
2568       fhPhiLam0LowE->SetYTitle("#lambda_{0}^{2}");
2569       fhPhiLam0LowE->SetXTitle("#phi");
2570       outputContainer->Add(fhPhiLam0LowE);
2571       
2572       fhEtaLam0HighE  = new TH2F ("hEtaLam0HighE","#eta vs #lambda_{0}^{2}, E > 2 GeV", netabins,etamin,etamax, ssbins,ssmin,ssmax);
2573       fhEtaLam0HighE->SetYTitle("#lambda_{0}^{2}");
2574       fhEtaLam0HighE->SetXTitle("#eta");
2575       outputContainer->Add(fhEtaLam0HighE);
2576       
2577       fhPhiLam0HighE  = new TH2F ("hPhiLam0HighE","#phi vs #lambda_{0}^{2}, E > 2 GeV", nphibins,phimin,phimax, ssbins,ssmin,ssmax);
2578       fhPhiLam0HighE->SetYTitle("#lambda_{0}^{2}");
2579       fhPhiLam0HighE->SetXTitle("#phi");
2580       outputContainer->Add(fhPhiLam0HighE);
2581       
2582       fhLam1Lam0LowE  = new TH2F ("hLam1Lam0LowE","#lambda_{0}^{2} vs #lambda_{1}^{2} in cluster of E < 2 GeV",  ssbins,ssmin,ssmax, ssbins,ssmin,ssmax);
2583       fhLam1Lam0LowE->SetYTitle("#lambda_{0}^{2}");
2584       fhLam1Lam0LowE->SetXTitle("#lambda_{1}^{2}");
2585       outputContainer->Add(fhLam1Lam0LowE);
2586       
2587       fhLam1Lam0HighE  = new TH2F ("hLam1Lam0HighE","#lambda_{0}^{2} vs #lambda_{1}^{2} in cluster of E > 2 GeV",  ssbins,ssmin,ssmax, ssbins,ssmin,ssmax);
2588       fhLam1Lam0HighE->SetYTitle("#lambda_{0}^{2}");
2589       fhLam1Lam0HighE->SetXTitle("#lambda_{1}^{2}");
2590       outputContainer->Add(fhLam1Lam0HighE);
2591       
2592       fhLam0DispLowE  = new TH2F ("hLam0DispLowE","#lambda_{0}^{2} vs dispersion^{2} in cluster of E < 2 GeV",  ssbins,ssmin,ssmax, ssbins,ssmin,ssmax);
2593       fhLam0DispLowE->SetXTitle("#lambda_{0}^{2}");
2594       fhLam0DispLowE->SetYTitle("D^{2}");
2595       outputContainer->Add(fhLam0DispLowE);
2596       
2597       fhLam0DispHighE  = new TH2F ("hLam0DispHighE","#lambda_{0}^{2} vs dispersion^{2} in cluster of E > 2 GeV",  ssbins,ssmin,ssmax, ssbins,ssmin,ssmax);
2598       fhLam0DispHighE->SetXTitle("#lambda_{0}^{2}");
2599       fhLam0DispHighE->SetYTitle("D^{2}");
2600       outputContainer->Add(fhLam0DispHighE);
2601       
2602       fhDispLam1LowE  = new TH2F ("hDispLam1LowE","Dispersion^{2} vs #lambda_{1}^{2} in cluster of E < 2 GeV",  ssbins,ssmin,ssmax, ssbins,ssmin,ssmax);
2603       fhDispLam1LowE->SetXTitle("D^{2}");
2604       fhDispLam1LowE->SetYTitle("#lambda_{1}^{2}");
2605       outputContainer->Add(fhDispLam1LowE);
2606       
2607       fhDispLam1HighE  = new TH2F ("hDispLam1HighE","Dispersion^{2} vs #lambda_{1^{2}} in cluster of E > 2 GeV",  ssbins,ssmin,ssmax, ssbins,ssmin,ssmax);
2608       fhDispLam1HighE->SetXTitle("D^{2}");
2609       fhDispLam1HighE->SetYTitle("#lambda_{1}^{2}");
2610       outputContainer->Add(fhDispLam1HighE);
2611       
2612       if(fCalorimeter == "EMCAL")
2613       {
2614         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);
2615         fhDispEtaE->SetXTitle("E (GeV)");
2616         fhDispEtaE->SetYTitle("#sigma^{2}_{#eta #eta}");
2617         outputContainer->Add(fhDispEtaE);
2618         
2619         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);
2620         fhDispPhiE->SetXTitle("E (GeV)");
2621         fhDispPhiE->SetYTitle("#sigma^{2}_{#phi #phi}");
2622         outputContainer->Add(fhDispPhiE);
2623         
2624         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);
2625         fhSumEtaE->SetXTitle("E (GeV)");
2626         fhSumEtaE->SetYTitle("#delta^{2}_{#eta #eta}");
2627         outputContainer->Add(fhSumEtaE);
2628         
2629         fhSumPhiE  = new TH2F ("hSumPhiE","#delta^{2}_{#phi #phi} = #Sigma w_{i}(#phi_{i})^{2}/ #Sigma w_{i} - <#phi>^{2} vs E",
2630                                nptbins,ptmin,ptmax, ssbins,ssmin,ssmax);
2631         fhSumPhiE->SetXTitle("E (GeV)");
2632         fhSumPhiE->SetYTitle("#delta^{2}_{#phi #phi}");
2633         outputContainer->Add(fhSumPhiE);
2634         
2635         fhSumEtaPhiE  = new TH2F ("hSumEtaPhiE","#delta^{2}_{#eta #phi} = #Sigma w_{i}(#phi_{i} #eta_{i} ) / #Sigma w_{i} - <#phi><#eta> vs E",
2636                                   nptbins,ptmin,ptmax, 2*ssbins,-ssmax,ssmax);
2637         fhSumEtaPhiE->SetXTitle("E (GeV)");
2638         fhSumEtaPhiE->SetYTitle("#delta^{2}_{#eta #phi}");
2639         outputContainer->Add(fhSumEtaPhiE);
2640         
2641         fhDispEtaPhiDiffE  = new TH2F ("hDispEtaPhiDiffE","#sigma^{2}_{#phi #phi} - #sigma^{2}_{#eta #eta} vs E",
2642                                        nptbins,ptmin,ptmax,200, -10,10);
2643         fhDispEtaPhiDiffE->SetXTitle("E (GeV)");
2644         fhDispEtaPhiDiffE->SetYTitle("#sigma^{2}_{#phi #phi}-#sigma^{2}_{#eta #eta}");
2645         outputContainer->Add(fhDispEtaPhiDiffE);
2646         
2647         fhSphericityE  = new TH2F ("hSphericityE","(#sigma^{2}_{#phi #phi} - #sigma^{2}_{#eta #eta}) / (#sigma^{2}_{#eta #eta} + #sigma^{2}_{#phi #phi}) vs E",
2648                                    nptbins,ptmin,ptmax, 200, -1,1);
2649         fhSphericityE->SetXTitle("E (GeV)");
2650         fhSphericityE->SetYTitle("s = (#sigma^{2}_{#phi #phi} - #sigma^{2}_{#eta #eta}) / (#sigma^{2}_{#eta #eta} + #sigma^{2}_{#phi #phi})");
2651         outputContainer->Add(fhSphericityE);
2652         
2653         fhDispSumEtaDiffE  = new TH2F ("hDispSumEtaDiffE","#sigma^{2}_{#eta #eta} - #delta^{2}_{#eta #eta} / average vs E",  nptbins,ptmin,ptmax, 200,-0.01,0.01);
2654         fhDispSumEtaDiffE->SetXTitle("E (GeV)");
2655         fhDispSumEtaDiffE->SetYTitle("#sigma^{2}_{#eta #eta} - #delta^{2}_{#eta #eta} / average");
2656         outputContainer->Add(fhDispSumEtaDiffE);
2657         
2658         fhDispSumPhiDiffE  = new TH2F ("hDispSumPhiDiffE","#sigma^{2}_{#phi #phi} - #delta^{2}_{#phi #phi} / average vs E",  nptbins,ptmin,ptmax, 200,-0.01,0.01);
2659         fhDispSumPhiDiffE->SetXTitle("E (GeV)");
2660         fhDispSumPhiDiffE->SetYTitle("#sigma^{2}_{#phi #phi} - #delta^{2}_{#phi #phi} / average");
2661         outputContainer->Add(fhDispSumPhiDiffE);
2662         
2663         for(Int_t i = 0; i < 7; i++)
2664         {
2665           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]),
2666                                           ssbins,ssmin,ssmax , ssbins,ssmin,ssmax);
2667           fhDispEtaDispPhi[i]->SetXTitle("#sigma^{2}_{#eta #eta}");
2668           fhDispEtaDispPhi[i]->SetYTitle("#sigma^{2}_{#phi #phi}");
2669           outputContainer->Add(fhDispEtaDispPhi[i]);
2670           
2671           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]),
2672                                           ssbins,ssmin,ssmax , ssbins,ssmin,ssmax);
2673           fhLambda0DispEta[i]->SetXTitle("#lambda^{2}_{0}");
2674           fhLambda0DispEta[i]->SetYTitle("#sigma^{2}_{#eta #eta}");
2675           outputContainer->Add(fhLambda0DispEta[i]);
2676           
2677           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]),
2678                                           ssbins,ssmin,ssmax , ssbins,ssmin,ssmax);
2679           fhLambda0DispPhi[i]->SetXTitle("#lambda^{2}_{0}");
2680           fhLambda0DispPhi[i]->SetYTitle("#sigma^{2}_{#phi #phi}");
2681           outputContainer->Add(fhLambda0DispPhi[i]);
2682         }
2683       }
2684     }
2685   } // Shower shape
2686   
2687   // Track Matching
2688   
2689   if(fFillTMHisto)
2690   {
2691     TString cutTM [] = {"NoCut",""};
2692     
2693     for(Int_t i = 0; i < 2; i++)
2694     {
2695       fhTrackMatchedDEta[i]  = new TH2F
2696       (Form("hTrackMatchedDEta%s",cutTM[i].Data()),
2697        Form("d#eta of cluster-track vs cluster energy, %s",cutTM[i].Data()),
2698        nptbins,ptmin,ptmax,nresetabins,resetamin,resetamax);
2699       fhTrackMatchedDEta[i]->SetYTitle("d#eta");
2700       fhTrackMatchedDEta[i]->SetXTitle("E_{cluster} (GeV)");
2701       
2702       fhTrackMatchedDPhi[i]  = new TH2F
2703       (Form("hTrackMatchedDPhi%s",cutTM[i].Data()),
2704        Form("d#phi of cluster-track vs cluster energy, %s",cutTM[i].Data()),
2705        nptbins,ptmin,ptmax,nresphibins,resphimin,resphimax);
2706       fhTrackMatchedDPhi[i]->SetYTitle("d#phi (rad)");
2707       fhTrackMatchedDPhi[i]->SetXTitle("E_{cluster} (GeV)");
2708       
2709       fhTrackMatchedDEtaDPhi[i]  = new TH2F
2710       (Form("hTrackMatchedDEtaDPhi%s",cutTM[i].Data()),
2711        Form("d#eta vs d#phi of cluster-track vs cluster energy, %s",cutTM[i].Data()),
2712        nresetabins,resetamin,resetamax,nresphibins,resphimin,resphimax);
2713       fhTrackMatchedDEtaDPhi[i]->SetYTitle("d#phi (rad)");
2714       fhTrackMatchedDEtaDPhi[i]->SetXTitle("d#eta");
2715       
2716       fhTrackMatchedDEtaPos[i]  = new TH2F
2717       (Form("hTrackMatchedDEtaPos%s",cutTM[i].Data()),
2718        Form("d#eta of cluster-track vs cluster energy, %s",cutTM[i].Data()),
2719        nptbins,ptmin,ptmax,nresetabins,resetamin,resetamax);
2720       fhTrackMatchedDEtaPos[i]->SetYTitle("d#eta");
2721       fhTrackMatchedDEtaPos[i]->SetXTitle("E_{cluster} (GeV)");
2722       
2723       fhTrackMatchedDPhiPos[i]  = new TH2F
2724       (Form("hTrackMatchedDPhiPos%s",cutTM[i].Data()),
2725        Form("d#phi of cluster-track vs cluster energy, %s",cutTM[i].Data()),
2726        nptbins,ptmin,ptmax,nresphibins,resphimin,resphimax);
2727       fhTrackMatchedDPhiPos[i]->SetYTitle("d#phi (rad)");
2728       fhTrackMatchedDPhiPos[i]->SetXTitle("E_{cluster} (GeV)");
2729       
2730       fhTrackMatchedDEtaDPhiPos[i]  = new TH2F
2731       (Form("hTrackMatchedDEtaDPhiPos%s",cutTM[i].Data()),
2732        Form("d#eta vs d#phi of cluster-track vs cluster energy, %s",cutTM[i].Data()),
2733        nresetabins,resetamin,resetamax,nresphibins,resphimin,resphimax);
2734       fhTrackMatchedDEtaDPhiPos[i]->SetYTitle("d#phi (rad)");
2735       fhTrackMatchedDEtaDPhiPos[i]->SetXTitle("d#eta");
2736       
2737       fhTrackMatchedDEtaNeg[i]  = new TH2F
2738       (Form("hTrackMatchedDEtaNeg%s",cutTM[i].Data()),
2739        Form("d#eta of cluster-track vs cluster energy, %s",cutTM[i].Data()),
2740        nptbins,ptmin,ptmax,nresetabins,resetamin,resetamax);
2741       fhTrackMatchedDEtaNeg[i]->SetYTitle("d#eta");
2742       fhTrackMatchedDEtaNeg[i]->SetXTitle("E_{cluster} (GeV)");
2743       
2744       fhTrackMatchedDPhiNeg[i]  = new TH2F
2745       (Form("hTrackMatchedDPhiNeg%s",cutTM[i].Data()),
2746        Form("d#phi of cluster-track vs cluster energy, %s",cutTM[i].Data()),
2747        nptbins,ptmin,ptmax,nresphibins,resphimin,resphimax);
2748       fhTrackMatchedDPhiNeg[i]->SetYTitle("d#phi (rad)");
2749       fhTrackMatchedDPhiNeg[i]->SetXTitle("E_{cluster} (GeV)");
2750       
2751       fhTrackMatchedDEtaDPhiNeg[i]  = new TH2F
2752       (Form("hTrackMatchedDEtaDPhiNeg%s",cutTM[i].Data()),
2753        Form("d#eta vs d#phi of cluster-track vs cluster energy, %s",cutTM[i].Data()),
2754        nresetabins,resetamin,resetamax,nresphibins,resphimin,resphimax);
2755       fhTrackMatchedDEtaDPhiNeg[i]->SetYTitle("d#phi (rad)");
2756       fhTrackMatchedDEtaDPhiNeg[i]->SetXTitle("d#eta");
2757       
2758       fhdEdx[i]  = new TH2F (Form("hdEdx%s",cutTM[i].Data()),Form("matched track <dE/dx> vs cluster E, %s",cutTM[i].Data()),
2759                              nptbins,ptmin,ptmax,ndedxbins, dedxmin, dedxmax);
2760       fhdEdx[i]->SetXTitle("E (GeV)");
2761       fhdEdx[i]->SetYTitle("<dE/dx>");
2762       
2763       fhEOverP[i]  = new TH2F (Form("hEOverP%s",cutTM[i].Data()),Form("matched track E/p vs cluster E, %s",cutTM[i].Data()),
2764                                nptbins,ptmin,ptmax,nPoverEbins,pOverEmin,pOverEmax);
2765       fhEOverP[i]->SetXTitle("E (GeV)");
2766       fhEOverP[i]->SetYTitle("E/p");
2767       
2768       outputContainer->Add(fhTrackMatchedDEta[i]) ;
2769       outputContainer->Add(fhTrackMatchedDPhi[i]) ;
2770       outputContainer->Add(fhTrackMatchedDEtaDPhi[i]) ;
2771       outputContainer->Add(fhTrackMatchedDEtaPos[i]) ;
2772       outputContainer->Add(fhTrackMatchedDPhiPos[i]) ;
2773       outputContainer->Add(fhTrackMatchedDEtaDPhiPos[i]) ;
2774       outputContainer->Add(fhTrackMatchedDEtaNeg[i]) ;
2775       outputContainer->Add(fhTrackMatchedDPhiNeg[i]) ;
2776       outputContainer->Add(fhTrackMatchedDEtaDPhiNeg[i]) ;
2777       outputContainer->Add(fhdEdx[i]);
2778       outputContainer->Add(fhEOverP[i]);
2779       
2780       if(fCalorimeter=="EMCAL")
2781       {
2782         fhTrackMatchedDEtaTRD[i]  = new TH2F
2783         (Form("hTrackMatchedDEtaTRD%s",cutTM[i].Data()),
2784          Form("d#eta of cluster-track vs cluster energy, SM behind TRD, %s",cutTM[i].Data()),
2785          nptbins,ptmin,ptmax,nresetabins,resetamin,resetamax);
2786         fhTrackMatchedDEtaTRD[i]->SetYTitle("d#eta");
2787         fhTrackMatchedDEtaTRD[i]->SetXTitle("E_{cluster} (GeV)");
2788         
2789         fhTrackMatchedDPhiTRD[i]  = new TH2F
2790         (Form("hTrackMatchedDPhiTRD%s",cutTM[i].Data()),
2791          Form("d#phi of cluster-track vs cluster energy, SM behing TRD, %s",cutTM[i].Data()),
2792          nptbins,ptmin,ptmax,nresphibins,resphimin,resphimax);
2793         fhTrackMatchedDPhiTRD[i]->SetYTitle("d#phi (rad)");
2794         fhTrackMatchedDPhiTRD[i]->SetXTitle("E_{cluster} (GeV)");
2795         
2796         fhEOverPTRD[i]  = new TH2F
2797         (Form("hEOverPTRD%s",cutTM[i].Data()),
2798          Form("matched track E/p vs cluster E, behind TRD, %s",cutTM[i].Data()),
2799          nptbins,ptmin,ptmax,nPoverEbins,pOverEmin,pOverEmax);
2800         fhEOverPTRD[i]->SetXTitle("E (GeV)");
2801         fhEOverPTRD[i]->SetYTitle("E/p");
2802         
2803         outputContainer->Add(fhTrackMatchedDEtaTRD[i]) ;
2804         outputContainer->Add(fhTrackMatchedDPhiTRD[i]) ;
2805         outputContainer->Add(fhEOverPTRD[i]);
2806       }
2807       
2808       if(IsDataMC())
2809       {
2810         fhTrackMatchedDEtaMCNoOverlap[i]  = new TH2F
2811         (Form("hTrackMatchedDEtaMCNoOverlap%s",cutTM[i].Data()),
2812          Form("d#eta of cluster-track vs cluster energy, no other MC particles overlap %s",cutTM[i].Data()),
2813          nptbins,ptmin,ptmax,nresetabins,resetamin,resetamax);
2814         fhTrackMatchedDEtaMCNoOverlap[i]->SetYTitle("d#eta");
2815         fhTrackMatchedDEtaMCNoOverlap[i]->SetXTitle("E_{cluster} (GeV)");
2816         
2817         fhTrackMatchedDPhiMCNoOverlap[i]  = new TH2F
2818         (Form("hTrackMatchedDPhiMCNoOverlap%s",cutTM[i].Data()),
2819          Form("d#phi of cluster-track vs cluster energy, no other MC particles overlap %s",cutTM[i].Data()),
2820          nptbins,ptmin,ptmax,nresphibins,resphimin,resphimax);
2821         fhTrackMatchedDPhiMCNoOverlap[i]->SetYTitle("d#phi (rad)");
2822         fhTrackMatchedDPhiMCNoOverlap[i]->SetXTitle("E_{cluster} (GeV)");
2823         
2824         outputContainer->Add(fhTrackMatchedDEtaMCNoOverlap[i]) ;
2825         outputContainer->Add(fhTrackMatchedDPhiMCNoOverlap[i]) ;
2826         fhTrackMatchedDEtaMCOverlap[i]  = new TH2F
2827         (Form("hTrackMatchedDEtaMCOverlap%s",cutTM[i].Data()),
2828          Form("d#eta of cluster-track vs cluster energy, several MC particles overlap %s",cutTM[i].Data()),
2829          nptbins,ptmin,ptmax,nresetabins,resetamin,resetamax);
2830         fhTrackMatchedDEtaMCOverlap[i]->SetYTitle("d#eta");
2831         fhTrackMatchedDEtaMCOverlap[i]->SetXTitle("E_{cluster} (GeV)");
2832         
2833         fhTrackMatchedDPhiMCOverlap[i]  = new TH2F
2834         (Form("hTrackMatchedDPhiMCOverlap%s",cutTM[i].Data()),
2835          Form("d#phi of cluster-track vs cluster energy, several MC particles overlap %s",cutTM[i].Data()),
2836          nptbins,ptmin,ptmax,nresphibins,resphimin,resphimax);
2837         fhTrackMatchedDPhiMCOverlap[i]->SetYTitle("d#phi (rad)");
2838         fhTrackMatchedDPhiMCOverlap[i]->SetXTitle("E_{cluster} (GeV)");
2839         
2840         outputContainer->Add(fhTrackMatchedDEtaMCOverlap[i]) ;
2841         outputContainer->Add(fhTrackMatchedDPhiMCOverlap[i]) ;
2842         
2843         fhTrackMatchedDEtaMCConversion[i]  = new TH2F
2844         (Form("hTrackMatchedDEtaMCConversion%s",cutTM[i].Data()),
2845          Form("d#eta of cluster-track vs cluster energy, no other MC particles overlap appart from conversions %s",cutTM[i].Data()),
2846          nptbins,ptmin,ptmax,nresetabins,resetamin,resetamax);
2847         fhTrackMatchedDEtaMCConversion[i]->SetYTitle("d#eta");
2848         fhTrackMatchedDEtaMCConversion[i]->SetXTitle("E_{cluster} (GeV)");
2849         
2850         fhTrackMatchedDPhiMCConversion[i]  = new TH2F
2851         (Form("hTrackMatchedDPhiMCConversion%s",cutTM[i].Data()),
2852          Form("d#phi of cluster-track vs cluster energy, no other MC particles overlap appart from conversions %s",cutTM[i].Data()),
2853          nptbins,ptmin,ptmax,nresphibins,resphimin,resphimax);
2854         fhTrackMatchedDPhiMCConversion[i]->SetYTitle("d#phi (rad)");
2855         fhTrackMatchedDPhiMCConversion[i]->SetXTitle("E_{cluster} (GeV)");
2856         
2857         outputContainer->Add(fhTrackMatchedDEtaMCConversion[i]) ;
2858         outputContainer->Add(fhTrackMatchedDPhiMCConversion[i]) ;
2859         
2860         fhTrackMatchedMCParticle[i]  = new TH2F
2861         (Form("hTrackMatchedMCParticle%s",cutTM[i].Data()),
2862          Form("Origin of particle vs energy %s",cutTM[i].Data()),
2863          nptbins,ptmin,ptmax,8,0,8);
2864         fhTrackMatchedMCParticle[i]->SetXTitle("E (GeV)");
2865         //fhTrackMatchedMCParticle[i]->SetYTitle("Particle type");
2866         
2867         fhTrackMatchedMCParticle[i]->GetYaxis()->SetBinLabel(1 ,"Photon");
2868         fhTrackMatchedMCParticle[i]->GetYaxis()->SetBinLabel(2 ,"Electron");
2869         fhTrackMatchedMCParticle[i]->GetYaxis()->SetBinLabel(3 ,"Meson Merged");
2870         fhTrackMatchedMCParticle[i]->GetYaxis()->SetBinLabel(4 ,"Rest");
2871         fhTrackMatchedMCParticle[i]->GetYaxis()->SetBinLabel(5 ,"Conv. Photon");
2872         fhTrackMatchedMCParticle[i]->GetYaxis()->SetBinLabel(6 ,"Conv. Electron");
2873         fhTrackMatchedMCParticle[i]->GetYaxis()->SetBinLabel(7 ,"Conv. Merged");
2874         fhTrackMatchedMCParticle[i]->GetYaxis()->SetBinLabel(8 ,"Conv. Rest");
2875         
2876         outputContainer->Add(fhTrackMatchedMCParticle[i]);
2877       }
2878     }
2879   }
2880   
2881   if(fFillPileUpHistograms)
2882   {
2883     
2884     TString pileUpName[] = {"SPD","EMCAL","SPDOrEMCAL","SPDAndEMCAL","SPDAndNotEMCAL","EMCALAndNotSPD","NotSPDAndNotEMCAL"} ;
2885     
2886     for(Int_t i = 0 ; i < 7 ; i++)
2887     {
2888       fhPtPileUp[i]  = new TH1F(Form("hPtPileUp%s",pileUpName[i].Data()),
2889                                 Form("Cluster  p_{T} distribution, %s Pile-Up event",pileUpName[i].Data()), nptbins,ptmin,ptmax);
2890       fhPtPileUp[i]->SetXTitle("p_{T} (GeV/c)");
2891       outputContainer->Add(fhPtPileUp[i]);
2892       
2893       fhPtChargedPileUp[i]  = new TH1F(Form("hPtChargedPileUp%s",pileUpName[i].Data()),
2894                                        Form("Charged clusters p_{T} distribution, %s Pile-Up event",pileUpName[i].Data()), nptbins,ptmin,ptmax);
2895       fhPtChargedPileUp[i]->SetXTitle("p_{T} (GeV/c)");
2896       outputContainer->Add(fhPtChargedPileUp[i]);
2897       
2898       fhPtPhotonPileUp[i]  = new TH1F(Form("hPtPhotonPileUp%s",pileUpName[i].Data()),
2899                                       Form("Selected photon p_{T} distribution, %s Pile-Up event",pileUpName[i].Data()), nptbins,ptmin,ptmax);
2900       fhPtPhotonPileUp[i]->SetXTitle("p_{T} (GeV/c)");
2901       outputContainer->Add(fhPtPhotonPileUp[i]);
2902       
2903       
2904       fhClusterEFracLongTimePileUp[i]  = new TH2F(Form("hClusterEFracLongTimePileUp%s",pileUpName[i].Data()),
2905                                                   Form("Cluster E vs fraction of cluster energy from large T cells, %s Pile-Up event",pileUpName[i].Data()),
2906                                                   nptbins,ptmin,ptmax,200,0,1);
2907       fhClusterEFracLongTimePileUp[i]->SetXTitle("E (GeV)");
2908       fhClusterEFracLongTimePileUp[i]->SetYTitle("E(large time) / E");
2909       outputContainer->Add(fhClusterEFracLongTimePileUp[i]);
2910       
2911       fhClusterCellTimePileUp[i]  = new TH2F(Form("hClusterCellTimePileUp%s",pileUpName[i].Data()),
2912                                              Form("Cluster E vs cell time in cluster, %s Pile-Up event",pileUpName[i].Data()),
2913                                              nptbins,ptmin,ptmax,ntimebins,timemin,timemax);
2914       fhClusterCellTimePileUp[i]->SetXTitle("E (GeV)");
2915       fhClusterCellTimePileUp[i]->SetYTitle("t_{cell} (ns)");
2916       outputContainer->Add(fhClusterCellTimePileUp[i]);
2917       
2918       fhClusterTimeDiffPileUp[i]  = new TH2F(Form("hClusterTimeDiffPileUp%s",pileUpName[i].Data()),
2919                                              Form("Cluster E vs t_{max}-t_{cell} in cluster, %s Pile-Up event",pileUpName[i].Data()),
2920                                              nptbins,ptmin,ptmax,400,-200,200);
2921       fhClusterTimeDiffPileUp[i]->SetXTitle("E (GeV)");
2922       fhClusterTimeDiffPileUp[i]->SetYTitle("t_{max}-t_{cell} (ns)");
2923       outputContainer->Add(fhClusterTimeDiffPileUp[i]);
2924       
2925       fhClusterTimeDiffChargedPileUp[i]  = new TH2F(Form("hClusterTimeDiffChargedPileUp%s",pileUpName[i].Data()),
2926                                                     Form("Charged clusters E vs t_{max}-t_{cell} in cluster, %s Pile-Up event",pileUpName[i].Data()),
2927                                                     nptbins,ptmin,ptmax,400,-200,200);
2928       fhClusterTimeDiffChargedPileUp[i]->SetXTitle("E (GeV)");
2929       fhClusterTimeDiffChargedPileUp[i]->SetYTitle("t_{max}-t_{cell} (ns)");
2930       outputContainer->Add(fhClusterTimeDiffChargedPileUp[i]);
2931       
2932       fhClusterTimeDiffPhotonPileUp[i]  = new TH2F(Form("hClusterTimeDiffPhotonPileUp%s",pileUpName[i].Data()),
2933                                                    Form("Selected photon E vs t_{max}-t_{cell} in cluster, %s Pile-Up event",pileUpName[i].Data()),
2934                                                    nptbins,ptmin,ptmax,400,-200,200);
2935       fhClusterTimeDiffPhotonPileUp[i]->SetXTitle("E (GeV)");
2936       fhClusterTimeDiffPhotonPileUp[i]->SetYTitle("t_{max}-t_{cell} (ns)");
2937       outputContainer->Add(fhClusterTimeDiffPhotonPileUp[i]);
2938       
2939       fhLambda0PileUp[i]  = new TH2F(Form("hLambda0PileUp%s",pileUpName[i].Data()),
2940                                      Form("Cluster E vs #lambda^{2}_{0} in cluster, %s Pile-Up event",pileUpName[i].Data()),
2941                                      nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
2942       fhLambda0PileUp[i]->SetXTitle("E (GeV)");
2943       fhLambda0PileUp[i]->SetYTitle("#lambda^{2}_{0}");
2944       outputContainer->Add(fhLambda0PileUp[i]);
2945       
2946       fhLambda0ChargedPileUp[i]  = new TH2F(Form("hLambda0ChargedPileUp%s",pileUpName[i].Data()),
2947                                             Form("Charged clusters E vs #lambda^{2}_{0}in cluster, %s Pile-Up event",pileUpName[i].Data()), nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
2948       fhLambda0ChargedPileUp[i]->SetXTitle("E (GeV)");
2949       fhLambda0ChargedPileUp[i]->SetYTitle("#lambda^{2}_{0}");
2950       outputContainer->Add(fhLambda0ChargedPileUp[i]);
2951       
2952     }
2953     
2954     fhEtaPhiBC0  = new TH2F ("hEtaPhiBC0","eta-phi for clusters tof corresponding to BC=0",netabins,etamin,etamax, nphibins,phimin,phimax);
2955     fhEtaPhiBC0->SetXTitle("#eta ");
2956     fhEtaPhiBC0->SetYTitle("#phi (rad)");
2957     outputContainer->Add(fhEtaPhiBC0);
2958     
2959     fhEtaPhiBCPlus  = new TH2F ("hEtaPhiBCPlus","eta-phi for clusters tof corresponding to BC>0",netabins,etamin,etamax, nphibins,phimin,phimax);
2960     fhEtaPhiBCPlus->SetXTitle("#eta ");
2961     fhEtaPhiBCPlus->SetYTitle("#phi (rad)");
2962     outputContainer->Add(fhEtaPhiBCPlus);
2963     
2964     fhEtaPhiBCMinus  = new TH2F ("hEtaPhiBCMinus","eta-phi for clusters tof corresponding to BC<0",netabins,etamin,etamax, nphibins,phimin,phimax);
2965     fhEtaPhiBCMinus->SetXTitle("#eta ");
2966     fhEtaPhiBCMinus->SetYTitle("#phi (rad)");
2967     outputContainer->Add(fhEtaPhiBCMinus);
2968     
2969     fhEtaPhiBC0PileUpSPD  = new TH2F ("hEtaPhiBC0PileUpSPD","eta-phi for clusters tof corresponding to BC=0, SPD pile-up",netabins,etamin,etamax, nphibins,phimin,phimax);
2970     fhEtaPhiBC0PileUpSPD->SetXTitle("#eta ");
2971     fhEtaPhiBC0PileUpSPD->SetYTitle("#phi (rad)");
2972     outputContainer->Add(fhEtaPhiBC0PileUpSPD);
2973     
2974     fhEtaPhiBCPlusPileUpSPD  = new TH2F ("hEtaPhiBCPlusPileUpSPD","eta-phi for clusters tof corresponding to BC>0, SPD pile-up",netabins,etamin,etamax, nphibins,phimin,phimax);
2975     fhEtaPhiBCPlusPileUpSPD->SetXTitle("#eta ");
2976     fhEtaPhiBCPlusPileUpSPD->SetYTitle("#phi (rad)");
2977     outputContainer->Add(fhEtaPhiBCPlusPileUpSPD);
2978     
2979     fhEtaPhiBCMinusPileUpSPD  = new TH2F ("hEtaPhiBCMinusPileUpSPD","eta-phi for clusters tof corresponding to BC<0, SPD pile-up",netabins,etamin,etamax, nphibins,phimin,phimax);
2980     fhEtaPhiBCMinusPileUpSPD->SetXTitle("#eta ");
2981     fhEtaPhiBCMinusPileUpSPD->SetYTitle("#phi (rad)");
2982     outputContainer->Add(fhEtaPhiBCMinusPileUpSPD);
2983     
2984     fhTimePtNoCut  = new TH2F ("hTimePt_NoCut","time of cluster vs pT of clusters, no cut", nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
2985     fhTimePtNoCut->SetXTitle("p_{T} (GeV/c)");
2986     fhTimePtNoCut->SetYTitle("time (ns)");
2987     outputContainer->Add(fhTimePtNoCut);
2988     
2989     fhTimePtSPD  = new TH2F ("hTimePt_SPD","time of cluster vs pT of clusters, SPD cut", nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
2990     fhTimePtSPD->SetXTitle("p_{T} (GeV/c)");
2991     fhTimePtSPD->SetYTitle("time (ns)");
2992     outputContainer->Add(fhTimePtSPD);
2993     
2994     fhTimePtPhotonNoCut  = new TH2F ("hTimePtPhoton_NoCut","time of photon cluster vs pT of clusters, no cut", nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
2995     fhTimePtPhotonNoCut->SetXTitle("p_{T} (GeV/c)");
2996     fhTimePtPhotonNoCut->SetYTitle("time (ns)");
2997     outputContainer->Add(fhTimePtPhotonNoCut);
2998     
2999     fhTimePtPhotonSPD  = new TH2F ("hTimePtPhoton_SPD","time of  photon cluster vs pT of clusters, SPD cut", nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
3000     fhTimePtPhotonSPD->SetXTitle("p_{T} (GeV/c)");
3001     fhTimePtPhotonSPD->SetYTitle("time (ns)");
3002     outputContainer->Add(fhTimePtPhotonSPD);
3003         
3004     fhTimeNPileUpVertSPD  = new TH2F ("hTime_NPileUpVertSPD","time of cluster vs N pile-up SPD vertex", ntimebins,timemin,timemax,20,0,20);
3005     fhTimeNPileUpVertSPD->SetYTitle("# vertex ");
3006     fhTimeNPileUpVertSPD->SetXTitle("time (ns)");
3007     outputContainer->Add(fhTimeNPileUpVertSPD);
3008     
3009     fhTimeNPileUpVertTrack  = new TH2F ("hTime_NPileUpVertTracks","time of cluster vs N pile-up Tracks vertex", ntimebins,timemin,timemax, 20,0,20 );
3010     fhTimeNPileUpVertTrack->SetYTitle("# vertex ");
3011     fhTimeNPileUpVertTrack->SetXTitle("time (ns)");
3012     outputContainer->Add(fhTimeNPileUpVertTrack);
3013     
3014     fhTimeNPileUpVertContributors  = new TH2F ("hTime_NPileUpVertContributors","time of cluster vs N constributors to pile-up SPD vertex", ntimebins,timemin,timemax,50,0,50);
3015     fhTimeNPileUpVertContributors->SetYTitle("# vertex ");
3016     fhTimeNPileUpVertContributors->SetXTitle("time (ns)");
3017     outputContainer->Add(fhTimeNPileUpVertContributors);
3018     
3019     fhTimePileUpMainVertexZDistance  = new TH2F ("hTime_PileUpMainVertexZDistance","time of cluster vs distance in Z pile-up SPD vertex - main SPD vertex",ntimebins,timemin,timemax,100,0,50);
3020     fhTimePileUpMainVertexZDistance->SetYTitle("distance Z (cm) ");
3021     fhTimePileUpMainVertexZDistance->SetXTitle("time (ns)");
3022     outputContainer->Add(fhTimePileUpMainVertexZDistance);
3023     
3024     fhTimePileUpMainVertexZDiamond  = new TH2F ("hTime_PileUpMainVertexZDiamond","time of cluster vs distance in Z pile-up SPD vertex - z diamond",ntimebins,timemin,timemax,100,0,50);
3025     fhTimePileUpMainVertexZDiamond->SetYTitle("diamond distance Z (cm) ");
3026     fhTimePileUpMainVertexZDiamond->SetXTitle("time (ns)");
3027     outputContainer->Add(fhTimePileUpMainVertexZDiamond);
3028     
3029     TString title[] = {"no |t diff| cut","|t diff|<20 ns","|t diff|>20 ns","|t diff|>40 ns"};
3030     TString name [] = {"TDiffNoCut","TDiffSmaller20ns","TDiffLarger20ns","TDiffLarger40ns"};
3031     for(Int_t i = 0; i < 4; i++)
3032     {
3033       fhClusterMultSPDPileUp[i] = new TH2F(Form("fhClusterMultSPDPileUp_%s", name[i].Data()),
3034                                            Form("Number of clusters per pile up event with E > 0.5 and %s respect cluster max vs cluster max E ",title[i].Data()),
3035                                            nptbins,ptmin,ptmax,100,0,100);
3036       fhClusterMultSPDPileUp[i]->SetYTitle("n clusters ");
3037       fhClusterMultSPDPileUp[i]->SetXTitle("E_{cluster max} (GeV)");
3038       outputContainer->Add(fhClusterMultSPDPileUp[i]) ;
3039       
3040       fhClusterMultNoPileUp[i] = new TH2F(Form("fhClusterMultNoPileUp_%s", name[i].Data()),
3041                                           Form("Number of clusters per non pile up event with E > 0.5 and %s respect cluster max vs cluster max E ",title[i].Data()),
3042                                           nptbins,ptmin,ptmax,100,0,100);
3043       fhClusterMultNoPileUp[i]->SetYTitle("n clusters ");
3044       fhClusterMultNoPileUp[i]->SetXTitle("E_{cluster max} (GeV)");
3045       outputContainer->Add(fhClusterMultNoPileUp[i]) ;
3046     }
3047     
3048     fhPtNPileUpSPDVtx  = new TH2F ("hPt_NPileUpVertSPD","pT of cluster vs N pile-up SPD vertex",
3049                                    nptbins,ptmin,ptmax,20,0,20);
3050     fhPtNPileUpSPDVtx->SetYTitle("# vertex ");
3051     fhPtNPileUpSPDVtx->SetXTitle("p_{T} (GeV/c)");
3052     outputContainer->Add(fhPtNPileUpSPDVtx);
3053           
3054     fhPtNPileUpTrkVtx  = new TH2F ("hPt_NPileUpVertTracks","pT of cluster vs N pile-up Tracks vertex",
3055                                    nptbins,ptmin,ptmax, 20,0,20 );
3056     fhPtNPileUpTrkVtx->SetYTitle("# vertex ");
3057     fhPtNPileUpTrkVtx->SetXTitle("p_{T} (GeV/c)");
3058     outputContainer->Add(fhPtNPileUpTrkVtx);
3059     
3060     fhPtNPileUpSPDVtxTimeCut  = new TH2F ("hPt_NPileUpVertSPD_TimeCut","pT of cluster vs N pile-up SPD vertex, |tof| < 25 ns",
3061                                           nptbins,ptmin,ptmax,20,0,20);
3062     fhPtNPileUpSPDVtxTimeCut->SetYTitle("# vertex ");
3063     fhPtNPileUpSPDVtxTimeCut->SetXTitle("p_{T} (GeV/c)");
3064     outputContainer->Add(fhPtNPileUpSPDVtxTimeCut);
3065           
3066     fhPtNPileUpTrkVtxTimeCut  = new TH2F ("hPt_NPileUpVertTracks_TimeCut","pT of cluster vs N pile-up Tracks vertex, |tof| < 25 ns",
3067                                           nptbins,ptmin,ptmax, 20,0,20 );
3068     fhPtNPileUpTrkVtxTimeCut->SetYTitle("# vertex ");
3069     fhPtNPileUpTrkVtxTimeCut->SetXTitle("p_{T} (GeV/c)");
3070     outputContainer->Add(fhPtNPileUpTrkVtxTimeCut);
3071     
3072                 fhPtNPileUpSPDVtxTimeCut2  = new TH2F ("hPt_NPileUpVertSPD_TimeCut2","pT of cluster vs N pile-up SPD vertex, -25 < tof < 75 ns",
3073                                            nptbins,ptmin,ptmax,20,0,20);
3074     fhPtNPileUpSPDVtxTimeCut2->SetYTitle("# vertex ");
3075     fhPtNPileUpSPDVtxTimeCut2->SetXTitle("p_{T} (GeV/c)");
3076     outputContainer->Add(fhPtNPileUpSPDVtxTimeCut2);
3077           
3078     fhPtNPileUpTrkVtxTimeCut2  = new TH2F ("hPt_NPileUpVertTracks_TimeCut2","pT of cluster vs N pile-up Tracks vertex, -25 < tof < 75 ns",
3079                                            nptbins,ptmin,ptmax, 20,0,20 );
3080     fhPtNPileUpTrkVtxTimeCut2->SetYTitle("# vertex ");
3081     fhPtNPileUpTrkVtxTimeCut2->SetXTitle("p_{T} (GeV/c)");
3082     outputContainer->Add(fhPtNPileUpTrkVtxTimeCut2);
3083     
3084     fhPtPhotonNPileUpSPDVtx  = new TH2F ("hPtPhoton_NPileUpVertSPD","pT of cluster vs N pile-up SPD vertex",
3085                                          nptbins,ptmin,ptmax,20,0,20);
3086     fhPtPhotonNPileUpSPDVtx->SetYTitle("# vertex ");
3087     fhPtPhotonNPileUpSPDVtx->SetXTitle("p_{T} (GeV/c)");
3088     outputContainer->Add(fhPtPhotonNPileUpSPDVtx);
3089           
3090     fhPtPhotonNPileUpTrkVtx  = new TH2F ("hPtPhoton_NPileUpVertTracks","pT of cluster vs N pile-up Tracks vertex",
3091                                          nptbins,ptmin,ptmax, 20,0,20 );
3092     fhPtPhotonNPileUpTrkVtx->SetYTitle("# vertex ");
3093     fhPtPhotonNPileUpTrkVtx->SetXTitle("p_{T} (GeV/c)");
3094     outputContainer->Add(fhPtPhotonNPileUpTrkVtx);
3095           
3096     fhPtPhotonNPileUpSPDVtxTimeCut  = new TH2F ("hPtPhoton_NPileUpVertSPD_TimeCut","pT of cluster vs N pile-up SPD vertex, |tof| < 25 ns",
3097                                                 nptbins,ptmin,ptmax,20,0,20);
3098     fhPtPhotonNPileUpSPDVtxTimeCut->SetYTitle("# vertex ");
3099     fhPtPhotonNPileUpSPDVtxTimeCut->SetXTitle("p_{T} (GeV/c)");
3100     outputContainer->Add(fhPtPhotonNPileUpSPDVtxTimeCut);
3101           
3102     fhPtPhotonNPileUpTrkVtxTimeCut  = new TH2F ("hPtPhoton_NPileUpVertTracks_TimeCut","pT of cluster vs N pile-up Tracks vertex, |tof| < 25 ns",
3103                                                 nptbins,ptmin,ptmax, 20,0,20 );
3104     fhPtPhotonNPileUpTrkVtxTimeCut->SetYTitle("# vertex ");
3105     fhPtPhotonNPileUpTrkVtxTimeCut->SetXTitle("p_{T} (GeV/c)");
3106     outputContainer->Add(fhPtPhotonNPileUpTrkVtxTimeCut);
3107           
3108     fhPtPhotonNPileUpSPDVtxTimeCut2  = new TH2F ("hPtPhoton_NPileUpVertSPD_TimeCut2","pT of cluster vs N pile-up SPD vertex, -25 < tof < 75 ns",
3109                                                  nptbins,ptmin,ptmax,20,0,20);
3110     fhPtPhotonNPileUpSPDVtxTimeCut2->SetYTitle("# vertex ");
3111     fhPtPhotonNPileUpSPDVtxTimeCut2->SetXTitle("p_{T} (GeV/c)");
3112     outputContainer->Add(fhPtPhotonNPileUpSPDVtxTimeCut2);
3113           
3114     fhPtPhotonNPileUpTrkVtxTimeCut2  = new TH2F ("hPtPhoton_NPileUpVertTracks_TimeCut2","pT of cluster vs N pile-up Tracks vertex,  -25 < tof < 75 ns",
3115                                                  nptbins,ptmin,ptmax, 20,0,20 );
3116     fhPtPhotonNPileUpTrkVtxTimeCut2->SetYTitle("# vertex ");
3117     fhPtPhotonNPileUpTrkVtxTimeCut2->SetXTitle("p_{T} (GeV/c)");
3118     outputContainer->Add(fhPtPhotonNPileUpTrkVtxTimeCut2);
3119     
3120   }
3121   
3122   if(IsDataMC())
3123   {
3124     TString ptype[] = { "#gamma", "#gamma_{#pi decay}","#gamma_{other decay}", "#pi^{0}","#eta",
3125       "e^{#pm}","#gamma->e^{#pm}","hadron?","Anti-N","Anti-P",
3126       "#gamma_{prompt}","#gamma_{fragmentation}","#gamma_{ISR}","String"                                } ;
3127     
3128     TString pname[] = { "Photon","PhotonPi0Decay","PhotonOtherDecay","Pi0","Eta","Electron",
3129       "Conversion", "Hadron", "AntiNeutron","AntiProton",
3130       "PhotonPrompt","PhotonFragmentation","PhotonISR","String" } ;
3131     
3132     for(Int_t i = 0; i < fNOriginHistograms; i++)
3133     {
3134       fhMCE[i]  = new TH1F(Form("hE_MC%s",pname[i].Data()),
3135                            Form("cluster from %s : E ",ptype[i].Data()),
3136                            nptbins,ptmin,ptmax);
3137       fhMCE[i]->SetXTitle("E (GeV)");
3138       outputContainer->Add(fhMCE[i]) ;
3139       
3140       fhMCPt[i]  = new TH1F(Form("hPt_MC%s",pname[i].Data()),
3141                             Form("cluster from %s : p_{T} ",ptype[i].Data()),
3142                             nptbins,ptmin,ptmax);
3143       fhMCPt[i]->SetXTitle("p_{T} (GeV/c)");
3144       outputContainer->Add(fhMCPt[i]) ;
3145       
3146       fhMCEta[i]  = new TH2F(Form("hEta_MC%s",pname[i].Data()),
3147                              Form("cluster from %s : #eta ",ptype[i].Data()),
3148                              nptbins,ptmin,ptmax,netabins,etamin,etamax);
3149       fhMCEta[i]->SetYTitle("#eta");
3150       fhMCEta[i]->SetXTitle("E (GeV)");
3151       outputContainer->Add(fhMCEta[i]) ;
3152       
3153       fhMCPhi[i]  = new TH2F(Form("hPhi_MC%s",pname[i].Data()),
3154                              Form("cluster from %s : #phi ",ptype[i].Data()),
3155                              nptbins,ptmin,ptmax,nphibins,phimin,phimax);
3156       fhMCPhi[i]->SetYTitle("#phi (rad)");
3157       fhMCPhi[i]->SetXTitle("E (GeV)");
3158       outputContainer->Add(fhMCPhi[i]) ;
3159       
3160       
3161       fhMCDeltaE[i]  = new TH2F (Form("hDeltaE_MC%s",pname[i].Data()),
3162                                  Form("MC - Reco E from %s",pname[i].Data()),
3163                                  nptbins,ptmin,ptmax, 200,-50,50);
3164       fhMCDeltaE[i]->SetYTitle("#Delta E (GeV)");
3165       fhMCDeltaE[i]->SetXTitle("E (GeV)");
3166       outputContainer->Add(fhMCDeltaE[i]);
3167       
3168       fhMCDeltaPt[i]  = new TH2F (Form("hDeltaPt_MC%s",pname[i].Data()),
3169                                   Form("MC - Reco p_{T} from %s",pname[i].Data()),
3170                                   nptbins,ptmin,ptmax, 200,-50,50);
3171       fhMCDeltaPt[i]->SetXTitle("p_{T,rec} (GeV/c)");
3172       fhMCDeltaPt[i]->SetYTitle("#Delta p_{T} (GeV/c)");
3173       outputContainer->Add(fhMCDeltaPt[i]);
3174       
3175       fhMC2E[i]  = new TH2F (Form("h2E_MC%s",pname[i].Data()),
3176                              Form("E distribution, reconstructed vs generated from %s",pname[i].Data()),
3177                              nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
3178       fhMC2E[i]->SetXTitle("E_{rec} (GeV)");
3179       fhMC2E[i]->SetYTitle("E_{gen} (GeV)");
3180       outputContainer->Add(fhMC2E[i]);
3181       
3182       fhMC2Pt[i]  = new TH2F (Form("h2Pt_MC%s",pname[i].Data()),
3183                               Form("p_T distribution, reconstructed vs generated from %s",pname[i].Data()),
3184                               nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
3185       fhMC2Pt[i]->SetXTitle("p_{T,rec} (GeV/c)");
3186       fhMC2Pt[i]->SetYTitle("p_{T,gen} (GeV/c)");
3187       outputContainer->Add(fhMC2Pt[i]);
3188       
3189       
3190     }
3191     
3192     TString pptype[] = { "#gamma", "#gamma_{#pi decay}","#gamma_{other decay}","hadron?",
3193       "#gamma_{prompt}","#gamma_{fragmentation}","#gamma_{ISR}"} ;
3194     
3195     TString ppname[] = { "Photon","PhotonPi0Decay","PhotonOtherDecay","Hadron",
3196       "PhotonPrompt","PhotonFragmentation","PhotonISR"} ;
3197     
3198     for(Int_t i = 0; i < fNPrimaryHistograms; i++)
3199     {
3200       fhEPrimMC[i]  = new TH1F(Form("hEPrim_MC%s",ppname[i].Data()),
3201                                Form("primary photon %s : E ",pptype[i].Data()),
3202                                nptbins,ptmin,ptmax);
3203       fhEPrimMC[i]->SetXTitle("E (GeV)");
3204       outputContainer->Add(fhEPrimMC[i]) ;
3205       
3206       fhPtPrimMC[i]  = new TH1F(Form("hPtPrim_MC%s",ppname[i].Data()),
3207                                 Form("primary photon %s : p_{T} ",pptype[i].Data()),
3208                                 nptbins,ptmin,ptmax);
3209       fhPtPrimMC[i]->SetXTitle("p_{T} (GeV/c)");
3210       outputContainer->Add(fhPtPrimMC[i]) ;
3211       
3212       fhYPrimMC[i]  = new TH2F(Form("hYPrim_MC%s",ppname[i].Data()),
3213                                Form("primary photon %s : Rapidity ",pptype[i].Data()),
3214                                nptbins,ptmin,ptmax,200,-2,2);
3215       fhYPrimMC[i]->SetYTitle("Rapidity");
3216       fhYPrimMC[i]->SetXTitle("E (GeV)");
3217       outputContainer->Add(fhYPrimMC[i]) ;
3218       
3219       fhEtaPrimMC[i]  = new TH2F(Form("hEtaPrim_MC%s",ppname[i].Data()),
3220                                Form("primary photon %s : #eta",pptype[i].Data()),
3221                                nptbins,ptmin,ptmax,200,-2,2);
3222       fhEtaPrimMC[i]->SetYTitle("#eta");
3223       fhEtaPrimMC[i]->SetXTitle("E (GeV)");
3224       outputContainer->Add(fhEtaPrimMC[i]) ;
3225       
3226       fhPhiPrimMC[i]  = new TH2F(Form("hPhiPrim_MC%s",ppname[i].Data()),
3227                                  Form("primary photon %s : #phi ",pptype[i].Data()),
3228                                  nptbins,ptmin,ptmax,nphibins,0,TMath::TwoPi());
3229       fhPhiPrimMC[i]->SetYTitle("#phi (rad)");
3230       fhPhiPrimMC[i]->SetXTitle("E (GeV)");
3231       outputContainer->Add(fhPhiPrimMC[i]) ;
3232       
3233       
3234       fhEPrimMCAcc[i]  = new TH1F(Form("hEPrimAcc_MC%s",ppname[i].Data()),
3235                                   Form("primary photon %s in acceptance: E ",pptype[i].Data()),
3236                                   nptbins,ptmin,ptmax);
3237       fhEPrimMCAcc[i]->SetXTitle("E (GeV)");
3238       outputContainer->Add(fhEPrimMCAcc[i]) ;
3239       
3240       fhPtPrimMCAcc[i]  = new TH1F(Form("hPtPrimAcc_MC%s",ppname[i].Data()),
3241                                    Form("primary photon %s in acceptance: p_{T} ",pptype[i].Data()),
3242                                    nptbins,ptmin,ptmax);
3243       fhPtPrimMCAcc[i]->SetXTitle("p_{T} (GeV/c)");
3244       outputContainer->Add(fhPtPrimMCAcc[i]) ;
3245       
3246       fhYPrimMCAcc[i]  = new TH2F(Form("hYPrimAcc_MC%s",ppname[i].Data()),
3247                                   Form("primary photon %s in acceptance: Rapidity ",pptype[i].Data()),
3248                                   nptbins,ptmin,ptmax,100,-1,1);
3249       fhYPrimMCAcc[i]->SetYTitle("Rapidity");
3250       fhYPrimMCAcc[i]->SetXTitle("E (GeV)");
3251       outputContainer->Add(fhYPrimMCAcc[i]) ;
3252
3253       fhEtaPrimMCAcc[i]  = new TH2F(Form("hEtaPrimAcc_MC%s",ppname[i].Data()),
3254                                   Form("primary photon %s in acceptance: #eta ",pptype[i].Data()),
3255                                   nptbins,ptmin,ptmax,netabins,etamin,etamax);
3256       fhEtaPrimMCAcc[i]->SetYTitle("#eta");
3257       fhEtaPrimMCAcc[i]->SetXTitle("E (GeV)");
3258       outputContainer->Add(fhEtaPrimMCAcc[i]) ;
3259       
3260       fhPhiPrimMCAcc[i]  = new TH2F(Form("hPhiPrimAcc_MC%s",ppname[i].Data()),
3261                                     Form("primary photon %s in acceptance: #phi ",pptype[i].Data()),
3262                                     nptbins,ptmin,ptmax,nphibins,phimin,phimax);
3263       fhPhiPrimMCAcc[i]->SetYTitle("#phi (rad)");
3264       fhPhiPrimMCAcc[i]->SetXTitle("E (GeV)");
3265       outputContainer->Add(fhPhiPrimMCAcc[i]) ;
3266       
3267     }
3268     
3269     if(fFillSSHistograms)
3270     {
3271       TString ptypess[] = { "#gamma","hadron?","#pi^{0}","#eta","#gamma->e^{#pm}","e^{#pm}"} ;
3272       
3273       TString pnamess[] = { "Photon","Hadron","Pi0","Eta","Conversion","Electron"} ;
3274       
3275       for(Int_t i = 0; i < 6; i++)
3276       {
3277         fhMCELambda0[i]  = new TH2F(Form("hELambda0_MC%s",pnamess[i].Data()),
3278                                     Form("cluster from %s : E vs #lambda_{0}^{2}",ptypess[i].Data()),
3279                                     nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
3280         fhMCELambda0[i]->SetYTitle("#lambda_{0}^{2}");
3281         fhMCELambda0[i]->SetXTitle("E (GeV)");
3282         outputContainer->Add(fhMCELambda0[i]) ;
3283         
3284         fhMCELambda1[i]  = new TH2F(Form("hELambda1_MC%s",pnamess[i].Data()),
3285                                     Form("cluster from %s : E vs #lambda_{1}^{2}",ptypess[i].Data()),
3286                                     nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
3287         fhMCELambda1[i]->SetYTitle("#lambda_{1}^{2}");
3288         fhMCELambda1[i]->SetXTitle("E (GeV)");
3289         outputContainer->Add(fhMCELambda1[i]) ;
3290         
3291         fhMCEDispersion[i]  = new TH2F(Form("hEDispersion_MC%s",pnamess[i].Data()),
3292                                        Form("cluster from %s : E vs dispersion^{2}",ptypess[i].Data()),
3293                                        nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
3294         fhMCEDispersion[i]->SetYTitle("D^{2}");
3295         fhMCEDispersion[i]->SetXTitle("E (GeV)");
3296         outputContainer->Add(fhMCEDispersion[i]) ;
3297         
3298         fhMCNCellsE[i]  = new TH2F (Form("hNCellsE_MC%s",pnamess[i].Data()),
3299                                     Form("# of cells in cluster from %s vs E of clusters",ptypess[i].Data()),
3300                                     nptbins,ptmin,ptmax, nbins,nmin,nmax);
3301         fhMCNCellsE[i]->SetXTitle("E (GeV)");
3302         fhMCNCellsE[i]->SetYTitle("# of cells in cluster");
3303         outputContainer->Add(fhMCNCellsE[i]);
3304         
3305         fhMCMaxCellDiffClusterE[i]  = new TH2F (Form("hMaxCellDiffClusterE_MC%s",pnamess[i].Data()),
3306                                                 Form("energy vs difference of cluster energy from %s - max cell energy / cluster energy, good clusters",ptypess[i].Data()),
3307                                                 nptbins,ptmin,ptmax, 500,0,1.);
3308         fhMCMaxCellDiffClusterE[i]->SetXTitle("E_{cluster} (GeV) ");
3309         fhMCMaxCellDiffClusterE[i]->SetYTitle("(E_{cluster} - E_{cell max})/ E_{cluster}");
3310         outputContainer->Add(fhMCMaxCellDiffClusterE[i]);
3311         
3312         if(!fFillOnlySimpleSSHisto)
3313         {
3314           fhMCLambda0vsClusterMaxCellDiffE0[i]  = new TH2F(Form("hLambda0vsClusterMaxCellDiffE0_MC%s",pnamess[i].Data()),
3315                                                            Form("cluster from %s : #lambda^{2}_{0} vs fraction of energy carried by max cell, E < 2 GeV",ptypess[i].Data()),
3316                                                            ssbins,ssmin,ssmax,500,0,1.);
3317           fhMCLambda0vsClusterMaxCellDiffE0[i]->SetXTitle("#lambda_{0}^{2}");
3318           fhMCLambda0vsClusterMaxCellDiffE0[i]->SetYTitle("(E_{cluster} - E_{cell max})/ E_{cluster}");
3319           outputContainer->Add(fhMCLambda0vsClusterMaxCellDiffE0[i]) ;
3320           
3321           fhMCLambda0vsClusterMaxCellDiffE2[i]  = new TH2F(Form("hLambda0vsClusterMaxCellDiffE2_MC%s",pnamess[i].Data()),
3322                                                            Form("cluster from %s : #lambda^{2}_{0} vs fraction of energy carried by max cell, 2< E < 6 GeV",ptypess[i].Data()),
3323                                                            ssbins,ssmin,ssmax,500,0,1.);
3324           fhMCLambda0vsClusterMaxCellDiffE2[i]->SetXTitle("#lambda_{0}^{2}");
3325           fhMCLambda0vsClusterMaxCellDiffE2[i]->SetYTitle("(E_{cluster} - E_{cell max})/ E_{cluster}");
3326           outputContainer->Add(fhMCLambda0vsClusterMaxCellDiffE2[i]) ;
3327           
3328           fhMCLambda0vsClusterMaxCellDiffE6[i]  = new TH2F(Form("hLambda0vsClusterMaxCellDiffE6_MC%s",pnamess[i].Data()),
3329                                                            Form("cluster from %s : #lambda^{2}_{0} vs fraction of energy carried by max cell, E > 6 GeV",ptypess[i].Data()),
3330                                                            ssbins,ssmin,ssmax,500,0,1.);
3331           fhMCLambda0vsClusterMaxCellDiffE6[i]->SetXTitle("#lambda_{0}^{2}");
3332           fhMCLambda0vsClusterMaxCellDiffE6[i]->SetYTitle("(E_{cluster} - E_{cell max})/ E_{cluster}");
3333           outputContainer->Add(fhMCLambda0vsClusterMaxCellDiffE6[i]) ;
3334           
3335           fhMCNCellsvsClusterMaxCellDiffE0[i]  = new TH2F(Form("hNCellsvsClusterMaxCellDiffE0_MC%s",pnamess[i].Data()),
3336                                                           Form("cluster from %s : N cells in cluster vs fraction of energy carried by max cell, E < 2 GeV",ptypess[i].Data()),
3337                                                           nbins/5,nmin,nmax/5,500,0,1.);
3338           fhMCNCellsvsClusterMaxCellDiffE0[i]->SetXTitle("N cells in cluster");
3339           fhMCNCellsvsClusterMaxCellDiffE0[i]->SetYTitle("(E_{cluster} - E_{cell max})/ E_{cluster}");
3340           outputContainer->Add(fhMCNCellsvsClusterMaxCellDiffE0[i]) ;
3341           
3342           fhMCNCellsvsClusterMaxCellDiffE2[i]  = new TH2F(Form("hNCellsvsClusterMaxCellDiffE2_MC%s",pnamess[i].Data()),
3343                                                           Form("cluster from %s : N cells in cluster vs fraction of energy carried by max cell, 2< E < 6 GeV",ptypess[i].Data()),
3344                                                           nbins/5,nmin,nmax/5,500,0,1.);
3345           fhMCNCellsvsClusterMaxCellDiffE2[i]->SetXTitle("N cells in cluster");
3346           fhMCNCellsvsClusterMaxCellDiffE2[i]->SetYTitle("(E_{cluster} - E_{cell max})/ E_{cluster}");
3347           outputContainer->Add(fhMCNCellsvsClusterMaxCellDiffE2[i]) ;
3348           
3349           fhMCNCellsvsClusterMaxCellDiffE6[i]  = new TH2F(Form("hNCellsvsClusterMaxCellDiffE6_MC%s",pnamess[i].Data()),
3350                                                           Form("cluster from %s : N cells in cluster vs fraction of energy carried by max cell, E > 6 GeV",ptypess[i].Data()),
3351                                                           nbins/5,nmin,nmax/5,500,0,1.);
3352           fhMCNCellsvsClusterMaxCellDiffE6[i]->SetXTitle("N cells in cluster");
3353           fhMCNCellsvsClusterMaxCellDiffE6[i]->SetYTitle("E (GeV)");
3354           outputContainer->Add(fhMCNCellsvsClusterMaxCellDiffE6[i]) ;
3355           
3356           if(fCalorimeter=="EMCAL")
3357           {
3358             fhMCEDispEta[i]  = new TH2F (Form("hEDispEtaE_MC%s",pnamess[i].Data()),
3359                                          Form("cluster from %s : #sigma^{2}_{#eta #eta} = #Sigma w_{i}(#eta_{i} - <#eta>)^{2}/ #Sigma w_{i} vs E",ptypess[i].Data()),
3360                                          nptbins,ptmin,ptmax, ssbins,ssmin,ssmax);
3361             fhMCEDispEta[i]->SetXTitle("E (GeV)");
3362             fhMCEDispEta[i]->SetYTitle("#sigma^{2}_{#eta #eta}");
3363             outputContainer->Add(fhMCEDispEta[i]);
3364             
3365             fhMCEDispPhi[i]  = new TH2F (Form("hEDispPhiE_MC%s",pnamess[i].Data()),
3366                                          Form("cluster from %s : #sigma^{2}_{#phi #phi} = #Sigma w_{i}(#phi_{i} - <#phi>)^{2} / #Sigma w_{i} vs E",ptypess[i].Data()),
3367                                          nptbins,ptmin,ptmax, ssbins,ssmin,ssmax);
3368             fhMCEDispPhi[i]->SetXTitle("E (GeV)");
3369             fhMCEDispPhi[i]->SetYTitle("#sigma^{2}_{#phi #phi}");
3370             outputContainer->Add(fhMCEDispPhi[i]);
3371             
3372             fhMCESumEtaPhi[i]  = new TH2F (Form("hESumEtaPhiE_MC%s",pnamess[i].Data()),
3373                                            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()),
3374                                            nptbins,ptmin,ptmax, 2*ssbins,-ssmax,ssmax);
3375             fhMCESumEtaPhi[i]->SetXTitle("E (GeV)");
3376             fhMCESumEtaPhi[i]->SetYTitle("#delta^{2}_{#eta #phi}");
3377             outputContainer->Add(fhMCESumEtaPhi[i]);
3378             
3379             fhMCEDispEtaPhiDiff[i]  = new TH2F (Form("hEDispEtaPhiDiffE_MC%s",pnamess[i].Data()),
3380                                                 Form("cluster from %s : #sigma^{2}_{#phi #phi} - #sigma^{2}_{#eta #eta} vs E",ptypess[i].Data()),
3381                                                 nptbins,ptmin,ptmax,200,-10,10);
3382             fhMCEDispEtaPhiDiff[i]->SetXTitle("E (GeV)");
3383             fhMCEDispEtaPhiDiff[i]->SetYTitle("#sigma^{2}_{#phi #phi}-#sigma^{2}_{#eta #eta}");
3384             outputContainer->Add(fhMCEDispEtaPhiDiff[i]);
3385             
3386             fhMCESphericity[i]  = new TH2F (Form("hESphericity_MC%s",pnamess[i].Data()),
3387                                             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()),
3388                                             nptbins,ptmin,ptmax, 200,-1,1);
3389             fhMCESphericity[i]->SetXTitle("E (GeV)");
3390             fhMCESphericity[i]->SetYTitle("s = (#sigma^{2}_{#phi #phi} - #sigma^{2}_{#eta #eta}) / (#sigma^{2}_{#eta #eta} + #sigma^{2}_{#phi #phi})");
3391             outputContainer->Add(fhMCESphericity[i]);
3392             
3393             for(Int_t ie = 0; ie < 7; ie++)
3394             {
3395               fhMCDispEtaDispPhi[ie][i] = new TH2F (Form("hMCDispEtaDispPhi_EBin%d_MC%s",ie,pnamess[i].Data()),
3396                                                     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]),
3397                                                     ssbins,ssmin,ssmax , ssbins,ssmin,ssmax);
3398               fhMCDispEtaDispPhi[ie][i]->SetXTitle("#sigma^{2}_{#eta #eta}");
3399               fhMCDispEtaDispPhi[ie][i]->SetYTitle("#sigma^{2}_{#phi #phi}");
3400               outputContainer->Add(fhMCDispEtaDispPhi[ie][i]);
3401               
3402               fhMCLambda0DispEta[ie][i] = new TH2F (Form("hMCLambda0DispEta_EBin%d_MC%s",ie,pnamess[i].Data()),
3403                                                     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]),
3404                                                     ssbins,ssmin,ssmax , ssbins,ssmin,ssmax);
3405               fhMCLambda0DispEta[ie][i]->SetXTitle("#lambda^{2}_{0}");
3406               fhMCLambda0DispEta[ie][i]->SetYTitle("#sigma^{2}_{#phi #phi}");
3407               outputContainer->Add(fhMCLambda0DispEta[ie][i]);
3408               
3409               fhMCLambd