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