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