]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWG4/PartCorrDep/AliAnaPhoton.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
open in the negative side default time cut in QA and Photon analysis, add time contro...
[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){
167 if(IsTrackMatched(calo)) {
168 if(GetDebug() > 2) printf("\t Reject track-matched clusters\n");
169 return kFALSE ;
170 }
171 else
172 if(GetDebug() > 2) printf(" Track-matching cut passed \n");
173 }// reject matched clusters
174
175 //.......................................
176 //Check Distance to Bad channel, set bit.
177 Double_t distBad=calo->GetDistanceToBadChannel() ; //Distance to bad channel
178 if(distBad < 0.) distBad=9999. ; //workout strange convension dist = -1. ;
179 if(distBad < fMinDist) {//In bad channel (PHOS cristal size 2.2x2.2 cm), EMCAL ( cell units )
180 return kFALSE ;
181 }
182 else if(GetDebug() > 2) printf("\t Bad channel cut passed %4.2f > %2.2f \n",distBad, fMinDist);
c4a7d28a 183
184 if(GetDebug() > 0)
185 printf("AliAnaPhoton::ClusterSelected() Current Event %d; After selection : E %2.2f, pT %2.2f, Ecl %2.2f, phi %2.2f, eta %2.2f\n",
186 GetReader()->GetEventNumber(),
187 mom.E(), mom.Pt(),calo->E(),mom.Phi()*TMath::RadToDeg(),mom.Eta());
188
189 //All checks passed, cluster selected
190 return kTRUE;
191
192}
193
3d5d5078 194//_____________________________________________________________
195void AliAnaPhoton::FillAcceptanceHistograms(){
196 //Fill acceptance histograms if MC data is available
197
198 if(GetReader()->ReadStack()){
199 AliStack * stack = GetMCStack();
200 if(stack){
201 for(Int_t i=0 ; i<stack->GetNtrack(); i++){
202 TParticle * prim = stack->Particle(i) ;
203 Int_t pdg = prim->GetPdgCode();
204 //printf("i %d, %s %d %s %d \n",i, stack->Particle(i)->GetName(), stack->Particle(i)->GetPdgCode(),
205 // prim->GetName(), prim->GetPdgCode());
206
207 if(pdg == 22){
208
209 // Get tag of this particle photon from fragmentation, decay, prompt ...
210 Int_t tag = GetMCAnalysisUtils()->CheckOrigin(i,GetReader(), 0);
211 if(!GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPhoton)){
212 //A conversion photon from a hadron, skip this kind of photon
213 // printf("AliAnaPhoton::FillAcceptanceHistograms() - not a photon, weird!: tag %d, conv %d, pi0 %d, hadron %d, electron %d, unk %d, muon %d,pion %d, proton %d, neutron %d, kaon %d, antiproton %d, antineutron %d\n",tag,
214 // GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCConversion),
215 // GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0),
216 // GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCOther),
217 // GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCElectron),
218 // GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCUnknown),
219 // GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCMuon),
220 // GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPion),
221 // GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCProton),
222 // GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCAntiNeutron),
223 // GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCKaon),
224 // GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCAntiProton),
225 // GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCAntiNeutron));
226
227 return;
228 }
229
230 //Get photon kinematics
231 if(prim->Energy() == TMath::Abs(prim->Pz())) continue ; //Protection against floating point exception
232
233 Double_t photonY = 0.5*TMath::Log((prim->Energy()-prim->Pz())/(prim->Energy()+prim->Pz())) ;
234 Double_t photonE = prim->Energy() ;
235 Double_t photonPt = prim->Pt() ;
236 Double_t photonPhi = TMath::RadToDeg()*prim->Phi() ;
237 if(photonPhi < 0) photonPhi+=TMath::TwoPi();
238 Double_t photonEta = prim->Eta() ;
239
240
241 //Check if photons hit the Calorimeter
242 TLorentzVector lv;
243 prim->Momentum(lv);
244 Bool_t inacceptance = kFALSE;
245 if (fCalorimeter == "PHOS"){
246 if(GetPHOSGeometry() && GetCaloUtils()->IsPHOSGeoMatrixSet()){
247 Int_t mod ;
248 Double_t x,z ;
249 if(GetPHOSGeometry()->ImpactOnEmc(prim,mod,z,x))
250 inacceptance = kTRUE;
251 if(GetDebug() > 2) printf("In %s Real acceptance? %d\n",fCalorimeter.Data(),inacceptance);
252 }
253 else{
254 if(GetFiducialCut()->IsInFiducialCut(lv,fCalorimeter))
255 inacceptance = kTRUE ;
256 if(GetDebug() > 2) printf("In %s fiducial cut acceptance? %d\n",fCalorimeter.Data(),inacceptance);
257 }
258 }
259 else if(fCalorimeter == "EMCAL" && GetCaloUtils()->IsEMCALGeoMatrixSet()){
260 if(GetEMCALGeometry()){
261
262 Int_t absID=0;
263
264 GetEMCALGeometry()->GetAbsCellIdFromEtaPhi(prim->Eta(),prim->Phi(),absID);
265
266 if( absID >= 0)
267 inacceptance = kTRUE;
268
269 // if(GetEMCALGeometry()->Impact(phot1) && GetEMCALGeometry()->Impact(phot2))
270 // inacceptance = kTRUE;
271 if(GetDebug() > 2) printf("In %s Real acceptance? %d\n",fCalorimeter.Data(),inacceptance);
272 }
273 else{
274 if(GetFiducialCut()->IsInFiducialCut(lv,fCalorimeter))
275 inacceptance = kTRUE ;
276 if(GetDebug() > 2) printf("In %s fiducial cut acceptance? %d\n",fCalorimeter.Data(),inacceptance);
277 }
278 } //In EMCAL
279
280 //Fill histograms
281
f66d95af 282 fhYPrimMC[mcPPhoton]->Fill(photonPt, photonY) ;
3d5d5078 283 if(TMath::Abs(photonY) < 1.0)
284 {
f66d95af 285 fhEPrimMC [mcPPhoton]->Fill(photonE ) ;
286 fhPtPrimMC [mcPPhoton]->Fill(photonPt) ;
287 fhPhiPrimMC[mcPPhoton]->Fill(photonE , photonPhi) ;
288 fhYPrimMC[mcPPhoton] ->Fill(photonE , photonEta) ;
3d5d5078 289 }
290 if(inacceptance){
f66d95af 291 fhEPrimMCAcc[mcPPhoton] ->Fill(photonE ) ;
292 fhPtPrimMCAcc[mcPPhoton] ->Fill(photonPt) ;
293 fhPhiPrimMCAcc[mcPPhoton]->Fill(photonE , photonPhi) ;
294 fhYPrimMCAcc[mcPPhoton] ->Fill(photonE , photonY) ;
3d5d5078 295 }//Accepted
296
297 //Origin of photon
f66d95af 298 if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPrompt) && fhEPrimMC[mcPPrompt])
3d5d5078 299 {
f66d95af 300 fhYPrimMC[mcPPrompt]->Fill(photonPt, photonY) ;
3d5d5078 301 if(TMath::Abs(photonY) < 1.0){
f66d95af 302 fhEPrimMC [mcPPrompt]->Fill(photonE ) ;
303 fhPtPrimMC [mcPPrompt]->Fill(photonPt) ;
304 fhPhiPrimMC[mcPPrompt]->Fill(photonE , photonPhi) ;
305 fhYPrimMC[mcPPrompt] ->Fill(photonE , photonEta) ;
3d5d5078 306 }
307 if(inacceptance){
f66d95af 308 fhEPrimMCAcc[mcPPrompt] ->Fill(photonE ) ;
309 fhPtPrimMCAcc[mcPPrompt] ->Fill(photonPt) ;
310 fhPhiPrimMCAcc[mcPPrompt]->Fill(photonE , photonPhi) ;
311 fhYPrimMCAcc[mcPPrompt] ->Fill(photonE , photonY) ;
3d5d5078 312 }//Accepted
313 }
f66d95af 314 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCFragmentation) && fhEPrimMC[mcPFragmentation])
3d5d5078 315 {
f66d95af 316 fhYPrimMC[mcPFragmentation]->Fill(photonPt, photonY) ;
3d5d5078 317 if(TMath::Abs(photonY) < 1.0){
f66d95af 318 fhEPrimMC [mcPFragmentation]->Fill(photonE ) ;
319 fhPtPrimMC [mcPFragmentation]->Fill(photonPt) ;
320 fhPhiPrimMC[mcPFragmentation]->Fill(photonE , photonPhi) ;
321 fhYPrimMC[mcPFragmentation] ->Fill(photonE , photonEta) ;
3d5d5078 322 }
323 if(inacceptance){
f66d95af 324 fhEPrimMCAcc[mcPFragmentation] ->Fill(photonE ) ;
325 fhPtPrimMCAcc[mcPFragmentation] ->Fill(photonPt) ;
326 fhPhiPrimMCAcc[mcPFragmentation]->Fill(photonE , photonPhi) ;
327 fhYPrimMCAcc[mcPFragmentation] ->Fill(photonE , photonY) ;
3d5d5078 328 }//Accepted
329 }
f66d95af 330 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCISR) && fhEPrimMC[mcPISR])
3d5d5078 331 {
f66d95af 332 fhYPrimMC[mcPISR]->Fill(photonPt, photonY) ;
3d5d5078 333 if(TMath::Abs(photonY) < 1.0){
f66d95af 334 fhEPrimMC [mcPISR]->Fill(photonE ) ;
335 fhPtPrimMC [mcPISR]->Fill(photonPt) ;
336 fhPhiPrimMC[mcPISR]->Fill(photonE , photonPhi) ;
337 fhYPrimMC[mcPISR]->Fill(photonE , photonEta) ;
3d5d5078 338 }
339 if(inacceptance){
f66d95af 340 fhEPrimMCAcc[mcPISR] ->Fill(photonE ) ;
341 fhPtPrimMCAcc[mcPISR] ->Fill(photonPt) ;
342 fhPhiPrimMCAcc[mcPISR]->Fill(photonE , photonPhi) ;
343 fhYPrimMCAcc[mcPISR] ->Fill(photonE , photonY) ;
3d5d5078 344 }//Accepted
345 }
f66d95af 346 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0Decay)&& fhEPrimMC[mcPPi0Decay])
3d5d5078 347 {
f66d95af 348 fhYPrimMC[mcPPi0Decay]->Fill(photonPt, photonY) ;
3d5d5078 349 if(TMath::Abs(photonY) < 1.0){
f66d95af 350 fhEPrimMC [mcPPi0Decay]->Fill(photonE ) ;
351 fhPtPrimMC [mcPPi0Decay]->Fill(photonPt) ;
352 fhPhiPrimMC[mcPPi0Decay]->Fill(photonE , photonPhi) ;
353 fhYPrimMC[mcPPi0Decay] ->Fill(photonE , photonEta) ;
3d5d5078 354 }
355 if(inacceptance){
f66d95af 356 fhEPrimMCAcc[mcPPi0Decay] ->Fill(photonE ) ;
357 fhPtPrimMCAcc[mcPPi0Decay] ->Fill(photonPt) ;
358 fhPhiPrimMCAcc[mcPPi0Decay]->Fill(photonE , photonPhi) ;
359 fhYPrimMCAcc[mcPPi0Decay] ->Fill(photonE , photonY) ;
3d5d5078 360 }//Accepted
361 }
f586f4aa 362 else if( (GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCEtaDecay) ||
363 GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCOtherDecay)) && fhEPrimMC[mcPOtherDecay])
3d5d5078 364 {
f66d95af 365 fhYPrimMC[mcPOtherDecay]->Fill(photonPt, photonY) ;
3d5d5078 366 if(TMath::Abs(photonY) < 1.0){
f66d95af 367 fhEPrimMC [mcPOtherDecay]->Fill(photonE ) ;
368 fhPtPrimMC [mcPOtherDecay]->Fill(photonPt) ;
369 fhPhiPrimMC[mcPOtherDecay]->Fill(photonE , photonPhi) ;
370 fhYPrimMC[mcPOtherDecay] ->Fill(photonE , photonEta) ;
3d5d5078 371 }
372 if(inacceptance){
f66d95af 373 fhEPrimMCAcc[mcPOtherDecay] ->Fill(photonE ) ;
374 fhPtPrimMCAcc[mcPOtherDecay] ->Fill(photonPt) ;
375 fhPhiPrimMCAcc[mcPOtherDecay]->Fill(photonE , photonPhi) ;
376 fhYPrimMCAcc[mcPOtherDecay] ->Fill(photonE , photonY) ;
3d5d5078 377 }//Accepted
378 }
f66d95af 379 else if(fhEPrimMC[mcPOther])
3d5d5078 380 {
f66d95af 381 fhYPrimMC[mcPOther]->Fill(photonPt, photonY) ;
3d5d5078 382 if(TMath::Abs(photonY) < 1.0){
f66d95af 383 fhEPrimMC [mcPOther]->Fill(photonE ) ;
384 fhPtPrimMC [mcPOther]->Fill(photonPt) ;
385 fhPhiPrimMC[mcPOther]->Fill(photonE , photonPhi) ;
386 fhYPrimMC[mcPOther] ->Fill(photonE , photonEta) ;
3d5d5078 387 }
388 if(inacceptance){
f66d95af 389 fhEPrimMCAcc[mcPOther] ->Fill(photonE ) ;
390 fhPtPrimMCAcc[mcPOther] ->Fill(photonPt) ;
391 fhPhiPrimMCAcc[mcPOther]->Fill(photonE , photonPhi) ;
392 fhYPrimMCAcc[mcPOther] ->Fill(photonE , photonY) ;
3d5d5078 393 }//Accepted
394 }//Other origin
395 }// Primary photon
396 }//loop on primaries
397 }//stack exists and data is MC
398 }//read stack
399 else if(GetReader()->ReadAODMCParticles()){
400 TClonesArray * mcparticles = GetReader()->GetAODMCParticles(0);
401 if(mcparticles){
402 Int_t nprim = mcparticles->GetEntriesFast();
403
404 for(Int_t i=0; i < nprim; i++)
405 {
406 AliAODMCParticle * prim = (AliAODMCParticle *) mcparticles->At(i);
407
408 Int_t pdg = prim->GetPdgCode();
409
410 if(pdg == 22){
411
412 // Get tag of this particle photon from fragmentation, decay, prompt ...
413 Int_t tag = GetMCAnalysisUtils()->CheckOrigin(i,GetReader(), 0);
414 if(!GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPhoton)){
415 //A conversion photon from a hadron, skip this kind of photon
416// printf("AliAnaPhoton::FillAcceptanceHistograms() - not a photon, weird!: tag %d, conv %d, pi0 %d, hadron %d, electron %d, unk %d, muon %d,pion %d, proton %d, neutron %d, kaon %d, antiproton %d, antineutron %d\n",tag,
417// GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCConversion),
418// GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0),
419// GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCOther),
420// GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCElectron),
421// GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCUnknown),
422// GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCMuon),
423// GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPion),
424// GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCProton),
425// GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCAntiNeutron),
426// GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCKaon),
427// GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCAntiProton),
428// GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCAntiNeutron));
429
430 return;
431 }
432
433 //Get photon kinematics
434 if(prim->E() == TMath::Abs(prim->Pz())) continue ; //Protection against floating point exception
435
436 Double_t photonY = 0.5*TMath::Log((prim->E()-prim->Pz())/(prim->E()+prim->Pz())) ;
437 Double_t photonE = prim->E() ;
438 Double_t photonPt = prim->Pt() ;
439 Double_t photonPhi = TMath::RadToDeg()*prim->Phi() ;
440 if(photonPhi < 0) photonPhi+=TMath::TwoPi();
441 Double_t photonEta = prim->Eta() ;
442
443 //Check if photons hit the Calorimeter
444 TLorentzVector lv;
445 lv.SetPxPyPzE(prim->Px(),prim->Py(),prim->Pz(),prim->E());
446 Bool_t inacceptance = kFALSE;
447 if (fCalorimeter == "PHOS"){
448 if(GetPHOSGeometry() && GetCaloUtils()->IsPHOSGeoMatrixSet()){
449 Int_t mod ;
450 Double_t x,z ;
451 Double_t vtx[]={prim->Xv(),prim->Yv(),prim->Zv()};
452 if(GetPHOSGeometry()->ImpactOnEmc(vtx, prim->Theta(),prim->Phi(),mod,z,x))
453 inacceptance = kTRUE;
454 if(GetDebug() > 2) printf("In %s Real acceptance? %d\n",fCalorimeter.Data(),inacceptance);
455 }
456 else{
457 if(GetFiducialCut()->IsInFiducialCut(lv,fCalorimeter))
458 inacceptance = kTRUE ;
459 if(GetDebug() > 2) printf("In %s fiducial cut acceptance? %d\n",fCalorimeter.Data(),inacceptance);
460 }
461 }
462 else if(fCalorimeter == "EMCAL" && GetCaloUtils()->IsEMCALGeoMatrixSet()){
463 if(GetEMCALGeometry()){
464
465 Int_t absID=0;
466
467 GetEMCALGeometry()->GetAbsCellIdFromEtaPhi(prim->Eta(),prim->Phi(),absID);
468
469 if( absID >= 0)
470 inacceptance = kTRUE;
471
472 if(GetDebug() > 2) printf("In %s Real acceptance? %d\n",fCalorimeter.Data(),inacceptance);
473 }
474 else{
475 if(GetFiducialCut()->IsInFiducialCut(lv,fCalorimeter))
476 inacceptance = kTRUE ;
477 if(GetDebug() > 2) printf("In %s fiducial cut acceptance? %d\n",fCalorimeter.Data(),inacceptance);
478 }
479 } //In EMCAL
480
481 //Fill histograms
482
f66d95af 483 fhYPrimMC[mcPPhoton]->Fill(photonPt, photonY) ;
3d5d5078 484 if(TMath::Abs(photonY) < 1.0)
485 {
f66d95af 486 fhEPrimMC [mcPPhoton]->Fill(photonE ) ;
487 fhPtPrimMC [mcPPhoton]->Fill(photonPt) ;
488 fhPhiPrimMC[mcPPhoton]->Fill(photonE , photonPhi) ;
489 fhYPrimMC[mcPPhoton]->Fill(photonE , photonEta) ;
3d5d5078 490 }
491 if(inacceptance){
f66d95af 492 fhEPrimMCAcc[mcPPhoton] ->Fill(photonE ) ;
493 fhPtPrimMCAcc[mcPPhoton] ->Fill(photonPt) ;
494 fhPhiPrimMCAcc[mcPPhoton]->Fill(photonE , photonPhi) ;
495 fhYPrimMCAcc[mcPPhoton] ->Fill(photonE , photonY) ;
3d5d5078 496 }//Accepted
497
498
499 //Origin of photon
f66d95af 500 if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPrompt) && fhEPrimMC[mcPPrompt])
3d5d5078 501 {
f66d95af 502 fhYPrimMC[mcPPrompt]->Fill(photonPt, photonY) ;
3d5d5078 503 if(TMath::Abs(photonY) < 1.0){
f66d95af 504 fhEPrimMC [mcPPrompt]->Fill(photonE ) ;
505 fhPtPrimMC [mcPPrompt]->Fill(photonPt) ;
506 fhPhiPrimMC[mcPPrompt]->Fill(photonE , photonPhi) ;
507 fhYPrimMC[mcPPrompt]->Fill(photonE , photonEta) ;
3d5d5078 508 }
509 if(inacceptance){
f66d95af 510 fhEPrimMCAcc[mcPPrompt] ->Fill(photonE ) ;
511 fhPtPrimMCAcc[mcPPrompt] ->Fill(photonPt) ;
512 fhPhiPrimMCAcc[mcPPrompt]->Fill(photonE , photonPhi) ;
513 fhYPrimMCAcc[mcPPrompt] ->Fill(photonE , photonY) ;
3d5d5078 514 }//Accepted
515 }
f66d95af 516 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCFragmentation) && fhEPrimMC[mcPFragmentation] )
3d5d5078 517 {
f66d95af 518 fhYPrimMC[mcPFragmentation]->Fill(photonPt, photonY) ;
3d5d5078 519 if(TMath::Abs(photonY) < 1.0){
f66d95af 520 fhEPrimMC [mcPFragmentation]->Fill(photonE ) ;
521 fhPtPrimMC [mcPFragmentation]->Fill(photonPt) ;
522 fhPhiPrimMC[mcPFragmentation]->Fill(photonE , photonPhi) ;
523 fhYPrimMC[mcPFragmentation]->Fill(photonE , photonEta) ;
3d5d5078 524 }
525 if(inacceptance){
f66d95af 526 fhEPrimMCAcc[mcPFragmentation] ->Fill(photonE ) ;
527 fhPtPrimMCAcc[mcPFragmentation] ->Fill(photonPt) ;
528 fhPhiPrimMCAcc[mcPFragmentation]->Fill(photonE , photonPhi) ;
529 fhYPrimMCAcc[mcPFragmentation] ->Fill(photonE , photonY) ;
3d5d5078 530 }//Accepted
531 }
f66d95af 532 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCISR) && fhEPrimMC[mcPISR])
3d5d5078 533 {
f66d95af 534 fhYPrimMC[mcPISR]->Fill(photonPt, photonY) ;
3d5d5078 535 if(TMath::Abs(photonY) < 1.0){
f66d95af 536 fhEPrimMC [mcPISR]->Fill(photonE ) ;
537 fhPtPrimMC [mcPISR]->Fill(photonPt) ;
538 fhPhiPrimMC[mcPISR]->Fill(photonE , photonPhi) ;
539 fhYPrimMC[mcPISR]->Fill(photonE , photonEta) ;
3d5d5078 540 }
541 if(inacceptance){
f66d95af 542 fhEPrimMCAcc[mcPISR] ->Fill(photonE ) ;
543 fhPtPrimMCAcc[mcPISR] ->Fill(photonPt) ;
544 fhPhiPrimMCAcc[mcPISR]->Fill(photonE , photonPhi) ;
545 fhYPrimMCAcc[mcPISR] ->Fill(photonE , photonY) ;
3d5d5078 546 }//Accepted
547 }
f66d95af 548 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0Decay)&& fhEPrimMC[mcPPi0Decay])
3d5d5078 549 {
f66d95af 550 fhYPrimMC[mcPPi0Decay]->Fill(photonPt, photonY) ;
3d5d5078 551 if(TMath::Abs(photonY) < 1.0){
f66d95af 552 fhEPrimMC [mcPPi0Decay]->Fill(photonE ) ;
553 fhPtPrimMC [mcPPi0Decay]->Fill(photonPt) ;
554 fhPhiPrimMC[mcPPi0Decay]->Fill(photonE , photonPhi) ;
555 fhYPrimMC[mcPPi0Decay]->Fill(photonE , photonEta) ;
3d5d5078 556 }
557 if(inacceptance){
f66d95af 558 fhEPrimMCAcc[mcPPi0Decay] ->Fill(photonE ) ;
559 fhPtPrimMCAcc[mcPPi0Decay] ->Fill(photonPt) ;
560 fhPhiPrimMCAcc[mcPPi0Decay]->Fill(photonE , photonPhi) ;
561 fhYPrimMCAcc[mcPPi0Decay] ->Fill(photonE , photonY) ;
3d5d5078 562 }//Accepted
563 }
f586f4aa 564 else if((GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCEtaDecay) ||
565 GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCOtherDecay) ) && fhEPrimMC[mcPOtherDecay])
3d5d5078 566 {
f66d95af 567 fhYPrimMC[mcPOtherDecay]->Fill(photonPt, photonY) ;
3d5d5078 568 if(TMath::Abs(photonY) < 1.0){
f66d95af 569 fhEPrimMC [mcPOtherDecay]->Fill(photonE ) ;
570 fhPtPrimMC [mcPOtherDecay]->Fill(photonPt) ;
571 fhPhiPrimMC[mcPOtherDecay]->Fill(photonE , photonPhi) ;
572 fhYPrimMC[mcPOtherDecay]->Fill(photonE , photonEta) ;
3d5d5078 573 }
574 if(inacceptance){
f66d95af 575 fhEPrimMCAcc[mcPOtherDecay] ->Fill(photonE ) ;
576 fhPtPrimMCAcc[mcPOtherDecay] ->Fill(photonPt) ;
577 fhPhiPrimMCAcc[mcPOtherDecay]->Fill(photonE , photonPhi) ;
578 fhYPrimMCAcc[mcPOtherDecay] ->Fill(photonE , photonY) ;
3d5d5078 579 }//Accepted
580 }
f66d95af 581 else if(fhEPrimMC[mcPOther])
3d5d5078 582 {
f66d95af 583 fhYPrimMC[mcPOther]->Fill(photonPt, photonY) ;
3d5d5078 584 if(TMath::Abs(photonY) < 1.0){
f66d95af 585 fhEPrimMC [mcPOther]->Fill(photonE ) ;
586 fhPtPrimMC [mcPOther]->Fill(photonPt) ;
587 fhPhiPrimMC[mcPOther]->Fill(photonE , photonPhi) ;
588 fhYPrimMC[mcPOther]->Fill(photonE , photonEta) ;
3d5d5078 589 }
590 if(inacceptance){
f66d95af 591 fhEPrimMCAcc[mcPOther] ->Fill(photonE ) ;
592 fhPtPrimMCAcc[mcPOther] ->Fill(photonPt) ;
593 fhPhiPrimMCAcc[mcPOther]->Fill(photonE , photonPhi) ;
594 fhYPrimMCAcc[mcPOther] ->Fill(photonE , photonY) ;
3d5d5078 595 }//Accepted
596 }//Other origin
597 }// Primary photon
598 }//loop on primaries
599
600 }//mc array exists and data is MC
601 } // read AOD MC
602}
521636d2 603
604//__________________________________________________________________
605void AliAnaPhoton::FillShowerShapeHistograms(AliVCluster* cluster, const Int_t mcTag){
606
607 //Fill cluster Shower Shape histograms
608
609 if(!fFillSSHistograms || GetMixedEvent()) return;
610
611 Float_t energy = cluster->E();
612 Int_t ncells = cluster->GetNCells();
521636d2 613 Float_t lambda0 = cluster->GetM02();
614 Float_t lambda1 = cluster->GetM20();
615 Float_t disp = cluster->GetDispersion()*cluster->GetDispersion();
616
617 TLorentzVector mom;
618 if(GetReader()->GetDataType() != AliCaloTrackReader::kMC){
619 cluster->GetMomentum(mom,GetVertex(0)) ;}//Assume that come from vertex in straight line
620 else{
621 Double_t vertex[]={0,0,0};
622 cluster->GetMomentum(mom,vertex) ;
623 }
624
625 Float_t eta = mom.Eta();
626 Float_t phi = mom.Phi();
627 if(phi < 0) phi+=TMath::TwoPi();
628
629 fhLam0E ->Fill(energy,lambda0);
630 fhLam1E ->Fill(energy,lambda1);
631 fhDispE ->Fill(energy,disp);
521636d2 632
633 if(fCalorimeter == "EMCAL" && GetModuleNumber(cluster) > 5){
634 fhLam0ETRD->Fill(energy,lambda0);
635 fhLam1ETRD->Fill(energy,lambda1);
636 fhDispETRD->Fill(energy,disp);
521636d2 637 }
638
639 if(energy < 2){
640 fhNCellsLam0LowE ->Fill(ncells,lambda0);
641 fhNCellsLam1LowE ->Fill(ncells,lambda1);
642 fhNCellsDispLowE ->Fill(ncells,disp);
521636d2 643
644 fhLam1Lam0LowE ->Fill(lambda1,lambda0);
645 fhLam0DispLowE ->Fill(lambda0,disp);
646 fhDispLam1LowE ->Fill(disp,lambda1);
647 fhEtaLam0LowE ->Fill(eta,lambda0);
648 fhPhiLam0LowE ->Fill(phi,lambda0);
649
521636d2 650 }
651 else {
652 fhNCellsLam0HighE ->Fill(ncells,lambda0);
653 fhNCellsLam1HighE ->Fill(ncells,lambda1);
654 fhNCellsDispHighE ->Fill(ncells,disp);
7c65ad18 655
521636d2 656 fhLam1Lam0HighE ->Fill(lambda1,lambda0);
657 fhLam0DispHighE ->Fill(lambda0,disp);
658 fhDispLam1HighE ->Fill(disp,lambda1);
659 fhEtaLam0HighE ->Fill(eta, lambda0);
660 fhPhiLam0HighE ->Fill(phi, lambda0);
521636d2 661 }
662
663 if(IsDataMC()){
3d5d5078 664
f66d95af 665 AliVCaloCells* cells = 0;
666 if(fCalorimeter == "EMCAL") cells = GetEMCALCells();
667 else cells = GetPHOSCells();
3d5d5078 668
669 //Fill histograms to check shape of embedded clusters
670 Float_t fraction = 0;
f66d95af 671 if(GetReader()->IsEmbeddedClusterSelectionOn()){//Only working for EMCAL
672
3d5d5078 673 Float_t clusterE = 0; // recalculate in case corrections applied.
674 Float_t cellE = 0;
675 for(Int_t icell = 0; icell < cluster->GetNCells(); icell++){
676 cellE = cells->GetCellAmplitude(cluster->GetCellAbsId(icell));
677 clusterE+=cellE;
678 fraction+=cellE*cluster->GetCellAmplitudeFraction(icell);
679 }
680
681 //Fraction of total energy due to the embedded signal
682 fraction/=clusterE;
683
684 if(GetDebug() > 1 ) printf("AliAnaPhoton::FillShowerShapeHistogram() - Energy fraction of embedded signal %2.3f, Energy %2.3f\n",fraction, clusterE);
685
686 fhEmbeddedSignalFractionEnergy->Fill(clusterE,fraction);
687
688 } // embedded fraction
689
f66d95af 690 // Get the fraction of the cluster energy that carries the cell with highest energy
691 Int_t absID =-1 ;
692 Float_t maxCellFraction = 0.;
693
694 absID = GetCaloUtils()->GetMaxEnergyCell(cells, cluster,maxCellFraction);
695
696 // Check the origin and fill histograms
521636d2 697 if( GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPhoton) &&
3d5d5078 698 !GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCConversion) &&
699 !GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPi0) &&
700 !GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCEta)){
701 fhMCELambda0[mcssPhoton] ->Fill(energy, lambda0);
3d5d5078 702 fhMCELambda1[mcssPhoton] ->Fill(energy, lambda1);
3d5d5078 703 fhMCEDispersion[mcssPhoton] ->Fill(energy, disp);
f66d95af 704 fhMCNCellsE[mcssPhoton] ->Fill(energy, ncells);
705 fhMCMaxCellDiffClusterE[mcssPhoton]->Fill(energy,maxCellFraction);
3d5d5078 706
f66d95af 707 if (energy < 2.){
708 fhMCLambda0vsClusterMaxCellDiffE0[mcssPhoton]->Fill(lambda0, maxCellFraction);
709 fhMCNCellsvsClusterMaxCellDiffE0[mcssPhoton] ->Fill(ncells, maxCellFraction);
710 }
711 else if(energy < 6.){
712 fhMCLambda0vsClusterMaxCellDiffE2[mcssPhoton]->Fill(lambda0, maxCellFraction);
713 fhMCNCellsvsClusterMaxCellDiffE2[mcssPhoton] ->Fill(ncells, maxCellFraction);
714 }
715 else{
716 fhMCLambda0vsClusterMaxCellDiffE6[mcssPhoton]->Fill(lambda0, maxCellFraction);
717 fhMCNCellsvsClusterMaxCellDiffE6[mcssPhoton] ->Fill(ncells, maxCellFraction);
718 }
3d5d5078 719
720 if(!GetReader()->IsEmbeddedClusterSelectionOn()){
721 //Check particle overlaps in cluster
722
723 //Compare the primary depositing more energy with the rest, if no photon/electron as comon ancestor (conversions), count as other particle
724 Int_t ancPDG = 0, ancStatus = -1;
725 TLorentzVector momentum; TVector3 prodVertex;
726 Int_t ancLabel = 0;
727 Int_t noverlaps = 1;
728 for (UInt_t ilab = 0; ilab < cluster->GetNLabels(); ilab++ ) {
729 ancLabel = GetMCAnalysisUtils()->CheckCommonAncestor(cluster->GetLabels()[0],cluster->GetLabels()[ilab], GetReader(),ancPDG,ancStatus,momentum,prodVertex);
730 if(ancPDG!=22 && TMath::Abs(ancPDG)!=11) noverlaps++;
731 }
732
733 if(noverlaps == 1){
734 fhMCPhotonELambda0NoOverlap ->Fill(energy, lambda0);
3d5d5078 735 }
736 else if(noverlaps == 2){
737 fhMCPhotonELambda0TwoOverlap ->Fill(energy, lambda0);
3d5d5078 738 }
739 else if(noverlaps > 2){
740 fhMCPhotonELambda0NOverlap ->Fill(energy, lambda0);
3d5d5078 741 }
742 else {
743 printf("AliAnaPhoton::FillShowerShapeHistogram() - n overlaps = %d!!", noverlaps);
744 }
745 }//No embedding
746
747 //Fill histograms to check shape of embedded clusters
748 if(GetReader()->IsEmbeddedClusterSelectionOn()){
749
750 if (fraction > 0.9)
751 {
752 fhEmbedPhotonELambda0FullSignal ->Fill(energy, lambda0);
3d5d5078 753 }
754 else if(fraction > 0.5)
755 {
756 fhEmbedPhotonELambda0MostlySignal ->Fill(energy, lambda0);
3d5d5078 757 }
758 else if(fraction > 0.1)
759 {
760 fhEmbedPhotonELambda0MostlyBkg ->Fill(energy, lambda0);
3d5d5078 761 }
762 else
763 {
764 fhEmbedPhotonELambda0FullBkg ->Fill(energy, lambda0);
3d5d5078 765 }
766 } // embedded
767
521636d2 768 }//photon no conversion
769 else if ( GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCElectron)){
3d5d5078 770 fhMCELambda0[mcssElectron] ->Fill(energy, lambda0);
3d5d5078 771 fhMCELambda1[mcssElectron] ->Fill(energy, lambda1);
3d5d5078 772 fhMCEDispersion[mcssElectron] ->Fill(energy, disp);
f66d95af 773 fhMCNCellsE[mcssElectron] ->Fill(energy, ncells);
774 fhMCMaxCellDiffClusterE[mcssElectron]->Fill(energy,maxCellFraction);
775
776 if (energy < 2.){
777 fhMCLambda0vsClusterMaxCellDiffE0[mcssElectron]->Fill(lambda0, maxCellFraction);
778 fhMCNCellsvsClusterMaxCellDiffE0[mcssElectron] ->Fill(ncells, maxCellFraction);
779 }
780 else if(energy < 6.){
781 fhMCLambda0vsClusterMaxCellDiffE2[mcssElectron]->Fill(lambda0, maxCellFraction);
782 fhMCNCellsvsClusterMaxCellDiffE2[mcssElectron] ->Fill(ncells, maxCellFraction);
783 }
784 else{
785 fhMCLambda0vsClusterMaxCellDiffE6[mcssElectron]->Fill(lambda0, maxCellFraction);
786 fhMCNCellsvsClusterMaxCellDiffE6[mcssElectron] ->Fill(ncells, maxCellFraction);
787 }
521636d2 788 }//electron
3d5d5078 789 else if ( GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPhoton) &&
521636d2 790 GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCConversion) ){
3d5d5078 791 fhMCELambda0[mcssConversion] ->Fill(energy, lambda0);
3d5d5078 792 fhMCELambda1[mcssConversion] ->Fill(energy, lambda1);
3d5d5078 793 fhMCEDispersion[mcssConversion] ->Fill(energy, disp);
f66d95af 794 fhMCNCellsE[mcssConversion] ->Fill(energy, ncells);
795 fhMCMaxCellDiffClusterE[mcssConversion]->Fill(energy,maxCellFraction);
796
797 if (energy < 2.){
798 fhMCLambda0vsClusterMaxCellDiffE0[mcssConversion]->Fill(lambda0, maxCellFraction);
799 fhMCNCellsvsClusterMaxCellDiffE0[mcssConversion] ->Fill(ncells, maxCellFraction);
800 }
801 else if(energy < 6.){
802 fhMCLambda0vsClusterMaxCellDiffE2[mcssConversion]->Fill(lambda0, maxCellFraction);
803 fhMCNCellsvsClusterMaxCellDiffE2[mcssConversion] ->Fill(ncells, maxCellFraction);
804 }
805 else{
806 fhMCLambda0vsClusterMaxCellDiffE6[mcssConversion]->Fill(lambda0, maxCellFraction);
807 fhMCNCellsvsClusterMaxCellDiffE6[mcssConversion] ->Fill(ncells, maxCellFraction);
808 }
3d5d5078 809
521636d2 810 }//conversion photon
811 else if ( GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPi0) ){
3d5d5078 812 fhMCELambda0[mcssPi0] ->Fill(energy, lambda0);
3d5d5078 813 fhMCELambda1[mcssPi0] ->Fill(energy, lambda1);
3d5d5078 814 fhMCEDispersion[mcssPi0] ->Fill(energy, disp);
f66d95af 815 fhMCNCellsE[mcssPi0] ->Fill(energy, ncells);
816 fhMCMaxCellDiffClusterE[mcssPi0]->Fill(energy,maxCellFraction);
817
818 if (energy < 2.){
819 fhMCLambda0vsClusterMaxCellDiffE0[mcssPi0]->Fill(lambda0, maxCellFraction);
820 fhMCNCellsvsClusterMaxCellDiffE0[mcssPi0] ->Fill(ncells, maxCellFraction);
821 }
822 else if(energy < 6.){
823 fhMCLambda0vsClusterMaxCellDiffE2[mcssPi0]->Fill(lambda0, maxCellFraction);
824 fhMCNCellsvsClusterMaxCellDiffE2[mcssPi0] ->Fill(ncells, maxCellFraction);
825 }
826 else{
827 fhMCLambda0vsClusterMaxCellDiffE6[mcssPi0]->Fill(lambda0, maxCellFraction);
828 fhMCNCellsvsClusterMaxCellDiffE6[mcssPi0] ->Fill(ncells, maxCellFraction);
829 }
3d5d5078 830
831 //Fill histograms to check shape of embedded clusters
832 if(GetReader()->IsEmbeddedClusterSelectionOn()){
833
834 if (fraction > 0.9)
835 {
836 fhEmbedPi0ELambda0FullSignal ->Fill(energy, lambda0);
3d5d5078 837 }
838 else if(fraction > 0.5)
839 {
840 fhEmbedPi0ELambda0MostlySignal ->Fill(energy, lambda0);
3d5d5078 841 }
842 else if(fraction > 0.1)
843 {
844 fhEmbedPi0ELambda0MostlyBkg ->Fill(energy, lambda0);
3d5d5078 845 }
846 else
847 {
848 fhEmbedPi0ELambda0FullBkg ->Fill(energy, lambda0);
3d5d5078 849 }
850 } // embedded
851
521636d2 852 }//pi0
3d5d5078 853 else if ( GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCEta) ){
854 fhMCELambda0[mcssEta] ->Fill(energy, lambda0);
3d5d5078 855 fhMCELambda1[mcssEta] ->Fill(energy, lambda1);
3d5d5078 856 fhMCEDispersion[mcssEta] ->Fill(energy, disp);
f66d95af 857 fhMCNCellsE[mcssEta] ->Fill(energy, ncells);
858 fhMCMaxCellDiffClusterE[mcssEta]->Fill(energy,maxCellFraction);
859
860 if (energy < 2.){
861 fhMCLambda0vsClusterMaxCellDiffE0[mcssEta]->Fill(lambda0, maxCellFraction);
862 fhMCNCellsvsClusterMaxCellDiffE0[mcssEta] ->Fill(ncells, maxCellFraction);
863 }
864 else if(energy < 6.){
865 fhMCLambda0vsClusterMaxCellDiffE2[mcssEta]->Fill(lambda0, maxCellFraction);
866 fhMCNCellsvsClusterMaxCellDiffE2[mcssEta] ->Fill(ncells, maxCellFraction);
867 }
868 else{
869 fhMCLambda0vsClusterMaxCellDiffE6[mcssEta]->Fill(lambda0, maxCellFraction);
870 fhMCNCellsvsClusterMaxCellDiffE6[mcssEta] ->Fill(ncells, maxCellFraction);
871 }
872
3d5d5078 873 }//eta
521636d2 874 else {
3d5d5078 875 fhMCELambda0[mcssOther] ->Fill(energy, lambda0);
3d5d5078 876 fhMCELambda1[mcssOther] ->Fill(energy, lambda1);
3d5d5078 877 fhMCEDispersion[mcssOther] ->Fill(energy, disp);
f66d95af 878 fhMCNCellsE[mcssOther] ->Fill(energy, ncells);
879 fhMCMaxCellDiffClusterE[mcssOther]->Fill(energy,maxCellFraction);
880
881 if (energy < 2.){
882 fhMCLambda0vsClusterMaxCellDiffE0[mcssOther]->Fill(lambda0, maxCellFraction);
883 fhMCNCellsvsClusterMaxCellDiffE0[mcssOther] ->Fill(ncells, maxCellFraction);
884 }
885 else if(energy < 6.){
886 fhMCLambda0vsClusterMaxCellDiffE2[mcssOther]->Fill(lambda0, maxCellFraction);
887 fhMCNCellsvsClusterMaxCellDiffE2[mcssOther] ->Fill(ncells, maxCellFraction);
888 }
889 else{
890 fhMCLambda0vsClusterMaxCellDiffE6[mcssOther]->Fill(lambda0, maxCellFraction);
891 fhMCNCellsvsClusterMaxCellDiffE6[mcssOther] ->Fill(ncells, maxCellFraction);
892 }
893
521636d2 894 }//other particles
895
896 }//MC data
897
898}
899
0c1383b5 900//________________________________________________________________________
901TObjString * AliAnaPhoton::GetAnalysisCuts()
902{
903 //Save parameters used for analysis
904 TString parList ; //this will be list of parameters used for this analysis.
5ae09196 905 const Int_t buffersize = 255;
906 char onePar[buffersize] ;
0c1383b5 907
5ae09196 908 snprintf(onePar,buffersize,"--- AliAnaPhoton ---\n") ;
0c1383b5 909 parList+=onePar ;
5ae09196 910 snprintf(onePar,buffersize,"Calorimeter: %s\n",fCalorimeter.Data()) ;
0c1383b5 911 parList+=onePar ;
5ae09196 912 snprintf(onePar,buffersize,"fMinDist =%2.2f (Minimal distance to bad channel to accept cluster) \n",fMinDist) ;
0c1383b5 913 parList+=onePar ;
5ae09196 914 snprintf(onePar,buffersize,"fMinDist2=%2.2f (Cuts on Minimal distance to study acceptance evaluation) \n",fMinDist2) ;
0c1383b5 915 parList+=onePar ;
5ae09196 916 snprintf(onePar,buffersize,"fMinDist3=%2.2f (One more cut on distance used for acceptance-efficiency study) \n",fMinDist3) ;
0c1383b5 917 parList+=onePar ;
5ae09196 918 snprintf(onePar,buffersize,"fRejectTrackMatch: %d\n",fRejectTrackMatch) ;
0c1383b5 919 parList+=onePar ;
920
921 //Get parameters set in base class.
922 parList += GetBaseParametersList() ;
923
924 //Get parameters set in PID class.
925 parList += GetCaloPID()->GetPIDParametersList() ;
926
927 //Get parameters set in FiducialCut class (not available yet)
928 //parlist += GetFidCut()->GetFidCutParametersList()
929
930 return new TObjString(parList) ;
931}
932
1c5acb87 933//________________________________________________________________________
934TList * AliAnaPhoton::GetCreateOutputObjects()
935{
477d6cee 936 // Create histograms to be saved in output file and
937 // store them in outputContainer
938 TList * outputContainer = new TList() ;
939 outputContainer->SetName("PhotonHistos") ;
4a745797 940
521636d2 941 Int_t nptbins = GetHistoPtBins(); Float_t ptmax = GetHistoPtMax(); Float_t ptmin = GetHistoPtMin();
942 Int_t nphibins = GetHistoPhiBins(); Float_t phimax = GetHistoPhiMax(); Float_t phimin = GetHistoPhiMin();
943 Int_t netabins = GetHistoEtaBins(); Float_t etamax = GetHistoEtaMax(); Float_t etamin = GetHistoEtaMin();
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
6639984f 1446//____________________________________________________________________________
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
1461}
1462
1c5acb87 1463//____________________________________________________________________________
1464void AliAnaPhoton::InitParameters()
1465{
1466
1467 //Initialize the parameters of the analysis.
a3aebfff 1468 AddToHistogramsName("AnaPhoton_");
521636d2 1469
6175da48 1470 fCalorimeter = "EMCAL" ;
1471 fMinDist = 2.;
1472 fMinDist2 = 4.;
1473 fMinDist3 = 5.;
1e86c71e 1474
4cf55759 1475 fTimeCutMin = -1;
1476 fTimeCutMax = 9999999;
6175da48 1477 fNCellsCut = 0;
2ac125bf 1478
1e86c71e 1479 fRejectTrackMatch = kTRUE ;
1e86c71e 1480
1c5acb87 1481}
1482
1483//__________________________________________________________________
1484void AliAnaPhoton::MakeAnalysisFillAOD()
1485{
f8006433 1486 //Do photon analysis and fill aods
f37fa8d2 1487
6175da48 1488 //Get the vertex
5025c139 1489 Double_t v[3] = {0,0,0}; //vertex ;
1490 GetReader()->GetVertex(v);
f8006433 1491
f37fa8d2 1492 //Select the Calorimeter of the photon
c8fe2783 1493 TObjArray * pl = 0x0;
477d6cee 1494 if(fCalorimeter == "PHOS")
be518ab0 1495 pl = GetPHOSClusters();
477d6cee 1496 else if (fCalorimeter == "EMCAL")
be518ab0 1497 pl = GetEMCALClusters();
5ae09196 1498
1499 if(!pl) {
1500 Info("MakeAnalysisFillAOD","TObjArray with %s clusters is NULL!\n",fCalorimeter.Data());
1501 return;
1502 }
521636d2 1503
6175da48 1504 //Init arrays, variables, get number of clusters
1e86c71e 1505 TLorentzVector mom, mom2 ;
1506 Int_t nCaloClusters = pl->GetEntriesFast();
20218aea 1507
6175da48 1508 if(GetDebug() > 0) printf("AliAnaPhoton::MakeAnalysisFillAOD() - input %s cluster entries %d\n", fCalorimeter.Data(), nCaloClusters);
521636d2 1509
6175da48 1510 //----------------------------------------------------
1511 // Fill AOD with PHOS/EMCAL AliAODPWG4Particle objects
1512 //----------------------------------------------------
1513 // Loop on clusters
1e86c71e 1514 for(Int_t icalo = 0; icalo < nCaloClusters; icalo++){
1515
0ae57829 1516 AliVCluster * calo = (AliVCluster*) (pl->At(icalo));
1517 //printf("calo %d, %f\n",icalo,calo->E());
521636d2 1518
f8006433 1519 //Get the index where the cluster comes, to retrieve the corresponding vertex
c8fe2783 1520 Int_t evtIndex = 0 ;
1521 if (GetMixedEvent()) {
1522 evtIndex=GetMixedEvent()->EventIndexForCaloCluster(calo->GetID()) ;
5025c139 1523 //Get the vertex and check it is not too large in z
96539743 1524 if(TMath::Abs(GetVertex(evtIndex)[2])> GetZvertexCut()) continue;
c8fe2783 1525 }
521636d2 1526
1527 //Cluster selection, not charged, with photon id and in fiducial cut
f8006433 1528 if(GetReader()->GetDataType() != AliCaloTrackReader::kMC){
1529 calo->GetMomentum(mom,GetVertex(evtIndex)) ;}//Assume that come from vertex in straight line
1530 else{
1531 Double_t vertex[]={0,0,0};
1532 calo->GetMomentum(mom,vertex) ;
1533 }
c8fe2783 1534
6175da48 1535 //--------------------------------------
1536 // Cluster selection
1537 //--------------------------------------
c4a7d28a 1538 if(!ClusterSelected(calo,mom)) continue;
6175da48 1539
1540 //----------------------------
1541 //Create AOD for analysis
1542 //----------------------------
1543 AliAODPWG4Particle aodph = AliAODPWG4Particle(mom);
1544
1545 //...............................................
1546 //Set the indeces of the original caloclusters (MC, ID), and calorimeter
1547 Int_t label = calo->GetLabel();
1548 aodph.SetLabel(label);
6175da48 1549 aodph.SetCaloLabel(calo->GetID(),-1);
1550 aodph.SetDetector(fCalorimeter);
c4a7d28a 1551 //printf("Index %d, Id %d, iaod %d\n",icalo, calo->GetID(),GetOutputAODBranch()->GetEntriesFast());
521636d2 1552
6175da48 1553 //...............................................
1554 //Set bad channel distance bit
c4a7d28a 1555 Double_t distBad=calo->GetDistanceToBadChannel() ; //Distance to bad channel
f37fa8d2 1556 if (distBad > fMinDist3) aodph.SetDistToBad(2) ;
477d6cee 1557 else if(distBad > fMinDist2) aodph.SetDistToBad(1) ;
f37fa8d2 1558 else aodph.SetDistToBad(0) ;
af7b3903 1559 //printf("DistBad %f Bit %d\n",distBad, aodph.DistToBad());
c8fe2783 1560
521636d2 1561 //--------------------------------------------------------------------------------------
1562 //Play with the MC stack if available
1563 //--------------------------------------------------------------------------------------
1564
1565 //Check origin of the candidates
1566 if(IsDataMC()){
1567 aodph.SetTag(GetMCAnalysisUtils()->CheckOrigin(calo->GetLabels(),calo->GetNLabels(),GetReader(), aodph.GetInputFileIndex()));
3d5d5078 1568
521636d2 1569 if(GetDebug() > 0)
1570 printf("AliAnaPhoton::MakeAnalysisFillAOD() - Origin of candidate, bit map %d\n",aodph.GetTag());
1571 }//Work with stack also
1572
1573 //--------------------------------------------------------------------------------------
1574 //Fill some shower shape histograms before PID is applied
1575 //--------------------------------------------------------------------------------------
1576
1577 FillShowerShapeHistograms(calo,aodph.GetTag());
6175da48 1578
1579 //-------------------------------------
f37fa8d2 1580 //PID selection or bit setting
6175da48 1581 //-------------------------------------
1582 // MC
477d6cee 1583 if(GetReader()->GetDataType() == AliCaloTrackReader::kMC){
f37fa8d2 1584 //Get most probable PID, check PID weights (in MC this option is mandatory)
21a4b1c0 1585 aodph.SetIdentifiedParticleType(GetCaloPID()->GetIdentifiedParticleType(fCalorimeter,calo->GetPID(),mom.E()));//PID with weights
1586 if(GetDebug() > 1) printf("AliAnaPhoton::MakeAnalysisFillAOD() - PDG of identified particle %d\n",aodph.GetIdentifiedParticleType());
f37fa8d2 1587 //If primary is not photon, skip it.
21a4b1c0 1588 if(aodph.GetIdentifiedParticleType() != AliCaloPID::kPhoton) continue ;
6175da48 1589 }
1590 //...............................................
1591 // Data, PID check on
477d6cee 1592 else if(IsCaloPIDOn()){
f37fa8d2 1593 //Get most probable PID, 2 options check PID weights
3d5d5078 1594 //or redo PID, recommended option for MCEal.
477d6cee 1595 if(!IsCaloPIDRecalculationOn())
21a4b1c0 1596 aodph.SetIdentifiedParticleType(GetCaloPID()->GetIdentifiedParticleType(fCalorimeter,calo->GetPID(),mom.E()));//PID with weights
477d6cee 1597 else
21a4b1c0 1598 aodph.SetIdentifiedParticleType(GetCaloPID()->GetIdentifiedParticleType(fCalorimeter,mom,calo));//PID recalculated
477d6cee 1599
21a4b1c0 1600 if(GetDebug() > 1) printf("AliAnaPhoton::MakeAnalysisFillAOD() - PDG of identified particle %d\n",aodph.GetIdentifiedParticleType());
477d6cee 1601
f37fa8d2 1602 //If cluster does not pass pid, not photon, skip it.
21a4b1c0 1603 if(aodph.GetIdentifiedParticleType() != AliCaloPID::kPhoton) continue ;
477d6cee 1604
1605 }
6175da48 1606 //...............................................
1607 // Data, PID check off
477d6cee 1608 else{
f37fa8d2 1609 //Set PID bits for later selection (AliAnaPi0 for example)
1610 //GetPDG already called in SetPIDBits.
f2ccb5b8 1611 GetCaloPID()->SetPIDBits(fCalorimeter,calo,&aodph, GetCaloUtils());
a3aebfff 1612 if(GetDebug() > 1) printf("AliAnaPhoton::MakeAnalysisFillAOD() - PID Bits set \n");
477d6cee 1613 }
1614
21a4b1c0 1615 if(GetDebug() > 1) printf("AliAnaPhoton::MakeAnalysisFillAOD() - Photon selection cuts passed: pT %3.2f, pdg %d\n",aodph.Pt(), aodph.GetIdentifiedParticleType());
c8fe2783 1616
f37fa8d2 1617 //Add AOD with photon object to aod branch
477d6cee 1618 AddAODParticle(aodph);
1619
1620 }//loop
5812a064 1621
f37fa8d2 1622 if(GetDebug() > 1) printf("AliAnaPhoton::MakeAnalysisFillAOD() End fill AODs, with %d entries \n",GetOutputAODBranch()->GetEntriesFast());
477d6cee 1623
1c5acb87 1624}
1625
1626//__________________________________________________________________
1627void AliAnaPhoton::MakeAnalysisFillHistograms()
1628{
6175da48 1629 //Fill histograms
f8006433 1630
6175da48 1631 //-------------------------------------------------------------------
577d9801 1632 // Access MC information in stack if requested, check that it exists.
521636d2 1633 AliStack * stack = 0x0;
1634 TParticle * primary = 0x0;
1635 TClonesArray * mcparticles = 0x0;
1636 AliAODMCParticle * aodprimary = 0x0;
1637
577d9801 1638 if(IsDataMC()){
521636d2 1639
577d9801 1640 if(GetReader()->ReadStack()){
1641 stack = GetMCStack() ;
1642 if(!stack) {
3d5d5078 1643 printf("AliAnaPhoton::MakeAnalysisFillHistograms() - Stack not available, is the MC handler called? STOP\n");
1644 abort();
577d9801 1645 }
f8006433 1646
577d9801 1647 }
1648 else if(GetReader()->ReadAODMCParticles()){
f8006433 1649
577d9801 1650 //Get the list of MC particles
521636d2 1651 mcparticles = GetReader()->GetAODMCParticles(0);
1652 if(!mcparticles && GetDebug() > 0) {
3d5d5078 1653 printf("AliAnaPhoton::MakeAnalysisFillHistograms() - Standard MCParticles not available!\n");
577d9801 1654 }
577d9801 1655 }
1656 }// is data and MC
521636d2 1657
6175da48 1658
1659 // Get vertex
2244659d 1660 Double_t v[3] = {0,0,0}; //vertex ;
1661 GetReader()->GetVertex(v);
6175da48 1662 //fhVertex->Fill(v[0],v[1],v[2]);
1663 if(TMath::Abs(v[2]) > GetZvertexCut()) return ; // done elsewhere for Single Event analysis, but there for mixed event
1664
1665 //----------------------------------
577d9801 1666 //Loop on stored AOD photons
1667 Int_t naod = GetOutputAODBranch()->GetEntriesFast();
577d9801 1668 if(GetDebug() > 0) printf("AliAnaPhoton::MakeAnalysisFillHistograms() - aod branch entries %d\n", naod);
521636d2 1669
577d9801 1670 for(Int_t iaod = 0; iaod < naod ; iaod++){
1671 AliAODPWG4Particle* ph = (AliAODPWG4Particle*) (GetOutputAODBranch()->At(iaod));
1672 Int_t pdg = ph->GetIdentifiedParticleType();
521636d2 1673
577d9801 1674 if(GetDebug() > 3)
1675 printf("AliAnaPhoton::MakeAnalysisFillHistograms() - PDG %d, MC TAG %d, Calorimeter %s\n", ph->GetIdentifiedParticleType(),ph->GetTag(), (ph->GetDetector()).Data()) ;
521636d2 1676
577d9801 1677 //If PID used, fill histos with photons in Calorimeter fCalorimeter
1678 if(IsCaloPIDOn() && pdg != AliCaloPID::kPhoton) continue;
1679 if(ph->GetDetector() != fCalorimeter) continue;
521636d2 1680
577d9801 1681 if(GetDebug() > 2)
1682 printf("AliAnaPhoton::MakeAnalysisFillHistograms() - ID Photon: pt %f, phi %f, eta %f\n", ph->Pt(),ph->Phi(),ph->Eta()) ;
521636d2 1683
6175da48 1684 //................................
577d9801 1685 //Fill photon histograms
1686 Float_t ptcluster = ph->Pt();
1687 Float_t phicluster = ph->Phi();
1688 Float_t etacluster = ph->Eta();
1689 Float_t ecluster = ph->E();
521636d2 1690
20218aea 1691 fhEPhoton ->Fill(ecluster);
577d9801 1692 fhPtPhoton ->Fill(ptcluster);
1693 fhPhiPhoton ->Fill(ptcluster,phicluster);
1694 fhEtaPhoton ->Fill(ptcluster,etacluster);
521636d2 1695 if (ecluster > 0.5) fhEtaPhiPhoton ->Fill(etacluster, phicluster);
20218aea 1696 else if(GetMinPt() < 0.5) fhEtaPhi05Photon->Fill(etacluster, phicluster);
1697
5812a064 1698
1699 //Get original cluster, to recover some information
1700 Int_t absID = 0;
1701 Float_t maxCellFraction = 0;
1702 AliVCaloCells* cells = 0;
1703 TObjArray * clusters = 0;
1704 if(fCalorimeter == "EMCAL"){
1705 cells = GetEMCALCells();
1706 clusters = GetEMCALClusters();
1707 }
1708 else{
1709 cells = GetPHOSCells();
1710 clusters = GetPHOSClusters();
6175da48 1711 }
20218aea 1712
5812a064 1713 Int_t iclus = -1;
1714 AliVCluster *cluster = FindCluster(clusters,ph->GetCaloLabel(0),iclus);
1715 absID = GetCaloUtils()->GetMaxEnergyCell(cells, cluster,maxCellFraction);
1716
f15c25da 1717 // Control histograms
5812a064 1718 fhMaxCellDiffClusterE->Fill(ph->E(),maxCellFraction);
1719 fhNCellsE ->Fill(ph->E(),cluster->GetNCells());
f15c25da 1720 fhTimeE ->Fill(ph->E(),cluster->GetTOF()*1.e9);
5812a064 1721
6175da48 1722 //.......................................
577d9801 1723 //Play with the MC data if available
1724 if(IsDataMC()){
521636d2 1725
3d5d5078 1726 FillAcceptanceHistograms();
1727
4c8f7c2e 1728 //....................................................................
1729 // Access MC information in stack if requested, check that it exists.
1730 Int_t label =ph->GetLabel();
1731 if(label < 0) {
1732 if(GetDebug() > 1) printf("AliAnaPhoton::MakeAnalysisFillHistograms() *** bad label ***: label %d \n", label);
1733 continue;
1734 }
1735
1736 Float_t eprim = 0;
1737 Float_t ptprim = 0;
1738 if(GetReader()->ReadStack()){
1739
1740 if(label >= stack->GetNtrack()) {
1741 if(GetDebug() > 2) printf("AliAnaPhoton::MakeAnalysisFillHistograms() *** large label ***: label %d, n tracks %d \n", label, stack->GetNtrack());
1742 continue ;
1743 }
1744
1745 primary = stack->Particle(label);
1746 if(!primary){
1747 printf("AliAnaPhoton::MakeAnalysisFillHistograms() *** no primary ***: label %d \n", label);
1748 continue;
1749 }
1750 eprim = primary->Energy();
1751 ptprim = primary->Pt();
1752
1753 }
1754 else if(GetReader()->ReadAODMCParticles()){
1755 //Check which is the input
1756 if(ph->GetInputFileIndex() == 0){
1757 if(!mcparticles) continue;
1758 if(label >= mcparticles->GetEntriesFast()) {
1759 if(GetDebug() > 2) printf("AliAnaPhoton::MakeAnalysisFillHistograms() *** large label ***: label %d, n tracks %d \n",
1760 label, mcparticles->GetEntriesFast());
1761 continue ;
1762 }
1763 //Get the particle
1764 aodprimary = (AliAODMCParticle*) mcparticles->At(label);
1765
1766 }
1767
1768 if(!aodprimary){
1769 printf("AliAnaPhoton::MakeAnalysisFillHistograms() *** no primary ***: label %d \n", label);
1770 continue;
1771 }
1772
1773 eprim = aodprimary->E();
1774 ptprim = aodprimary->Pt();
1775
1776 }
1777
577d9801 1778 Int_t tag =ph->GetTag();
521636d2 1779
f66d95af 1780 if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPhoton) && fhMCE[mcPhoton])
3d5d5078 1781 {
1782 fhMCE [mcPhoton] ->Fill(ecluster);
4c8f7c2e 1783 fhMCPt [mcPhoton] ->Fill(ptcluster);
1784 fhMCPhi[mcPhoton] ->Fill(ecluster,phicluster);
1785 fhMCEta[mcPhoton] ->Fill(ecluster,etacluster);
1786
1787 fhMC2E[mcPhoton] ->Fill(ecluster, eprim);
1788 fhMC2Pt[mcPhoton] ->Fill(ptcluster, ptprim);
d9105d92 1789 fhMCDeltaE[mcPhoton] ->Fill(ecluster,eprim-ecluster);
1790 fhMCDeltaPt[mcPhoton]->Fill(ptcluster,ptprim-ptcluster);
3d5d5078 1791
f66d95af 1792 if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCConversion) && fhMCE[mcConversion])
3d5d5078 1793 {
1794 fhMCE [mcConversion] ->Fill(ecluster);
4c8f7c2e 1795 fhMCPt [mcConversion] ->Fill(ptcluster);
1796 fhMCPhi[mcConversion] ->Fill(ecluster,phicluster);
1797 fhMCEta[mcConversion] ->Fill(ecluster,etacluster);
1798
1799 fhMC2E[mcConversion] ->Fill(ecluster, eprim);
1800 fhMC2Pt[mcConversion] ->Fill(ptcluster, ptprim);
d9105d92 1801 fhMCDeltaE[mcConversion] ->Fill(ecluster,eprim-ecluster);
1802 fhMCDeltaPt[mcConversion]->Fill(ptcluster,ptprim-ptcluster);
3d5d5078 1803 }
1804
f66d95af 1805 if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPrompt) && fhMCE[mcPrompt]){
3d5d5078 1806 fhMCE [mcPrompt] ->Fill(ecluster);
4c8f7c2e 1807 fhMCPt [mcPrompt] ->Fill(ptcluster);
1808 fhMCPhi[mcPrompt] ->Fill(ecluster,phicluster);
1809 fhMCEta[mcPrompt] ->Fill(ecluster,etacluster);
1810
1811 fhMC2E[mcPrompt] ->Fill(ecluster, eprim);
1812 fhMC2Pt[mcPrompt] ->Fill(ptcluster, ptprim);
d9105d92 1813 fhMCDeltaE[mcPrompt] ->Fill(ecluster,eprim-ecluster);
1814 fhMCDeltaPt[mcPrompt]->Fill(ptcluster,ptprim-ptcluster);
4c8f7c2e 1815
3d5d5078 1816 }
f66d95af 1817 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCFragmentation)&& fhMCE[mcFragmentation])
3d5d5078 1818 {
1819 fhMCE [mcFragmentation] ->Fill(ecluster);
4c8f7c2e 1820 fhMCPt [mcFragmentation] ->Fill(ptcluster);
1821 fhMCPhi[mcFragmentation] ->Fill(ecluster,phicluster);
1822 fhMCEta[mcFragmentation] ->Fill(ecluster,etacluster);
1823
1824 fhMC2E[mcFragmentation] ->Fill(ecluster, eprim);
1825 fhMC2Pt[mcFragmentation] ->Fill(ptcluster, ptprim);
d9105d92 1826 fhMCDeltaE[mcFragmentation] ->Fill(ecluster,eprim-ecluster);
1827 fhMCDeltaPt[mcFragmentation]->Fill(ptcluster,ptprim-ptcluster);
3d5d5078 1828
1829 }
f66d95af 1830 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCISR)&& fhMCE[mcISR])
3d5d5078 1831 {
1832 fhMCE [mcISR] ->Fill(ecluster);
4c8f7c2e 1833 fhMCPt [mcISR] ->Fill(ptcluster);
1834 fhMCPhi[mcISR] ->Fill(ecluster,phicluster);
1835 fhMCEta[mcISR] ->Fill(ecluster,etacluster);
1836
1837 fhMC2E[mcISR] ->Fill(ecluster, eprim);
1838 fhMC2Pt[mcISR] ->Fill(ptcluster, ptprim);
d9105d92 1839 fhMCDeltaE[mcISR] ->Fill(ecluster,eprim-ecluster);
1840 fhMCDeltaPt[mcISR]->Fill(ptcluster,ptprim-ptcluster);
4c8f7c2e 1841
3d5d5078 1842 }
1843 else if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0Decay) &&
f66d95af 1844 !GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0) && fhMCE[mcPi0Decay])
3d5d5078 1845 {
1846 fhMCE [mcPi0Decay] ->Fill(ecluster);
4c8f7c2e 1847 fhMCPt [mcPi0Decay] ->Fill(ptcluster);
1848 fhMCPhi[mcPi0Decay] ->Fill(ecluster,phicluster);
1849 fhMCEta[mcPi0Decay] ->Fill(ecluster,etacluster);
1850
1851 fhMC2E[mcPi0Decay] ->Fill(ecluster, eprim);
1852 fhMC2Pt[mcPi0Decay] ->Fill(ptcluster, ptprim);
d9105d92 1853 fhMCDeltaE[mcPi0Decay] ->Fill(ecluster,eprim-ecluster);
1854 fhMCDeltaPt[mcPi0Decay]->Fill(ptcluster,ptprim-ptcluster);
4c8f7c2e 1855
3d5d5078 1856 }
f586f4aa 1857 else if( (GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCEtaDecay) ||
1858 GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCOtherDecay) ) && fhMCE[mcOtherDecay])
3d5d5078 1859 {
1860 fhMCE [mcOtherDecay] ->Fill(ecluster);
4c8f7c2e 1861 fhMCPt [mcOtherDecay] ->Fill(ptcluster);
1862 fhMCPhi[mcOtherDecay] ->Fill(ecluster,phicluster);
1863 fhMCEta[mcOtherDecay] ->Fill(ecluster,etacluster);
1864
1865 fhMC2E[mcOtherDecay] ->Fill(ecluster, eprim);
1866 fhMC2Pt[mcOtherDecay] ->Fill(ptcluster, ptprim);
d9105d92 1867 fhMCDeltaE[mcOtherDecay] ->Fill(ecluster,eprim-ecluster);
1868 fhMCDeltaPt[mcOtherDecay]->Fill(ptcluster,ptprim-ptcluster);
4c8f7c2e 1869
3d5d5078 1870 }
f66d95af 1871 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0) && fhMCE [mcPi0])
3d5d5078 1872 {
1873 fhMCE [mcPi0] ->Fill(ecluster);
4c8f7c2e 1874 fhMCPt [mcPi0] ->Fill(ptcluster);
1875 fhMCPhi[mcPi0] ->Fill(ecluster,phicluster);
1876 fhMCEta[mcPi0] ->Fill(ecluster,etacluster);
1877
1878 fhMC2E[mcPi0] ->Fill(ecluster, eprim);
1879 fhMC2Pt[mcPi0] ->Fill(ptcluster, ptprim);
d9105d92 1880 fhMCDeltaE[mcPi0] ->Fill(ecluster,eprim-ecluster);
1881 fhMCDeltaPt[mcPi0]->Fill(ptcluster,ptprim-ptcluster);
4c8f7c2e 1882
f66d95af 1883 }
1884 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCEta) && fhMCE[mcEta])
1885 {
1886 fhMCE [mcEta] ->Fill(ecluster);
4c8f7c2e 1887 fhMCPt [mcEta] ->Fill(ptcluster);
1888 fhMCPhi[mcEta] ->Fill(ecluster,phicluster);
1889 fhMCEta[mcEta] ->Fill(ecluster,etacluster);
1890
1891 fhMC2E[mcEta] ->Fill(ecluster, eprim);
1892 fhMC2Pt[mcEta] ->Fill(ptcluster, ptprim);
d9105d92 1893 fhMCDeltaE[mcEta] ->Fill(ecluster,eprim-ecluster);
1894 fhMCDeltaPt[mcEta]->Fill(ptcluster,ptprim-ptcluster);
4c8f7c2e 1895
f66d95af 1896 }
3d5d5078 1897 }
f66d95af 1898 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCAntiNeutron) && fhMCE[mcAntiNeutron])
3d5d5078 1899 {
1900 fhMCE [mcAntiNeutron] ->Fill(ecluster);
4c8f7c2e 1901 fhMCPt [mcAntiNeutron] ->Fill(ptcluster);
1902 fhMCPhi[mcAntiNeutron] ->Fill(ecluster,phicluster);
1903 fhMCEta[mcAntiNeutron] ->Fill(ecluster,etacluster);
3d5d5078 1904
4c8f7c2e 1905 fhMC2E[mcAntiNeutron] ->Fill(ecluster, eprim);
1906 fhMC2Pt[mcAntiNeutron] ->Fill(ptcluster, ptprim);
d9105d92 1907 fhMCDeltaE[mcAntiNeutron] ->Fill(ecluster,eprim-ecluster);
1908 fhMCDeltaPt[mcAntiNeutron]->Fill(ptcluster,ptprim-ptcluster);
4c8f7c2e 1909
3d5d5078 1910 }
f66d95af 1911 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCAntiProton) && fhMCE[mcAntiProton])
3d5d5078 1912 {
1913 fhMCE [mcAntiProton] ->Fill(ecluster);
4c8f7c2e 1914 fhMCPt [mcAntiProton] ->Fill(ptcluster);
1915 fhMCPhi[mcAntiProton] ->Fill(ecluster,phicluster);
1916 fhMCEta[mcAntiProton] ->Fill(ecluster,etacluster);
4c8f7c2e 1917
1918 fhMC2E[mcAntiProton] ->Fill(ecluster, eprim);
1919 fhMC2Pt[mcAntiProton] ->Fill(ptcluster, ptprim);
d9105d92 1920 fhMCDeltaE[mcAntiProton] ->Fill(ecluster,eprim-ecluster);
1921 fhMCDeltaPt[mcAntiProton]->Fill(ecluster,ptprim-ptcluster);
4c8f7c2e 1922
3d5d5078 1923 }
f66d95af 1924 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCElectron) && fhMCE[mcElectron])
3d5d5078 1925 {
1926 fhMCE [mcElectron] ->Fill(ecluster);
4c8f7c2e 1927 fhMCPt [mcElectron] ->Fill(ptcluster);
1928 fhMCPhi[mcElectron] ->Fill(ecluster,phicluster);
1929 fhMCEta[mcElectron] ->Fill(ecluster,etacluster);
1930
1931 fhMC2E[mcElectron] ->Fill(ecluster, eprim);
1932 fhMC2Pt[mcElectron] ->Fill(ptcluster, ptprim);
d9105d92 1933 fhMCDeltaE[mcElectron] ->Fill(ecluster,eprim-ecluster);
1934 fhMCDeltaPt[mcElectron]->Fill(ecluster,ptprim-ptcluster);
3d5d5078 1935 }
f66d95af 1936 else if( fhMCE[mcOther]){
3d5d5078 1937 fhMCE [mcOther] ->Fill(ecluster);
4c8f7c2e 1938 fhMCPt [mcOther] ->Fill(ptcluster);
1939 fhMCPhi[mcOther] ->Fill(ecluster,phicluster);
1940 fhMCEta[mcOther] ->Fill(ecluster,etacluster);
521636d2 1941
4c8f7c2e 1942 fhMC2E[mcOther] ->Fill(ecluster, eprim);
1943 fhMC2Pt[mcOther] ->Fill(ptcluster, ptprim);
d9105d92 1944 fhMCDeltaE[mcOther] ->Fill(ecluster,eprim-ecluster);
1945 fhMCDeltaPt[mcOther]->Fill(ecluster,ptprim-ptcluster);
4c8f7c2e 1946
f8006433 1947 // printf(" AliAnaPhoton::MakeAnalysisFillHistograms() - Label %d, pT %2.3f Unknown, bits set: ",
1948 // ph->GetLabel(),ph->Pt());
1949 // for(Int_t i = 0; i < 20; i++) {
1950 // if(GetMCAnalysisUtils()->CheckTagBit(tag,i)) printf(" %d, ",i);
1951 // }
1952 // printf("\n");
1953
577d9801 1954 }
521636d2 1955
577d9801 1956 }//Histograms with MC
521636d2 1957
577d9801 1958 }// aod loop
521636d2 1959
1c5acb87 1960}
1961
1962
1963//__________________________________________________________________
1964void AliAnaPhoton::Print(const Option_t * opt) const
1965{
477d6cee 1966 //Print some relevant parameters set for the analysis
1967
1968 if(! opt)
1969 return;
1970
1971 printf("**** Print %s %s ****\n", GetName(), GetTitle() ) ;
1972 AliAnaPartCorrBaseClass::Print(" ");
a3aebfff 1973
477d6cee 1974 printf("Calorimeter = %s\n", fCalorimeter.Data()) ;
1975 printf("Min Distance to Bad Channel = %2.1f\n",fMinDist);
1976 printf("Min Distance to Bad Channel 2 = %2.1f\n",fMinDist2);
1977 printf("Min Distance to Bad Channel 3 = %2.1f\n",fMinDist3);
a3aebfff 1978 printf("Reject clusters with a track matched = %d\n",fRejectTrackMatch);
4cf55759 1979 printf("Time Cut: %3.1f < TOF < %3.1f\n", fTimeCutMin, fTimeCutMax);
2ac125bf 1980 printf("Number of cells in cluster is > %d \n", fNCellsCut);
477d6cee 1981 printf(" \n") ;
1c5acb87 1982
1983}