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