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