fixing coding violations
[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.),
f15c25da 55 fRejectTrackMatch(0), fTimeCutMin(-10000), fTimeCutMax(10000),
f66d95af 56 fNCellsCut(0), fFillSSHistograms(kFALSE),
57 fNOriginHistograms(8), fNPrimaryHistograms(4),
521636d2 58
c4a7d28a 59 // Histograms
f15c25da 60 fhNCellsE(0), fhMaxCellDiffClusterE(0), fhTimeE(0), // Control histograms
f66d95af 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){
49b5c49b 167 if(IsTrackMatched(calo,GetReader()->GetInputEvent())) {
c4a7d28a 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();
f15c25da 944 Int_t ssbins = GetHistoShowerShapeBins(); Float_t ssmax = GetHistoShowerShapeMax(); Float_t ssmin = GetHistoShowerShapeMin();
945 Int_t nbins = GetHistoNClusterCellBins(); Int_t nmax = GetHistoNClusterCellMax(); Int_t nmin = GetHistoNClusterCellMin();
946 Int_t ntimebins= GetHistoTimeBins(); Float_t timemax = GetHistoTimeMax(); Float_t timemin = GetHistoTimeMin();
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");
f15c25da 951 outputContainer->Add(fhNCellsE);
952
953 fhTimeE = new TH2F ("hTimeE","time of cluster vs E of clusters", nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
954 fhTimeE->SetXTitle("E (GeV)");
955 fhTimeE->SetYTitle("time (ns)");
956 outputContainer->Add(fhTimeE);
6175da48 957
f66d95af 958 fhMaxCellDiffClusterE = new TH2F ("hMaxCellDiffClusterE","energy vs difference of cluster energy - max cell energy / cluster energy, good clusters",
959 nptbins,ptmin,ptmax, 500,0,1.);
960 fhMaxCellDiffClusterE->SetXTitle("E_{cluster} (GeV) ");
961 fhMaxCellDiffClusterE->SetYTitle("(E_{cluster} - E_{cell max})/ E_{cluster}");
962 outputContainer->Add(fhMaxCellDiffClusterE);
963
20218aea 964 fhEPhoton = new TH1F("hEPhoton","Number of #gamma over calorimeter vs energy",nptbins,ptmin,ptmax);
965 fhEPhoton->SetYTitle("N");
966 fhEPhoton->SetXTitle("E_{#gamma}(GeV)");
967 outputContainer->Add(fhEPhoton) ;
968
969 fhPtPhoton = new TH1F("hPtPhoton","Number of #gamma over calorimeter vs p_{T}",nptbins,ptmin,ptmax);
477d6cee 970 fhPtPhoton->SetYTitle("N");
971 fhPtPhoton->SetXTitle("p_{T #gamma}(GeV/c)");
972 outputContainer->Add(fhPtPhoton) ;
973
974 fhPhiPhoton = new TH2F
20218aea 975 ("hPhiPhoton","#phi_{#gamma} vs p_{T}",nptbins,ptmin,ptmax,nphibins,phimin,phimax);
6175da48 976 fhPhiPhoton->SetYTitle("#phi (rad)");
477d6cee 977 fhPhiPhoton->SetXTitle("p_{T #gamma} (GeV/c)");
978 outputContainer->Add(fhPhiPhoton) ;
979
980 fhEtaPhoton = new TH2F
20218aea 981 ("hEtaPhoton","#eta_{#gamma} vs p_{T}",nptbins,ptmin,ptmax,netabins,etamin,etamax);
477d6cee 982 fhEtaPhoton->SetYTitle("#eta");
983 fhEtaPhoton->SetXTitle("p_{T #gamma} (GeV/c)");
984 outputContainer->Add(fhEtaPhoton) ;
985
6175da48 986 fhEtaPhiPhoton = new TH2F
987 ("hEtaPhiPhoton","#eta vs #phi",netabins,etamin,etamax,nphibins,phimin,phimax);
988 fhEtaPhiPhoton->SetYTitle("#phi (rad)");
989 fhEtaPhiPhoton->SetXTitle("#eta");
990 outputContainer->Add(fhEtaPhiPhoton) ;
20218aea 991 if(GetMinPt() < 0.5){
992 fhEtaPhi05Photon = new TH2F
993 ("hEtaPhi05Photon","#eta vs #phi, E > 0.5",netabins,etamin,etamax,nphibins,phimin,phimax);
994 fhEtaPhi05Photon->SetYTitle("#phi (rad)");
995 fhEtaPhi05Photon->SetXTitle("#eta");
996 outputContainer->Add(fhEtaPhi05Photon) ;
997 }
6175da48 998
521636d2 999 //Shower shape
1000 if(fFillSSHistograms){
1001
1002 fhLam0E = new TH2F ("hLam0E","#lambda_{0}^{2} vs E", nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
1003 fhLam0E->SetYTitle("#lambda_{0}^{2}");
1004 fhLam0E->SetXTitle("E (GeV)");
1005 outputContainer->Add(fhLam0E);
1006
1007 fhLam1E = new TH2F ("hLam1E","#lambda_{1}^{2} vs E", nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
1008 fhLam1E->SetYTitle("#lambda_{1}^{2}");
1009 fhLam1E->SetXTitle("E (GeV)");
1010 outputContainer->Add(fhLam1E);
1011
1012 fhDispE = new TH2F ("hDispE"," dispersion^{2} vs E", nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
1013 fhDispE->SetYTitle("D^{2}");
1014 fhDispE->SetXTitle("E (GeV) ");
1015 outputContainer->Add(fhDispE);
1016
521636d2 1017 if(fCalorimeter == "EMCAL"){
1018 fhLam0ETRD = new TH2F ("hLam0ETRD","#lambda_{0}^{2} vs E, EMCAL SM covered by TRD", nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
1019 fhLam0ETRD->SetYTitle("#lambda_{0}^{2}");
1020 fhLam0ETRD->SetXTitle("E (GeV)");
1021 outputContainer->Add(fhLam0ETRD);
1022
1023 fhLam1ETRD = new TH2F ("hLam1ETRD","#lambda_{1}^{2} vs E, EMCAL SM covered by TRD", nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
1024 fhLam1ETRD->SetYTitle("#lambda_{1}^{2}");
1025 fhLam1ETRD->SetXTitle("E (GeV)");
1026 outputContainer->Add(fhLam1ETRD);
1027
1028 fhDispETRD = new TH2F ("hDispETRD"," dispersion^{2} vs E, EMCAL SM covered by TRD", nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
1029 fhDispETRD->SetYTitle("Dispersion^{2}");
1030 fhDispETRD->SetXTitle("E (GeV) ");
1031 outputContainer->Add(fhDispETRD);
521636d2 1032 }
1033
d9105d92 1034 fhNCellsLam0LowE = new TH2F ("hNCellsLam0LowE","N_{cells} in cluster vs #lambda_{0}^{2}, E < 2 GeV", nbins,nmin, nmax, ssbins,ssmin,ssmax);
521636d2 1035 fhNCellsLam0LowE->SetXTitle("N_{cells}");
1036 fhNCellsLam0LowE->SetYTitle("#lambda_{0}^{2}");
1037 outputContainer->Add(fhNCellsLam0LowE);
1038
d9105d92 1039 fhNCellsLam0HighE = new TH2F ("hNCellsLam0HighE","N_{cells} in cluster vs #lambda_{0}^{2}, E > 2 GeV", nbins,nmin, nmax, ssbins,ssmin,ssmax);
521636d2 1040 fhNCellsLam0HighE->SetXTitle("N_{cells}");
1041 fhNCellsLam0HighE->SetYTitle("#lambda_{0}^{2}");
1042 outputContainer->Add(fhNCellsLam0HighE);
1043
d9105d92 1044 fhNCellsLam1LowE = new TH2F ("hNCellsLam1LowE","N_{cells} in cluster vs #lambda_{1}^{2}, E < 2 GeV", nbins,nmin, nmax, ssbins,ssmin,ssmax);
521636d2 1045 fhNCellsLam1LowE->SetXTitle("N_{cells}");
1046 fhNCellsLam1LowE->SetYTitle("#lambda_{0}^{2}");
1047 outputContainer->Add(fhNCellsLam1LowE);
1048
d9105d92 1049 fhNCellsLam1HighE = new TH2F ("hNCellsLam1HighE","N_{cells} in cluster vs #lambda_{1}^{2}, E > 2 GeV", nbins,nmin, nmax, ssbins,ssmin,ssmax);
521636d2 1050 fhNCellsLam1HighE->SetXTitle("N_{cells}");
1051 fhNCellsLam1HighE->SetYTitle("#lambda_{0}^{2}");
1052 outputContainer->Add(fhNCellsLam1HighE);
1053
d9105d92 1054 fhNCellsDispLowE = new TH2F ("hNCellsDispLowE","N_{cells} in cluster vs dispersion^{2}, E < 2 GeV", nbins,nmin, nmax, ssbins,ssmin,ssmax);
521636d2 1055 fhNCellsDispLowE->SetXTitle("N_{cells}");
1056 fhNCellsDispLowE->SetYTitle("D^{2}");
1057 outputContainer->Add(fhNCellsDispLowE);
1058
d9105d92 1059 fhNCellsDispHighE = new TH2F ("hNCellsDispHighE","N_{cells} in cluster vs dispersion^{2}, E < 2 GeV", nbins,nmin, nmax, ssbins,ssmin,ssmax);
521636d2 1060 fhNCellsDispHighE->SetXTitle("N_{cells}");
1061 fhNCellsDispHighE->SetYTitle("D^{2}");
1062 outputContainer->Add(fhNCellsDispHighE);
1063
521636d2 1064 fhEtaLam0LowE = new TH2F ("hEtaLam0LowE","#eta vs #lambda_{0}^{2}, E < 2 GeV", netabins,etamin,etamax, ssbins,ssmin,ssmax);
1065 fhEtaLam0LowE->SetYTitle("#lambda_{0}^{2}");
1066 fhEtaLam0LowE->SetXTitle("#eta");
1067 outputContainer->Add(fhEtaLam0LowE);
1068
1069 fhPhiLam0LowE = new TH2F ("hPhiLam0LowE","#phi vs #lambda_{0}^{2}, E < 2 GeV", nphibins,phimin,phimax, ssbins,ssmin,ssmax);
1070 fhPhiLam0LowE->SetYTitle("#lambda_{0}^{2}");
1071 fhPhiLam0LowE->SetXTitle("#phi");
1072 outputContainer->Add(fhPhiLam0LowE);
1073
1074 fhEtaLam0HighE = new TH2F ("hEtaLam0HighE","#eta vs #lambda_{0}^{2}, E > 2 GeV", netabins,etamin,etamax, ssbins,ssmin,ssmax);
1075 fhEtaLam0HighE->SetYTitle("#lambda_{0}^{2}");
1076 fhEtaLam0HighE->SetXTitle("#eta");
1077 outputContainer->Add(fhEtaLam0HighE);
1078
1079 fhPhiLam0HighE = new TH2F ("hPhiLam0HighE","#phi vs #lambda_{0}^{2}, E > 2 GeV", nphibins,phimin,phimax, ssbins,ssmin,ssmax);
1080 fhPhiLam0HighE->SetYTitle("#lambda_{0}^{2}");
1081 fhPhiLam0HighE->SetXTitle("#phi");
1082 outputContainer->Add(fhPhiLam0HighE);
1083
1084 fhLam1Lam0LowE = new TH2F ("hLam1Lam0LowE","#lambda_{0}^{2} vs #lambda_{1}^{2} in cluster of E < 2 GeV", ssbins,ssmin,ssmax, ssbins,ssmin,ssmax);
1085 fhLam1Lam0LowE->SetYTitle("#lambda_{0}^{2}");
1086 fhLam1Lam0LowE->SetXTitle("#lambda_{1}^{2}");
1087 outputContainer->Add(fhLam1Lam0LowE);
1088
1089 fhLam1Lam0HighE = new TH2F ("hLam1Lam0HighE","#lambda_{0}^{2} vs #lambda_{1}^{2} in cluster of E > 2 GeV", ssbins,ssmin,ssmax, ssbins,ssmin,ssmax);
1090 fhLam1Lam0HighE->SetYTitle("#lambda_{0}^{2}");
1091 fhLam1Lam0HighE->SetXTitle("#lambda_{1}^{2}");
1092 outputContainer->Add(fhLam1Lam0HighE);
1093
1094 fhLam0DispLowE = new TH2F ("hLam0DispLowE","#lambda_{0}^{2} vs dispersion^{2} in cluster of E < 2 GeV", ssbins,ssmin,ssmax, ssbins,ssmin,ssmax);
1095 fhLam0DispLowE->SetXTitle("#lambda_{0}^{2}");
1096 fhLam0DispLowE->SetYTitle("D^{2}");
1097 outputContainer->Add(fhLam0DispLowE);
1098
1099 fhLam0DispHighE = new TH2F ("hLam0DispHighE","#lambda_{0}^{2} vs dispersion^{2} in cluster of E > 2 GeV", ssbins,ssmin,ssmax, ssbins,ssmin,ssmax);
1100 fhLam0DispHighE->SetXTitle("#lambda_{0}^{2}");
1101 fhLam0DispHighE->SetYTitle("D^{2}");
1102 outputContainer->Add(fhLam0DispHighE);
1103
1104 fhDispLam1LowE = new TH2F ("hDispLam1LowE","Dispersion^{2} vs #lambda_{1}^{2} in cluster of E < 2 GeV", ssbins,ssmin,ssmax, ssbins,ssmin,ssmax);
1105 fhDispLam1LowE->SetXTitle("D^{2}");
1106 fhDispLam1LowE->SetYTitle("#lambda_{1}^{2}");
1107 outputContainer->Add(fhDispLam1LowE);
1108
1109 fhDispLam1HighE = new TH2F ("hDispLam1HighE","Dispersion^{2} vs #lambda_{1^{2}} in cluster of E > 2 GeV", ssbins,ssmin,ssmax, ssbins,ssmin,ssmax);
1110 fhDispLam1HighE->SetXTitle("D^{2}");
1111 fhDispLam1HighE->SetYTitle("#lambda_{1}^{2}");
1112 outputContainer->Add(fhDispLam1HighE);
1113
521636d2 1114 } // Shower shape
1115
6175da48 1116
477d6cee 1117 if(IsDataMC()){
123fc3bd 1118
f66d95af 1119 TString ptype[] = { "#gamma", "#gamma_{#pi decay}","#gamma_{other decay}", "#pi^{0}","#eta",
1120 "e^{#pm}","#gamma->e^{#pm}","hadron?","Anti-N","Anti-P",
1121 "#gamma_{prompt}","#gamma_{fragmentation}","#gamma_{ISR}","String" } ;
3d5d5078 1122
f66d95af 1123 TString pname[] = { "Photon","PhotonPi0Decay","PhotonOtherDecay","Pi0","Eta","Electron",
1124 "Conversion", "Hadron", "AntiNeutron","AntiProton",
1125 "PhotonPrompt","PhotonFragmentation","PhotonISR","String" } ;
521636d2 1126
f66d95af 1127 for(Int_t i = 0; i < fNOriginHistograms; i++){
521636d2 1128
3d5d5078 1129 fhMCE[i] = new TH1F(Form("hE_MC%s",pname[i].Data()),
521636d2 1130 Form("cluster from %s : E ",ptype[i].Data()),
1131 nptbins,ptmin,ptmax);
3d5d5078 1132 fhMCE[i]->SetXTitle("E (GeV)");
1133 outputContainer->Add(fhMCE[i]) ;
521636d2 1134
4c8f7c2e 1135 fhMCPt[i] = new TH1F(Form("hPt_MC%s",pname[i].Data()),
521636d2 1136 Form("cluster from %s : p_{T} ",ptype[i].Data()),
1137 nptbins,ptmin,ptmax);
4c8f7c2e 1138 fhMCPt[i]->SetXTitle("p_{T} (GeV/c)");
1139 outputContainer->Add(fhMCPt[i]) ;
521636d2 1140
4c8f7c2e 1141 fhMCEta[i] = new TH2F(Form("hEta_MC%s",pname[i].Data()),
521636d2 1142 Form("cluster from %s : #eta ",ptype[i].Data()),
1143 nptbins,ptmin,ptmax,netabins,etamin,etamax);
4c8f7c2e 1144 fhMCEta[i]->SetYTitle("#eta");
1145 fhMCEta[i]->SetXTitle("E (GeV)");
1146 outputContainer->Add(fhMCEta[i]) ;
521636d2 1147
4c8f7c2e 1148 fhMCPhi[i] = new TH2F(Form("hPhi_MC%s",pname[i].Data()),
521636d2 1149 Form("cluster from %s : #phi ",ptype[i].Data()),
1150 nptbins,ptmin,ptmax,nphibins,phimin,phimax);
4c8f7c2e 1151 fhMCPhi[i]->SetYTitle("#phi (rad)");
1152 fhMCPhi[i]->SetXTitle("E (GeV)");
1153 outputContainer->Add(fhMCPhi[i]) ;
1154
1155
d9105d92 1156 fhMCDeltaE[i] = new TH2F (Form("hDeltaE_MC%s",pname[i].Data()),
1157 Form("MC - Reco E from %s",pname[i].Data()),
1158 nptbins,ptmin,ptmax, 200,-50,50);
4c8f7c2e 1159 fhMCDeltaE[i]->SetXTitle("#Delta E (GeV)");
1160 outputContainer->Add(fhMCDeltaE[i]);
1161
d9105d92 1162 fhMCDeltaPt[i] = new TH2F (Form("hDeltaPt_MC%s",pname[i].Data()),
1163 Form("MC - Reco p_{T} from %s",pname[i].Data()),
1164 nptbins,ptmin,ptmax, 200,-50,50);
4c8f7c2e 1165 fhMCDeltaPt[i]->SetXTitle("#Delta p_{T} (GeV/c)");
1166 outputContainer->Add(fhMCDeltaPt[i]);
d9105d92 1167
4c8f7c2e 1168 fhMC2E[i] = new TH2F (Form("h2E_MC%s",pname[i].Data()),
1169 Form("E distribution, reconstructed vs generated from %s",pname[i].Data()),
1170 nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
1171 fhMC2E[i]->SetXTitle("E_{rec} (GeV)");
1172 fhMC2E[i]->SetYTitle("E_{gen} (GeV)");
1173 outputContainer->Add(fhMC2E[i]);
1174
1175 fhMC2Pt[i] = new TH2F (Form("h2Pt_MC%s",pname[i].Data()),
1176 Form("p_T distribution, reconstructed vs generated from %s",pname[i].Data()),
1177 nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
1178 fhMC2Pt[i]->SetXTitle("p_{T,rec} (GeV/c)");
1179 fhMC2Pt[i]->SetYTitle("p_{T,gen} (GeV/c)");
1180 outputContainer->Add(fhMC2Pt[i]);
1181
521636d2 1182
1183 }
3d5d5078 1184
f66d95af 1185 TString pptype[] = { "#gamma", "#gamma_{#pi decay}","#gamma_{other decay}","hadron?",
1186 "#gamma_{prompt}","#gamma_{fragmentation}","#gamma_{ISR}"} ;
1187
1188 TString ppname[] = { "Photon","PhotonPi0Decay","PhotonOtherDecay","Hadron",
1189 "PhotonPrompt","PhotonFragmentation","PhotonISR"} ;
1190
1191 for(Int_t i = 0; i < fNPrimaryHistograms; i++){
1192 fhEPrimMC[i] = new TH1F(Form("hEPrim_MC%s",ppname[i].Data()),
1193 Form("primary photon %s : E ",pptype[i].Data()),
3d5d5078 1194 nptbins,ptmin,ptmax);
1195 fhEPrimMC[i]->SetXTitle("E (GeV)");
1196 outputContainer->Add(fhEPrimMC[i]) ;
1197
f66d95af 1198 fhPtPrimMC[i] = new TH1F(Form("hPtPrim_MC%s",ppname[i].Data()),
1199 Form("primary photon %s : p_{T} ",pptype[i].Data()),
3d5d5078 1200 nptbins,ptmin,ptmax);
1201 fhPtPrimMC[i]->SetXTitle("p_{T} (GeV/c)");
1202 outputContainer->Add(fhPtPrimMC[i]) ;
1203
f66d95af 1204 fhYPrimMC[i] = new TH2F(Form("hYPrim_MC%s",ppname[i].Data()),
1205 Form("primary photon %s : Rapidity ",pptype[i].Data()),
3d5d5078 1206 nptbins,ptmin,ptmax,800,-8,8);
1207 fhYPrimMC[i]->SetYTitle("Rapidity");
1208 fhYPrimMC[i]->SetXTitle("E (GeV)");
1209 outputContainer->Add(fhYPrimMC[i]) ;
1210
f66d95af 1211 fhPhiPrimMC[i] = new TH2F(Form("hPhiPrim_MC%s",ppname[i].Data()),
1212 Form("primary photon %s : #phi ",pptype[i].Data()),
3d5d5078 1213 nptbins,ptmin,ptmax,nphibins,phimin,phimax);
1214 fhPhiPrimMC[i]->SetYTitle("#phi (rad)");
1215 fhPhiPrimMC[i]->SetXTitle("E (GeV)");
1216 outputContainer->Add(fhPhiPrimMC[i]) ;
1217
1218
f66d95af 1219 fhEPrimMCAcc[i] = new TH1F(Form("hEPrimAcc_MC%s",ppname[i].Data()),
1220 Form("primary photon %s in acceptance: E ",pptype[i].Data()),
3d5d5078 1221 nptbins,ptmin,ptmax);
1222 fhEPrimMCAcc[i]->SetXTitle("E (GeV)");
1223 outputContainer->Add(fhEPrimMCAcc[i]) ;
1224
f66d95af 1225 fhPtPrimMCAcc[i] = new TH1F(Form("hPtPrimAcc_MC%s",ppname[i].Data()),
1226 Form("primary photon %s in acceptance: p_{T} ",pptype[i].Data()),
3d5d5078 1227 nptbins,ptmin,ptmax);
1228 fhPtPrimMCAcc[i]->SetXTitle("p_{T} (GeV/c)");
1229 outputContainer->Add(fhPtPrimMCAcc[i]) ;
1230
f66d95af 1231 fhYPrimMCAcc[i] = new TH2F(Form("hYPrimAcc_MC%s",ppname[i].Data()),
1232 Form("primary photon %s in acceptance: Rapidity ",pptype[i].Data()),
3d5d5078 1233 nptbins,ptmin,ptmax,100,-1,1);
1234 fhYPrimMCAcc[i]->SetYTitle("Rapidity");
1235 fhYPrimMCAcc[i]->SetXTitle("E (GeV)");
1236 outputContainer->Add(fhYPrimMCAcc[i]) ;
1237
f66d95af 1238 fhPhiPrimMCAcc[i] = new TH2F(Form("hPhiPrimAcc_MC%s",ppname[i].Data()),
1239 Form("primary photon %s in acceptance: #phi ",pptype[i].Data()),
3d5d5078 1240 nptbins,ptmin,ptmax,nphibins,phimin,phimax);
1241 fhPhiPrimMCAcc[i]->SetYTitle("#phi (rad)");
1242 fhPhiPrimMCAcc[i]->SetXTitle("E (GeV)");
1243 outputContainer->Add(fhPhiPrimMCAcc[i]) ;
1244
1245 }
1246
521636d2 1247 if(fFillSSHistograms){
1248
3d5d5078 1249 TString ptypess[] = { "#gamma","hadron?","#pi^{0}","#eta","#gamma->e^{#pm}","e^{#pm}"} ;
1250
1251 TString pnamess[] = { "Photon","Hadron","Pi0","Eta","Conversion","Electron"} ;
1252
1253 for(Int_t i = 0; i < 6; i++){
521636d2 1254
3d5d5078 1255 fhMCELambda0[i] = new TH2F(Form("hELambda0_MC%s",pnamess[i].Data()),
1256 Form("cluster from %s : E vs #lambda_{0}^{2}",ptypess[i].Data()),
521636d2 1257 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
3d5d5078 1258 fhMCELambda0[i]->SetYTitle("#lambda_{0}^{2}");
1259 fhMCELambda0[i]->SetXTitle("E (GeV)");
1260 outputContainer->Add(fhMCELambda0[i]) ;
521636d2 1261
3d5d5078 1262 fhMCELambda1[i] = new TH2F(Form("hELambda1_MC%s",pnamess[i].Data()),
1263 Form("cluster from %s : E vs #lambda_{1}^{2}",ptypess[i].Data()),
521636d2 1264 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
3d5d5078 1265 fhMCELambda1[i]->SetYTitle("#lambda_{1}^{2}");
1266 fhMCELambda1[i]->SetXTitle("E (GeV)");
1267 outputContainer->Add(fhMCELambda1[i]) ;
7c65ad18 1268
3d5d5078 1269 fhMCEDispersion[i] = new TH2F(Form("hEDispersion_MC%s",pnamess[i].Data()),
1270 Form("cluster from %s : E vs dispersion^{2}",ptypess[i].Data()),
521636d2 1271 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
3d5d5078 1272 fhMCEDispersion[i]->SetYTitle("D^{2}");
1273 fhMCEDispersion[i]->SetXTitle("E (GeV)");
1274 outputContainer->Add(fhMCEDispersion[i]) ;
7c65ad18 1275
f66d95af 1276 fhMCNCellsE[i] = new TH2F (Form("hNCellsE_MC%s",pnamess[i].Data()),
1277 Form("# of cells in cluster from %s vs E of clusters",ptypess[i].Data()),
1278 nptbins,ptmin,ptmax, nbins,nmin,nmax);
1279 fhMCNCellsE[i]->SetXTitle("E (GeV)");
1280 fhMCNCellsE[i]->SetYTitle("# of cells in cluster");
1281 outputContainer->Add(fhMCNCellsE[i]);
1282
1283 fhMCMaxCellDiffClusterE[i] = new TH2F (Form("hMaxCellDiffClusterE_MC%s",pnamess[i].Data()),
1284 Form("energy vs difference of cluster energy from %s - max cell energy / cluster energy, good clusters",ptypess[i].Data()),
1285 nptbins,ptmin,ptmax, 500,0,1.);
1286 fhMCMaxCellDiffClusterE[i]->SetXTitle("E_{cluster} (GeV) ");
1287 fhMCMaxCellDiffClusterE[i]->SetYTitle("(E_{cluster} - E_{cell max})/ E_{cluster}");
1288 outputContainer->Add(fhMCMaxCellDiffClusterE[i]);
1289
1290 fhMCLambda0vsClusterMaxCellDiffE0[i] = new TH2F(Form("hLambda0vsClusterMaxCellDiffE0_MC%s",pnamess[i].Data()),
1291 Form("cluster from %s : #lambda^{2}_{0} vs fraction of energy carried by max cell, E < 2 GeV",ptypess[i].Data()),
1292 ssbins,ssmin,ssmax,500,0,1.);
1293 fhMCLambda0vsClusterMaxCellDiffE0[i]->SetXTitle("#lambda_{0}^{2}");
1294 fhMCLambda0vsClusterMaxCellDiffE0[i]->SetYTitle("(E_{cluster} - E_{cell max})/ E_{cluster}");
1295 outputContainer->Add(fhMCLambda0vsClusterMaxCellDiffE0[i]) ;
1296
1297 fhMCLambda0vsClusterMaxCellDiffE2[i] = new TH2F(Form("hLambda0vsClusterMaxCellDiffE2_MC%s",pnamess[i].Data()),
1298 Form("cluster from %s : #lambda^{2}_{0} vs fraction of energy carried by max cell, 2< E < 6 GeV",ptypess[i].Data()),
1299 ssbins,ssmin,ssmax,500,0,1.);
1300 fhMCLambda0vsClusterMaxCellDiffE2[i]->SetXTitle("#lambda_{0}^{2}");
1301 fhMCLambda0vsClusterMaxCellDiffE2[i]->SetYTitle("(E_{cluster} - E_{cell max})/ E_{cluster}");
1302 outputContainer->Add(fhMCLambda0vsClusterMaxCellDiffE2[i]) ;
1303
1304 fhMCLambda0vsClusterMaxCellDiffE6[i] = new TH2F(Form("hLambda0vsClusterMaxCellDiffE6_MC%s",pnamess[i].Data()),
1305 Form("cluster from %s : #lambda^{2}_{0} vs fraction of energy carried by max cell, E > 6 GeV",ptypess[i].Data()),
1306 ssbins,ssmin,ssmax,500,0,1.);
1307 fhMCLambda0vsClusterMaxCellDiffE6[i]->SetXTitle("#lambda_{0}^{2}");
1308 fhMCLambda0vsClusterMaxCellDiffE6[i]->SetYTitle("(E_{cluster} - E_{cell max})/ E_{cluster}");
1309 outputContainer->Add(fhMCLambda0vsClusterMaxCellDiffE6[i]) ;
1310
1311 fhMCNCellsvsClusterMaxCellDiffE0[i] = new TH2F(Form("hNCellsvsClusterMaxCellDiffE0_MC%s",pnamess[i].Data()),
1312 Form("cluster from %s : N cells in cluster vs fraction of energy carried by max cell, E < 2 GeV",ptypess[i].Data()),
1313 nbins/5,nmin,nmax/5,500,0,1.);
1314 fhMCNCellsvsClusterMaxCellDiffE0[i]->SetXTitle("N cells in cluster");
1315 fhMCNCellsvsClusterMaxCellDiffE0[i]->SetYTitle("(E_{cluster} - E_{cell max})/ E_{cluster}");
1316 outputContainer->Add(fhMCNCellsvsClusterMaxCellDiffE0[i]) ;
1317
1318 fhMCNCellsvsClusterMaxCellDiffE2[i] = new TH2F(Form("hNCellsvsClusterMaxCellDiffE2_MC%s",pnamess[i].Data()),
1319 Form("cluster from %s : N cells in cluster vs fraction of energy carried by max cell, 2< E < 6 GeV",ptypess[i].Data()),
1320 nbins/5,nmin,nmax/5,500,0,1.);
1321 fhMCNCellsvsClusterMaxCellDiffE2[i]->SetXTitle("N cells in cluster");
1322 fhMCNCellsvsClusterMaxCellDiffE2[i]->SetYTitle("(E_{cluster} - E_{cell max})/ E_{cluster}");
1323 outputContainer->Add(fhMCNCellsvsClusterMaxCellDiffE2[i]) ;
1324
1325 fhMCNCellsvsClusterMaxCellDiffE6[i] = new TH2F(Form("hNCellsvsClusterMaxCellDiffE6_MC%s",pnamess[i].Data()),
1326 Form("cluster from %s : N cells in cluster vs fraction of energy carried by max cell, E > 6 GeV",ptypess[i].Data()),
1327 nbins/5,nmin,nmax/5,500,0,1.);
1328 fhMCNCellsvsClusterMaxCellDiffE6[i]->SetXTitle("N cells in cluster");
1329 fhMCNCellsvsClusterMaxCellDiffE6[i]->SetYTitle("E (GeV)");
1330 outputContainer->Add(fhMCNCellsvsClusterMaxCellDiffE6[i]) ;
1331
3d5d5078 1332 }// loop
1333
1334 if(!GetReader()->IsEmbeddedClusterSelectionOn())
1335 {
1336 fhMCPhotonELambda0NoOverlap = new TH2F("hELambda0_MCPhoton_NoOverlap",
1337 "cluster from Photon : E vs #lambda_{0}^{2}",
1338 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
1339 fhMCPhotonELambda0NoOverlap->SetYTitle("#lambda_{0}^{2}");
1340 fhMCPhotonELambda0NoOverlap->SetXTitle("E (GeV)");
1341 outputContainer->Add(fhMCPhotonELambda0NoOverlap) ;
1342
3d5d5078 1343 fhMCPhotonELambda0TwoOverlap = new TH2F("hELambda0_MCPhoton_TwoOverlap",
1344 "cluster from Photon : E vs #lambda_{0}^{2}",
1345 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
1346 fhMCPhotonELambda0TwoOverlap->SetYTitle("#lambda_{0}^{2}");
1347 fhMCPhotonELambda0TwoOverlap->SetXTitle("E (GeV)");
1348 outputContainer->Add(fhMCPhotonELambda0TwoOverlap) ;
1349
3d5d5078 1350 fhMCPhotonELambda0NOverlap = new TH2F("hELambda0_MCPhoton_NOverlap",
1351 "cluster from Photon : E vs #lambda_{0}^{2}",
1352 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
1353 fhMCPhotonELambda0NOverlap->SetYTitle("#lambda_{0}^{2}");
1354 fhMCPhotonELambda0NOverlap->SetXTitle("E (GeV)");
1355 outputContainer->Add(fhMCPhotonELambda0NOverlap) ;
521636d2 1356
3d5d5078 1357 } //No embedding
1358
1359 //Fill histograms to check shape of embedded clusters
1360 if(GetReader()->IsEmbeddedClusterSelectionOn())
1361 {
1362
1363 fhEmbeddedSignalFractionEnergy = new TH2F("hEmbeddedSignal_FractionEnergy",
1364 "Energy Fraction of embedded signal versus cluster energy",
1365 nptbins,ptmin,ptmax,100,0.,1.);
1366 fhEmbeddedSignalFractionEnergy->SetYTitle("Fraction");
1367 fhEmbeddedSignalFractionEnergy->SetXTitle("E (GeV)");
1368 outputContainer->Add(fhEmbeddedSignalFractionEnergy) ;
1369
1370 fhEmbedPhotonELambda0FullSignal = new TH2F("hELambda0_EmbedPhoton_FullSignal",
1371 "cluster from Photon embedded with more than 90% energy in cluster : E vs #lambda_{0}^{2}",
1372 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
1373 fhEmbedPhotonELambda0FullSignal->SetYTitle("#lambda_{0}^{2}");
1374 fhEmbedPhotonELambda0FullSignal->SetXTitle("E (GeV)");
1375 outputContainer->Add(fhEmbedPhotonELambda0FullSignal) ;
7c65ad18 1376
3d5d5078 1377 fhEmbedPhotonELambda0MostlySignal = new TH2F("hELambda0_EmbedPhoton_MostlySignal",
1378 "cluster from Photon embedded with 50% to 90% energy in cluster : E vs #lambda_{0}^{2}",
1379 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
1380 fhEmbedPhotonELambda0MostlySignal->SetYTitle("#lambda_{0}^{2}");
1381 fhEmbedPhotonELambda0MostlySignal->SetXTitle("E (GeV)");
1382 outputContainer->Add(fhEmbedPhotonELambda0MostlySignal) ;
1383
3d5d5078 1384 fhEmbedPhotonELambda0MostlyBkg = new TH2F("hELambda0_EmbedPhoton_MostlyBkg",
1385 "cluster from Photon embedded with 10% to 50% energy in cluster : E vs #lambda_{0}^{2}",
1386 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
1387 fhEmbedPhotonELambda0MostlyBkg->SetYTitle("#lambda_{0}^{2}");
1388 fhEmbedPhotonELambda0MostlyBkg->SetXTitle("E (GeV)");
1389 outputContainer->Add(fhEmbedPhotonELambda0MostlyBkg) ;
7c65ad18 1390
3d5d5078 1391 fhEmbedPhotonELambda0FullBkg = new TH2F("hELambda0_EmbedPhoton_FullBkg",
1392 "cluster from Photonm embedded with 0% to 10% energy in cluster : E vs #lambda_{0}^{2}",
1393 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
1394 fhEmbedPhotonELambda0FullBkg->SetYTitle("#lambda_{0}^{2}");
1395 fhEmbedPhotonELambda0FullBkg->SetXTitle("E (GeV)");
1396 outputContainer->Add(fhEmbedPhotonELambda0FullBkg) ;
1397
3d5d5078 1398 fhEmbedPi0ELambda0FullSignal = new TH2F("hELambda0_EmbedPi0_FullSignal",
1399 "cluster from Pi0 embedded with more than 90% energy in cluster : E vs #lambda_{0}^{2}",
1400 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
1401 fhEmbedPi0ELambda0FullSignal->SetYTitle("#lambda_{0}^{2}");
1402 fhEmbedPi0ELambda0FullSignal->SetXTitle("E (GeV)");
1403 outputContainer->Add(fhEmbedPi0ELambda0FullSignal) ;
7c65ad18 1404
3d5d5078 1405 fhEmbedPi0ELambda0MostlySignal = new TH2F("hELambda0_EmbedPi0_MostlySignal",
1406 "cluster from Pi0 embedded with 50% to 90% energy in cluster : E vs #lambda_{0}^{2}",
1407 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
1408 fhEmbedPi0ELambda0MostlySignal->SetYTitle("#lambda_{0}^{2}");
1409 fhEmbedPi0ELambda0MostlySignal->SetXTitle("E (GeV)");
1410 outputContainer->Add(fhEmbedPi0ELambda0MostlySignal) ;
1411
3d5d5078 1412 fhEmbedPi0ELambda0MostlyBkg = new TH2F("hELambda0_EmbedPi0_MostlyBkg",
1413 "cluster from Pi0 embedded with 10% to 50% energy in cluster : E vs #lambda_{0}^{2}",
1414 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
1415 fhEmbedPi0ELambda0MostlyBkg->SetYTitle("#lambda_{0}^{2}");
1416 fhEmbedPi0ELambda0MostlyBkg->SetXTitle("E (GeV)");
1417 outputContainer->Add(fhEmbedPi0ELambda0MostlyBkg) ;
1418
3d5d5078 1419 fhEmbedPi0ELambda0FullBkg = new TH2F("hELambda0_EmbedPi0_FullBkg",
1420 "cluster from Pi0 embedded with 0% to 10% energy in cluster : E vs #lambda_{0}^{2}",
1421 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
1422 fhEmbedPi0ELambda0FullBkg->SetYTitle("#lambda_{0}^{2}");
1423 fhEmbedPi0ELambda0FullBkg->SetXTitle("E (GeV)");
1424 outputContainer->Add(fhEmbedPi0ELambda0FullBkg) ;
7c65ad18 1425
3d5d5078 1426 }// embedded histograms
1427
521636d2 1428
1429 }// Fill SS MC histograms
1430
477d6cee 1431 }//Histos with MC
0c1383b5 1432
d39cba7e 1433 //Store calo PID histograms
05d0d05d 1434 if(fRejectTrackMatch){
d39cba7e 1435 TList * caloPIDHistos = GetCaloPID()->GetCreateOutputObjects() ;
1436 for(Int_t i = 0; i < caloPIDHistos->GetEntries(); i++) {
1437 outputContainer->Add(caloPIDHistos->At(i)) ;
1438 }
1439 delete caloPIDHistos;
1440 }
1441
477d6cee 1442 return outputContainer ;
1443
1c5acb87 1444}
1445
1446//____________________________________________________________________________
6639984f 1447void AliAnaPhoton::Init()
1448{
1449
1450 //Init
1451 //Do some checks
1e86c71e 1452 if(fCalorimeter == "PHOS" && !GetReader()->IsPHOSSwitchedOn() && NewOutputAOD()){
591cc579 1453 printf("AliAnaPhoton::Init() - !!STOP: You want to use PHOS in analysis but it is not read!! \n!!Check the configuration file!!\n");
6639984f 1454 abort();
1455 }
1e86c71e 1456 else if(fCalorimeter == "EMCAL" && !GetReader()->IsEMCALSwitchedOn() && NewOutputAOD()){
591cc579 1457 printf("AliAnaPhoton::Init() - !!STOP: You want to use EMCAL in analysis but it is not read!! \n!!Check the configuration file!!\n");
6639984f 1458 abort();
1459 }
1460
49b5c49b 1461 if(GetReader()->GetDataType() == AliCaloTrackReader::kMC) GetCaloPID()->SwitchOnBayesian();
1462
6639984f 1463}
1464
6639984f 1465//____________________________________________________________________________
1c5acb87 1466void AliAnaPhoton::InitParameters()
1467{
1468
1469 //Initialize the parameters of the analysis.
a3aebfff 1470 AddToHistogramsName("AnaPhoton_");
521636d2 1471
6175da48 1472 fCalorimeter = "EMCAL" ;
1473 fMinDist = 2.;
1474 fMinDist2 = 4.;
1475 fMinDist3 = 5.;
1e86c71e 1476
4cf55759 1477 fTimeCutMin = -1;
1478 fTimeCutMax = 9999999;
6175da48 1479 fNCellsCut = 0;
2ac125bf 1480
1e86c71e 1481 fRejectTrackMatch = kTRUE ;
1e86c71e 1482
1c5acb87 1483}
1484
1485//__________________________________________________________________
1486void AliAnaPhoton::MakeAnalysisFillAOD()
1487{
f8006433 1488 //Do photon analysis and fill aods
f37fa8d2 1489
6175da48 1490 //Get the vertex
5025c139 1491 Double_t v[3] = {0,0,0}; //vertex ;
1492 GetReader()->GetVertex(v);
f8006433 1493
f37fa8d2 1494 //Select the Calorimeter of the photon
c8fe2783 1495 TObjArray * pl = 0x0;
477d6cee 1496 if(fCalorimeter == "PHOS")
be518ab0 1497 pl = GetPHOSClusters();
477d6cee 1498 else if (fCalorimeter == "EMCAL")
be518ab0 1499 pl = GetEMCALClusters();
5ae09196 1500
1501 if(!pl) {
1502 Info("MakeAnalysisFillAOD","TObjArray with %s clusters is NULL!\n",fCalorimeter.Data());
1503 return;
1504 }
521636d2 1505
6175da48 1506 //Init arrays, variables, get number of clusters
1e86c71e 1507 TLorentzVector mom, mom2 ;
1508 Int_t nCaloClusters = pl->GetEntriesFast();
20218aea 1509
6175da48 1510 if(GetDebug() > 0) printf("AliAnaPhoton::MakeAnalysisFillAOD() - input %s cluster entries %d\n", fCalorimeter.Data(), nCaloClusters);
521636d2 1511
6175da48 1512 //----------------------------------------------------
1513 // Fill AOD with PHOS/EMCAL AliAODPWG4Particle objects
1514 //----------------------------------------------------
1515 // Loop on clusters
1e86c71e 1516 for(Int_t icalo = 0; icalo < nCaloClusters; icalo++){
1517
0ae57829 1518 AliVCluster * calo = (AliVCluster*) (pl->At(icalo));
1519 //printf("calo %d, %f\n",icalo,calo->E());
521636d2 1520
f8006433 1521 //Get the index where the cluster comes, to retrieve the corresponding vertex
c8fe2783 1522 Int_t evtIndex = 0 ;
1523 if (GetMixedEvent()) {
1524 evtIndex=GetMixedEvent()->EventIndexForCaloCluster(calo->GetID()) ;
5025c139 1525 //Get the vertex and check it is not too large in z
96539743 1526 if(TMath::Abs(GetVertex(evtIndex)[2])> GetZvertexCut()) continue;
c8fe2783 1527 }
521636d2 1528
1529 //Cluster selection, not charged, with photon id and in fiducial cut
f8006433 1530 if(GetReader()->GetDataType() != AliCaloTrackReader::kMC){
1531 calo->GetMomentum(mom,GetVertex(evtIndex)) ;}//Assume that come from vertex in straight line
1532 else{
1533 Double_t vertex[]={0,0,0};
1534 calo->GetMomentum(mom,vertex) ;
1535 }
c8fe2783 1536
6175da48 1537 //--------------------------------------
1538 // Cluster selection
1539 //--------------------------------------
c4a7d28a 1540 if(!ClusterSelected(calo,mom)) continue;
6175da48 1541
1542 //----------------------------
1543 //Create AOD for analysis
1544 //----------------------------
1545 AliAODPWG4Particle aodph = AliAODPWG4Particle(mom);
1546
1547 //...............................................
1548 //Set the indeces of the original caloclusters (MC, ID), and calorimeter
1549 Int_t label = calo->GetLabel();
1550 aodph.SetLabel(label);
6175da48 1551 aodph.SetCaloLabel(calo->GetID(),-1);
1552 aodph.SetDetector(fCalorimeter);
c4a7d28a 1553 //printf("Index %d, Id %d, iaod %d\n",icalo, calo->GetID(),GetOutputAODBranch()->GetEntriesFast());
521636d2 1554
6175da48 1555 //...............................................
1556 //Set bad channel distance bit
c4a7d28a 1557 Double_t distBad=calo->GetDistanceToBadChannel() ; //Distance to bad channel
f37fa8d2 1558 if (distBad > fMinDist3) aodph.SetDistToBad(2) ;
477d6cee 1559 else if(distBad > fMinDist2) aodph.SetDistToBad(1) ;
f37fa8d2 1560 else aodph.SetDistToBad(0) ;
af7b3903 1561 //printf("DistBad %f Bit %d\n",distBad, aodph.DistToBad());
c8fe2783 1562
521636d2 1563 //--------------------------------------------------------------------------------------
1564 //Play with the MC stack if available
1565 //--------------------------------------------------------------------------------------
1566
1567 //Check origin of the candidates
1568 if(IsDataMC()){
1569 aodph.SetTag(GetMCAnalysisUtils()->CheckOrigin(calo->GetLabels(),calo->GetNLabels(),GetReader(), aodph.GetInputFileIndex()));
3d5d5078 1570
521636d2 1571 if(GetDebug() > 0)
1572 printf("AliAnaPhoton::MakeAnalysisFillAOD() - Origin of candidate, bit map %d\n",aodph.GetTag());
1573 }//Work with stack also
1574
1575 //--------------------------------------------------------------------------------------
1576 //Fill some shower shape histograms before PID is applied
1577 //--------------------------------------------------------------------------------------
1578
1579 FillShowerShapeHistograms(calo,aodph.GetTag());
6175da48 1580
1581 //-------------------------------------
f37fa8d2 1582 //PID selection or bit setting
6175da48 1583 //-------------------------------------
49b5c49b 1584
6175da48 1585 //...............................................
1586 // Data, PID check on
49b5c49b 1587 if(IsCaloPIDOn()){
1588 // Get most probable PID, 2 options check bayesian PID weights or redo PID
1589 // By default, redo PID
1590
1591 aodph.SetIdentifiedParticleType(GetCaloPID()->GetIdentifiedParticleType(fCalorimeter,mom,calo));
477d6cee 1592
21a4b1c0 1593 if(GetDebug() > 1) printf("AliAnaPhoton::MakeAnalysisFillAOD() - PDG of identified particle %d\n",aodph.GetIdentifiedParticleType());
477d6cee 1594
f37fa8d2 1595 //If cluster does not pass pid, not photon, skip it.
21a4b1c0 1596 if(aodph.GetIdentifiedParticleType() != AliCaloPID::kPhoton) continue ;
477d6cee 1597
1598 }
6175da48 1599 //...............................................
1600 // Data, PID check off
477d6cee 1601 else{
f37fa8d2 1602 //Set PID bits for later selection (AliAnaPi0 for example)
49b5c49b 1603 //GetIdentifiedParticleType already called in SetPIDBits.
1604
1605 GetCaloPID()->SetPIDBits(fCalorimeter,calo,&aodph, GetCaloUtils(),GetReader()->GetInputEvent());
1606
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
f15c25da 1712 // Control histograms
5812a064 1713 fhMaxCellDiffClusterE->Fill(ph->E(),maxCellFraction);
1714 fhNCellsE ->Fill(ph->E(),cluster->GetNCells());
f15c25da 1715 fhTimeE ->Fill(ph->E(),cluster->GetTOF()*1.e9);
5812a064 1716
6175da48 1717 //.......................................
577d9801 1718 //Play with the MC data if available
1719 if(IsDataMC()){
521636d2 1720
3d5d5078 1721 FillAcceptanceHistograms();
1722
4c8f7c2e 1723 //....................................................................
1724 // Access MC information in stack if requested, check that it exists.
1725 Int_t label =ph->GetLabel();
1726 if(label < 0) {
1727 if(GetDebug() > 1) printf("AliAnaPhoton::MakeAnalysisFillHistograms() *** bad label ***: label %d \n", label);
1728 continue;
1729 }
1730
1731 Float_t eprim = 0;
1732 Float_t ptprim = 0;
1733 if(GetReader()->ReadStack()){
1734
1735 if(label >= stack->GetNtrack()) {
1736 if(GetDebug() > 2) printf("AliAnaPhoton::MakeAnalysisFillHistograms() *** large label ***: label %d, n tracks %d \n", label, stack->GetNtrack());
1737 continue ;
1738 }
1739
1740 primary = stack->Particle(label);
1741 if(!primary){
1742 printf("AliAnaPhoton::MakeAnalysisFillHistograms() *** no primary ***: label %d \n", label);
1743 continue;
1744 }
1745 eprim = primary->Energy();
1746 ptprim = primary->Pt();
1747
1748 }
1749 else if(GetReader()->ReadAODMCParticles()){
1750 //Check which is the input
1751 if(ph->GetInputFileIndex() == 0){
1752 if(!mcparticles) continue;
1753 if(label >= mcparticles->GetEntriesFast()) {
1754 if(GetDebug() > 2) printf("AliAnaPhoton::MakeAnalysisFillHistograms() *** large label ***: label %d, n tracks %d \n",
1755 label, mcparticles->GetEntriesFast());
1756 continue ;
1757 }
1758 //Get the particle
1759 aodprimary = (AliAODMCParticle*) mcparticles->At(label);
1760
1761 }
1762
1763 if(!aodprimary){
1764 printf("AliAnaPhoton::MakeAnalysisFillHistograms() *** no primary ***: label %d \n", label);
1765 continue;
1766 }
1767
1768 eprim = aodprimary->E();
1769 ptprim = aodprimary->Pt();
1770
1771 }
1772
577d9801 1773 Int_t tag =ph->GetTag();
521636d2 1774
f66d95af 1775 if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPhoton) && fhMCE[mcPhoton])
3d5d5078 1776 {
1777 fhMCE [mcPhoton] ->Fill(ecluster);
4c8f7c2e 1778 fhMCPt [mcPhoton] ->Fill(ptcluster);
1779 fhMCPhi[mcPhoton] ->Fill(ecluster,phicluster);
1780 fhMCEta[mcPhoton] ->Fill(ecluster,etacluster);
1781
1782 fhMC2E[mcPhoton] ->Fill(ecluster, eprim);
1783 fhMC2Pt[mcPhoton] ->Fill(ptcluster, ptprim);
d9105d92 1784 fhMCDeltaE[mcPhoton] ->Fill(ecluster,eprim-ecluster);
1785 fhMCDeltaPt[mcPhoton]->Fill(ptcluster,ptprim-ptcluster);
3d5d5078 1786
f66d95af 1787 if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCConversion) && fhMCE[mcConversion])
3d5d5078 1788 {
1789 fhMCE [mcConversion] ->Fill(ecluster);
4c8f7c2e 1790 fhMCPt [mcConversion] ->Fill(ptcluster);
1791 fhMCPhi[mcConversion] ->Fill(ecluster,phicluster);
1792 fhMCEta[mcConversion] ->Fill(ecluster,etacluster);
1793
1794 fhMC2E[mcConversion] ->Fill(ecluster, eprim);
1795 fhMC2Pt[mcConversion] ->Fill(ptcluster, ptprim);
d9105d92 1796 fhMCDeltaE[mcConversion] ->Fill(ecluster,eprim-ecluster);
1797 fhMCDeltaPt[mcConversion]->Fill(ptcluster,ptprim-ptcluster);
3d5d5078 1798 }
1799
f66d95af 1800 if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPrompt) && fhMCE[mcPrompt]){
3d5d5078 1801 fhMCE [mcPrompt] ->Fill(ecluster);
4c8f7c2e 1802 fhMCPt [mcPrompt] ->Fill(ptcluster);
1803 fhMCPhi[mcPrompt] ->Fill(ecluster,phicluster);
1804 fhMCEta[mcPrompt] ->Fill(ecluster,etacluster);
1805
1806 fhMC2E[mcPrompt] ->Fill(ecluster, eprim);
1807 fhMC2Pt[mcPrompt] ->Fill(ptcluster, ptprim);
d9105d92 1808 fhMCDeltaE[mcPrompt] ->Fill(ecluster,eprim-ecluster);
1809 fhMCDeltaPt[mcPrompt]->Fill(ptcluster,ptprim-ptcluster);
4c8f7c2e 1810
3d5d5078 1811 }
f66d95af 1812 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCFragmentation)&& fhMCE[mcFragmentation])
3d5d5078 1813 {
1814 fhMCE [mcFragmentation] ->Fill(ecluster);
4c8f7c2e 1815 fhMCPt [mcFragmentation] ->Fill(ptcluster);
1816 fhMCPhi[mcFragmentation] ->Fill(ecluster,phicluster);
1817 fhMCEta[mcFragmentation] ->Fill(ecluster,etacluster);
1818
1819 fhMC2E[mcFragmentation] ->Fill(ecluster, eprim);
1820 fhMC2Pt[mcFragmentation] ->Fill(ptcluster, ptprim);
d9105d92 1821 fhMCDeltaE[mcFragmentation] ->Fill(ecluster,eprim-ecluster);
1822 fhMCDeltaPt[mcFragmentation]->Fill(ptcluster,ptprim-ptcluster);
3d5d5078 1823
1824 }
f66d95af 1825 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCISR)&& fhMCE[mcISR])
3d5d5078 1826 {
1827 fhMCE [mcISR] ->Fill(ecluster);
4c8f7c2e 1828 fhMCPt [mcISR] ->Fill(ptcluster);
1829 fhMCPhi[mcISR] ->Fill(ecluster,phicluster);
1830 fhMCEta[mcISR] ->Fill(ecluster,etacluster);
1831
1832 fhMC2E[mcISR] ->Fill(ecluster, eprim);
1833 fhMC2Pt[mcISR] ->Fill(ptcluster, ptprim);
d9105d92 1834 fhMCDeltaE[mcISR] ->Fill(ecluster,eprim-ecluster);
1835 fhMCDeltaPt[mcISR]->Fill(ptcluster,ptprim-ptcluster);
4c8f7c2e 1836
3d5d5078 1837 }
1838 else if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0Decay) &&
f66d95af 1839 !GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0) && fhMCE[mcPi0Decay])
3d5d5078 1840 {
1841 fhMCE [mcPi0Decay] ->Fill(ecluster);
4c8f7c2e 1842 fhMCPt [mcPi0Decay] ->Fill(ptcluster);
1843 fhMCPhi[mcPi0Decay] ->Fill(ecluster,phicluster);
1844 fhMCEta[mcPi0Decay] ->Fill(ecluster,etacluster);
1845
1846 fhMC2E[mcPi0Decay] ->Fill(ecluster, eprim);
1847 fhMC2Pt[mcPi0Decay] ->Fill(ptcluster, ptprim);
d9105d92 1848 fhMCDeltaE[mcPi0Decay] ->Fill(ecluster,eprim-ecluster);
1849 fhMCDeltaPt[mcPi0Decay]->Fill(ptcluster,ptprim-ptcluster);
4c8f7c2e 1850
3d5d5078 1851 }
f586f4aa 1852 else if( (GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCEtaDecay) ||
1853 GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCOtherDecay) ) && fhMCE[mcOtherDecay])
3d5d5078 1854 {
1855 fhMCE [mcOtherDecay] ->Fill(ecluster);
4c8f7c2e 1856 fhMCPt [mcOtherDecay] ->Fill(ptcluster);
1857 fhMCPhi[mcOtherDecay] ->Fill(ecluster,phicluster);
1858 fhMCEta[mcOtherDecay] ->Fill(ecluster,etacluster);
1859
1860 fhMC2E[mcOtherDecay] ->Fill(ecluster, eprim);
1861 fhMC2Pt[mcOtherDecay] ->Fill(ptcluster, ptprim);
d9105d92 1862 fhMCDeltaE[mcOtherDecay] ->Fill(ecluster,eprim-ecluster);
1863 fhMCDeltaPt[mcOtherDecay]->Fill(ptcluster,ptprim-ptcluster);
4c8f7c2e 1864
3d5d5078 1865 }
f66d95af 1866 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0) && fhMCE [mcPi0])
3d5d5078 1867 {
1868 fhMCE [mcPi0] ->Fill(ecluster);
4c8f7c2e 1869 fhMCPt [mcPi0] ->Fill(ptcluster);
1870 fhMCPhi[mcPi0] ->Fill(ecluster,phicluster);
1871 fhMCEta[mcPi0] ->Fill(ecluster,etacluster);
1872
1873 fhMC2E[mcPi0] ->Fill(ecluster, eprim);
1874 fhMC2Pt[mcPi0] ->Fill(ptcluster, ptprim);
d9105d92 1875 fhMCDeltaE[mcPi0] ->Fill(ecluster,eprim-ecluster);
1876 fhMCDeltaPt[mcPi0]->Fill(ptcluster,ptprim-ptcluster);
4c8f7c2e 1877
f66d95af 1878 }
1879 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCEta) && fhMCE[mcEta])
1880 {
1881 fhMCE [mcEta] ->Fill(ecluster);
4c8f7c2e 1882 fhMCPt [mcEta] ->Fill(ptcluster);
1883 fhMCPhi[mcEta] ->Fill(ecluster,phicluster);
1884 fhMCEta[mcEta] ->Fill(ecluster,etacluster);
1885
1886 fhMC2E[mcEta] ->Fill(ecluster, eprim);
1887 fhMC2Pt[mcEta] ->Fill(ptcluster, ptprim);
d9105d92 1888 fhMCDeltaE[mcEta] ->Fill(ecluster,eprim-ecluster);
1889 fhMCDeltaPt[mcEta]->Fill(ptcluster,ptprim-ptcluster);
4c8f7c2e 1890
f66d95af 1891 }
3d5d5078 1892 }
f66d95af 1893 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCAntiNeutron) && fhMCE[mcAntiNeutron])
3d5d5078 1894 {
1895 fhMCE [mcAntiNeutron] ->Fill(ecluster);
4c8f7c2e 1896 fhMCPt [mcAntiNeutron] ->Fill(ptcluster);
1897 fhMCPhi[mcAntiNeutron] ->Fill(ecluster,phicluster);
1898 fhMCEta[mcAntiNeutron] ->Fill(ecluster,etacluster);
3d5d5078 1899
4c8f7c2e 1900 fhMC2E[mcAntiNeutron] ->Fill(ecluster, eprim);
1901 fhMC2Pt[mcAntiNeutron] ->Fill(ptcluster, ptprim);
d9105d92 1902 fhMCDeltaE[mcAntiNeutron] ->Fill(ecluster,eprim-ecluster);
1903 fhMCDeltaPt[mcAntiNeutron]->Fill(ptcluster,ptprim-ptcluster);
4c8f7c2e 1904
3d5d5078 1905 }
f66d95af 1906 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCAntiProton) && fhMCE[mcAntiProton])
3d5d5078 1907 {
1908 fhMCE [mcAntiProton] ->Fill(ecluster);
4c8f7c2e 1909 fhMCPt [mcAntiProton] ->Fill(ptcluster);
1910 fhMCPhi[mcAntiProton] ->Fill(ecluster,phicluster);
1911 fhMCEta[mcAntiProton] ->Fill(ecluster,etacluster);
4c8f7c2e 1912
1913 fhMC2E[mcAntiProton] ->Fill(ecluster, eprim);
1914 fhMC2Pt[mcAntiProton] ->Fill(ptcluster, ptprim);
d9105d92 1915 fhMCDeltaE[mcAntiProton] ->Fill(ecluster,eprim-ecluster);
1916 fhMCDeltaPt[mcAntiProton]->Fill(ecluster,ptprim-ptcluster);
4c8f7c2e 1917
3d5d5078 1918 }
f66d95af 1919 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCElectron) && fhMCE[mcElectron])
3d5d5078 1920 {
1921 fhMCE [mcElectron] ->Fill(ecluster);
4c8f7c2e 1922 fhMCPt [mcElectron] ->Fill(ptcluster);
1923 fhMCPhi[mcElectron] ->Fill(ecluster,phicluster);
1924 fhMCEta[mcElectron] ->Fill(ecluster,etacluster);
1925
1926 fhMC2E[mcElectron] ->Fill(ecluster, eprim);
1927 fhMC2Pt[mcElectron] ->Fill(ptcluster, ptprim);
d9105d92 1928 fhMCDeltaE[mcElectron] ->Fill(ecluster,eprim-ecluster);
1929 fhMCDeltaPt[mcElectron]->Fill(ecluster,ptprim-ptcluster);
3d5d5078 1930 }
f66d95af 1931 else if( fhMCE[mcOther]){
3d5d5078 1932 fhMCE [mcOther] ->Fill(ecluster);
4c8f7c2e 1933 fhMCPt [mcOther] ->Fill(ptcluster);
1934 fhMCPhi[mcOther] ->Fill(ecluster,phicluster);
1935 fhMCEta[mcOther] ->Fill(ecluster,etacluster);
521636d2 1936
4c8f7c2e 1937 fhMC2E[mcOther] ->Fill(ecluster, eprim);
1938 fhMC2Pt[mcOther] ->Fill(ptcluster, ptprim);
d9105d92 1939 fhMCDeltaE[mcOther] ->Fill(ecluster,eprim-ecluster);
1940 fhMCDeltaPt[mcOther]->Fill(ecluster,ptprim-ptcluster);
4c8f7c2e 1941
f8006433 1942 // printf(" AliAnaPhoton::MakeAnalysisFillHistograms() - Label %d, pT %2.3f Unknown, bits set: ",
1943 // ph->GetLabel(),ph->Pt());
1944 // for(Int_t i = 0; i < 20; i++) {
1945 // if(GetMCAnalysisUtils()->CheckTagBit(tag,i)) printf(" %d, ",i);
1946 // }
1947 // printf("\n");
1948
577d9801 1949 }
521636d2 1950
577d9801 1951 }//Histograms with MC
521636d2 1952
577d9801 1953 }// aod loop
521636d2 1954
1c5acb87 1955}
1956
1957
1958//__________________________________________________________________
1959void AliAnaPhoton::Print(const Option_t * opt) const
1960{
477d6cee 1961 //Print some relevant parameters set for the analysis
1962
1963 if(! opt)
1964 return;
1965
1966 printf("**** Print %s %s ****\n", GetName(), GetTitle() ) ;
1967 AliAnaPartCorrBaseClass::Print(" ");
a3aebfff 1968
477d6cee 1969 printf("Calorimeter = %s\n", fCalorimeter.Data()) ;
1970 printf("Min Distance to Bad Channel = %2.1f\n",fMinDist);
1971 printf("Min Distance to Bad Channel 2 = %2.1f\n",fMinDist2);
1972 printf("Min Distance to Bad Channel 3 = %2.1f\n",fMinDist3);
a3aebfff 1973 printf("Reject clusters with a track matched = %d\n",fRejectTrackMatch);
4cf55759 1974 printf("Time Cut: %3.1f < TOF < %3.1f\n", fTimeCutMin, fTimeCutMax);
2ac125bf 1975 printf("Number of cells in cluster is > %d \n", fNCellsCut);
477d6cee 1976 printf(" \n") ;
1c5acb87 1977
1978}