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