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