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