]> 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 more correlation histograms to study exotic clusters, remove others, fill MIP...
[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>
477d6cee 34//#include <Riostream.h>
123fc3bd 35#include "TParticle.h"
6175da48 36#include "TDatabasePDG.h"
1c5acb87 37
38// --- Analysis system ---
39#include "AliAnaPhoton.h"
40#include "AliCaloTrackReader.h"
123fc3bd 41#include "AliStack.h"
1c5acb87 42#include "AliCaloPID.h"
6639984f 43#include "AliMCAnalysisUtils.h"
ff45398a 44#include "AliFiducialCut.h"
0ae57829 45#include "AliVCluster.h"
591cc579 46#include "AliAODMCParticle.h"
c8fe2783 47#include "AliMixedEvent.h"
48
1c5acb87 49
50ClassImp(AliAnaPhoton)
51
52//____________________________________________________________________________
521636d2 53AliAnaPhoton::AliAnaPhoton() :
54 AliAnaPartCorrBaseClass(), fCalorimeter(""),
55 fMinDist(0.), fMinDist2(0.), fMinDist3(0.),
56 fRejectTrackMatch(0), fTimeCutMin(-1), fTimeCutMax(999999),
57 fNCellsCut(0), fFillSSHistograms(kFALSE),
58 fCheckConversion(kFALSE), fRemoveConvertedPair(kFALSE),
59 fAddConvertedPairsToAOD(kFALSE),
60 fMassCut(0), fConvAsymCut(1.), fConvDEtaCut(2.),
61 fConvDPhiMinCut(-1.), fConvDPhiMaxCut(7.),
62
c4a7d28a 63 // Histograms
521636d2 64 fhNCellsE(0), fhEPhoton(0), fhPtPhoton(0),
65 fhPhiPhoton(0), fhEtaPhoton(0),
66 fhEtaPhiPhoton(0), fhEtaPhi05Photon(0),
67
c4a7d28a 68 // Conversion histograms
521636d2 69 fhPtPhotonConv(0), fhEtaPhiPhotonConv(0), fhEtaPhi05PhotonConv(0),
70 fhConvDeltaEta(0), fhConvDeltaPhi(0), fhConvDeltaEtaPhi(0),
71 fhConvAsym(0), fhConvPt(0),
72 fhConvDistEta(0), fhConvDistEn(0), fhConvDistMass(0),
73 fhConvDistEtaCutEta(0), fhConvDistEnCutEta(0), fhConvDistMassCutEta(0),
74 fhConvDistEtaCutMass(0), fhConvDistEnCutMass(0),
75 fhConvDistEtaCutAsy(0), fhConvDistEnCutAsy(0),
76
77 //Shower shape histograms
78 fhDispE(0), fhLam0E(0), fhLam1E(0),
79 fhdDispE(0), fhdLam0E(0), fhdLam1E(0),
80 fhDispETRD(0), fhLam0ETRD(0), fhLam1ETRD(0),
81 fhdDispETRD(0), fhdLam0ETRD(0), fhdLam1ETRD(0),
82
83 fhNCellsLam0LowE(0), fhNCellsLam1LowE(0), fhNCellsDispLowE(0),
84 fhNCellsLam0HighE(0), fhNCellsLam1HighE(0), fhNCellsDispHighE(0),
85 fhNCellsdLam0LowE(0), fhNCellsdLam1LowE(0), fhNCellsdDispLowE(0),
86 fhNCellsdLam0HighE(0), fhNCellsdLam1HighE(0), fhNCellsdDispHighE(0),
87
88 fhEtaLam0LowE(0), fhPhiLam0LowE(0),
89 fhEtaLam0HighE(0), fhPhiLam0HighE(0),
90 fhLam0DispLowE(0), fhLam0DispHighE(0),
91 fhLam1Lam0LowE(0), fhLam1Lam0HighE(0),
92 fhDispLam1LowE(0), fhDispLam1HighE(0),
93 fhEtadLam0LowE(0), fhPhidLam0LowE(0),
94 fhEtadLam0HighE(0), fhPhidLam0HighE(0),
95 fhdLam0dDispLowE(0), fhdLam0dDispHighE(0),
96 fhdLam1dLam0LowE(0), fhdLam1dLam0HighE(0),
97 fhdDispdLam1LowE(0), fhdDispdLam1HighE(0),
98
99 //MC histograms
100 fhDeltaE(0), fhDeltaPt(0),
101 fhRatioE(0), fhRatioPt(0),
102 fh2E(0), fh2Pt(0),
103
104 // Conversion MC histograms
105 fhPtConversionTagged(0), fhPtAntiNeutronTagged(0),
106 fhPtAntiProtonTagged(0), fhPtUnknownTagged(0),
107 fhEtaPhiConversion(0), fhEtaPhi05Conversion(0),
108
109 fhConvDeltaEtaMCConversion(0), fhConvDeltaPhiMCConversion(0), fhConvDeltaEtaPhiMCConversion(0),
110 fhConvAsymMCConversion(0), fhConvPtMCConversion(0),
111 fhConvDispersionMCConversion(0), fhConvM02MCConversion(0),
112
113 fhConvDeltaEtaMCAntiNeutron(0), fhConvDeltaPhiMCAntiNeutron(0), fhConvDeltaEtaPhiMCAntiNeutron(0),
114 fhConvAsymMCAntiNeutron(0), fhConvPtMCAntiNeutron(0),
115 fhConvDispersionMCAntiNeutron(0), fhConvM02MCAntiNeutron(0),
116 fhConvDeltaEtaMCAntiProton(0), fhConvDeltaPhiMCAntiProton(0), fhConvDeltaEtaPhiMCAntiProton(0),
117 fhConvAsymMCAntiProton(0), fhConvPtMCAntiProton(0),
118 fhConvDispersionMCAntiProton(0), fhConvM02MCAntiProton(0),
119 fhConvDeltaEtaMCString(0), fhConvDeltaPhiMCString(0), fhConvDeltaEtaPhiMCString(0),
120 fhConvAsymMCString(0), fhConvPtMCString(0),
121 fhConvDispersionMCString(0), fhConvM02MCString(0),
3d5d5078 122 fhConvDistMCConversion(0), fhConvDistMCConversionCuts(0),
123
124 // Photon SS MC histograms
125 fhMCPhotonELambda0NoOverlap(0), fhMCPhotonELambda0TwoOverlap(0), fhMCPhotonELambda0NOverlap(0),
126 fhMCPhotonEdLambda0NoOverlap(0), fhMCPhotonEdLambda0TwoOverlap(0), fhMCPhotonEdLambda0NOverlap(0),
127
128 //Embedding
129 fhEmbeddedSignalFractionEnergy(0),
130 fhEmbedPhotonELambda0FullSignal(0), fhEmbedPhotonEdLambda0FullSignal(0),
131 fhEmbedPhotonELambda0MostlySignal(0), fhEmbedPhotonEdLambda0MostlySignal(0),
132 fhEmbedPhotonELambda0MostlyBkg(0), fhEmbedPhotonEdLambda0MostlyBkg(0),
133 fhEmbedPhotonELambda0FullBkg(0), fhEmbedPhotonEdLambda0FullBkg(0),
134 fhEmbedPi0ELambda0FullSignal(0), fhEmbedPi0EdLambda0FullSignal(0),
135 fhEmbedPi0ELambda0MostlySignal(0), fhEmbedPi0EdLambda0MostlySignal(0),
136 fhEmbedPi0ELambda0MostlyBkg(0), fhEmbedPi0EdLambda0MostlyBkg(0),
137 fhEmbedPi0ELambda0FullBkg(0), fhEmbedPi0EdLambda0FullBkg(0)
1c5acb87 138{
139 //default ctor
140
521636d2 141 for(Int_t i = 0; i < 12; i++){
3d5d5078 142 fhPtMC [i] = 0;
143 fhMCE [i] = 0;
521636d2 144 fhPhiMC[i] = 0;
145 fhEtaMC[i] = 0;
146 }
147
3d5d5078 148 for(Int_t i = 0; i < 7; i++){
149 fhPtPrimMC [i] = 0;
150 fhEPrimMC [i] = 0;
151 fhPhiPrimMC[i] = 0;
152 fhYPrimMC [i] = 0;
153
154 fhPtPrimMCAcc [i] = 0;
155 fhEPrimMCAcc [i] = 0;
156 fhPhiPrimMCAcc[i] = 0;
157 fhYPrimMCAcc [i] = 0;
158 }
159
160 for(Int_t i = 0; i < 6; i++){
161 fhMCELambda0 [i] = 0;
162 fhMCELambda1 [i] = 0;
163 fhMCEDispersion [i] = 0;
164 fhMCEdLambda0 [i] = 0;
165 fhMCEdLambda1 [i] = 0;
166 fhMCEdDispersion[i] = 0;
521636d2 167 }
168
1c5acb87 169 //Initialize parameters
170 InitParameters();
171
5ae09196 172}//____________________________________________________________________________
1c5acb87 173AliAnaPhoton::~AliAnaPhoton()
174{
175 //dtor
176
177}
178
c4a7d28a 179//__________________________________________________________________
180Bool_t AliAnaPhoton::ClusterSelected(AliVCluster* calo, TLorentzVector mom)
181{
182 //Select clusters if they pass different cuts
183 if(GetDebug() > 2)
184 printf("AliAnaPhoton::ClusterSelected() Current Event %d; Before selection : E %2.2f, pT %2.2f, Ecl %2.2f, phi %2.2f, eta %2.2f\n",
185 GetReader()->GetEventNumber(),
186 mom.E(), mom.Pt(),calo->E(),mom.Phi()*TMath::RadToDeg(),mom.Eta());
187
188 //.......................................
189 //If too small or big energy, skip it
190 if(mom.E() < GetMinEnergy() || mom.E() > GetMaxEnergy() ) return kFALSE ;
191 if(GetDebug() > 2) printf("\t Cluster %d Pass E Cut \n",calo->GetID());
192
193 //.......................................
194 // TOF cut, BE CAREFUL WITH THIS CUT
195 Double_t tof = calo->GetTOF()*1e9;
196 if(tof < fTimeCutMin || tof > fTimeCutMax) return kFALSE;
197 if(GetDebug() > 2) printf("\t Cluster %d Pass Time Cut \n",calo->GetID());
198
199 //.......................................
200 if(calo->GetNCells() <= fNCellsCut && GetReader()->GetDataType() != AliCaloTrackReader::kMC) return kFALSE;
201 if(GetDebug() > 2) printf("\t Cluster %d Pass NCell Cut \n",calo->GetID());
202
203 //.......................................
204 //Check acceptance selection
205 if(IsFiducialCutOn()){
206 Bool_t in = GetFiducialCut()->IsInFiducialCut(mom,fCalorimeter) ;
207 if(! in ) return kFALSE ;
208 }
209 if(GetDebug() > 2) printf("Fiducial cut passed \n");
210
211 //.......................................
212 //Skip matched clusters with tracks
213 if(fRejectTrackMatch){
214 if(IsTrackMatched(calo)) {
215 if(GetDebug() > 2) printf("\t Reject track-matched clusters\n");
216 return kFALSE ;
217 }
218 else
219 if(GetDebug() > 2) printf(" Track-matching cut passed \n");
220 }// reject matched clusters
221
222 //.......................................
223 //Check Distance to Bad channel, set bit.
224 Double_t distBad=calo->GetDistanceToBadChannel() ; //Distance to bad channel
225 if(distBad < 0.) distBad=9999. ; //workout strange convension dist = -1. ;
226 if(distBad < fMinDist) {//In bad channel (PHOS cristal size 2.2x2.2 cm), EMCAL ( cell units )
227 return kFALSE ;
228 }
229 else if(GetDebug() > 2) printf("\t Bad channel cut passed %4.2f > %2.2f \n",distBad, fMinDist);
230 //printf("Cluster %d Pass Bad Dist Cut \n",icalo);
231
232 if(GetDebug() > 0)
233 printf("AliAnaPhoton::ClusterSelected() Current Event %d; After selection : E %2.2f, pT %2.2f, Ecl %2.2f, phi %2.2f, eta %2.2f\n",
234 GetReader()->GetEventNumber(),
235 mom.E(), mom.Pt(),calo->E(),mom.Phi()*TMath::RadToDeg(),mom.Eta());
236
237 //All checks passed, cluster selected
238 return kTRUE;
239
240}
241
3d5d5078 242//_____________________________________________________________
243void AliAnaPhoton::FillAcceptanceHistograms(){
244 //Fill acceptance histograms if MC data is available
245
246 if(GetReader()->ReadStack()){
247 AliStack * stack = GetMCStack();
248 if(stack){
249 for(Int_t i=0 ; i<stack->GetNtrack(); i++){
250 TParticle * prim = stack->Particle(i) ;
251 Int_t pdg = prim->GetPdgCode();
252 //printf("i %d, %s %d %s %d \n",i, stack->Particle(i)->GetName(), stack->Particle(i)->GetPdgCode(),
253 // prim->GetName(), prim->GetPdgCode());
254
255 if(pdg == 22){
256
257 // Get tag of this particle photon from fragmentation, decay, prompt ...
258 Int_t tag = GetMCAnalysisUtils()->CheckOrigin(i,GetReader(), 0);
259 if(!GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPhoton)){
260 //A conversion photon from a hadron, skip this kind of photon
261 // 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,
262 // GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCConversion),
263 // GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0),
264 // GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCOther),
265 // GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCElectron),
266 // GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCUnknown),
267 // GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCMuon),
268 // GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPion),
269 // GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCProton),
270 // GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCAntiNeutron),
271 // GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCKaon),
272 // GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCAntiProton),
273 // GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCAntiNeutron));
274
275 return;
276 }
277
278 //Get photon kinematics
279 if(prim->Energy() == TMath::Abs(prim->Pz())) continue ; //Protection against floating point exception
280
281 Double_t photonY = 0.5*TMath::Log((prim->Energy()-prim->Pz())/(prim->Energy()+prim->Pz())) ;
282 Double_t photonE = prim->Energy() ;
283 Double_t photonPt = prim->Pt() ;
284 Double_t photonPhi = TMath::RadToDeg()*prim->Phi() ;
285 if(photonPhi < 0) photonPhi+=TMath::TwoPi();
286 Double_t photonEta = prim->Eta() ;
287
288
289 //Check if photons hit the Calorimeter
290 TLorentzVector lv;
291 prim->Momentum(lv);
292 Bool_t inacceptance = kFALSE;
293 if (fCalorimeter == "PHOS"){
294 if(GetPHOSGeometry() && GetCaloUtils()->IsPHOSGeoMatrixSet()){
295 Int_t mod ;
296 Double_t x,z ;
297 if(GetPHOSGeometry()->ImpactOnEmc(prim,mod,z,x))
298 inacceptance = kTRUE;
299 if(GetDebug() > 2) printf("In %s Real acceptance? %d\n",fCalorimeter.Data(),inacceptance);
300 }
301 else{
302 if(GetFiducialCut()->IsInFiducialCut(lv,fCalorimeter))
303 inacceptance = kTRUE ;
304 if(GetDebug() > 2) printf("In %s fiducial cut acceptance? %d\n",fCalorimeter.Data(),inacceptance);
305 }
306 }
307 else if(fCalorimeter == "EMCAL" && GetCaloUtils()->IsEMCALGeoMatrixSet()){
308 if(GetEMCALGeometry()){
309
310 Int_t absID=0;
311
312 GetEMCALGeometry()->GetAbsCellIdFromEtaPhi(prim->Eta(),prim->Phi(),absID);
313
314 if( absID >= 0)
315 inacceptance = kTRUE;
316
317 // if(GetEMCALGeometry()->Impact(phot1) && GetEMCALGeometry()->Impact(phot2))
318 // inacceptance = kTRUE;
319 if(GetDebug() > 2) printf("In %s Real acceptance? %d\n",fCalorimeter.Data(),inacceptance);
320 }
321 else{
322 if(GetFiducialCut()->IsInFiducialCut(lv,fCalorimeter))
323 inacceptance = kTRUE ;
324 if(GetDebug() > 2) printf("In %s fiducial cut acceptance? %d\n",fCalorimeter.Data(),inacceptance);
325 }
326 } //In EMCAL
327
328 //Fill histograms
329
330 fhYPrimMC[mcPhoton]->Fill(photonPt, photonY) ;
331 if(TMath::Abs(photonY) < 1.0)
332 {
333 fhEPrimMC [mcPhoton]->Fill(photonE ) ;
334 fhPtPrimMC [mcPhoton]->Fill(photonPt) ;
335 fhPhiPrimMC[mcPhoton]->Fill(photonE , photonPhi) ;
336 fhYPrimMC[mcPhoton] ->Fill(photonE , photonEta) ;
337 }
338 if(inacceptance){
339 fhEPrimMCAcc[mcPhoton] ->Fill(photonE ) ;
340 fhPtPrimMCAcc[mcPhoton] ->Fill(photonPt) ;
341 fhPhiPrimMCAcc[mcPhoton]->Fill(photonE , photonPhi) ;
342 fhYPrimMCAcc[mcPhoton] ->Fill(photonE , photonY) ;
343 }//Accepted
344
345 //Origin of photon
346 if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPrompt))
347 {
348 fhYPrimMC[mcPrompt]->Fill(photonPt, photonY) ;
349 if(TMath::Abs(photonY) < 1.0){
350 fhEPrimMC [mcPrompt]->Fill(photonE ) ;
351 fhPtPrimMC [mcPrompt]->Fill(photonPt) ;
352 fhPhiPrimMC[mcPrompt]->Fill(photonE , photonPhi) ;
353 fhYPrimMC[mcPrompt] ->Fill(photonE , photonEta) ;
354 }
355 if(inacceptance){
356 fhEPrimMCAcc[mcPrompt] ->Fill(photonE ) ;
357 fhPtPrimMCAcc[mcPrompt] ->Fill(photonPt) ;
358 fhPhiPrimMCAcc[mcPrompt]->Fill(photonE , photonPhi) ;
359 fhYPrimMCAcc[mcPrompt] ->Fill(photonE , photonY) ;
360 }//Accepted
361 }
362 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCFragmentation))
363 {
364 fhYPrimMC[mcFragmentation]->Fill(photonPt, photonY) ;
365 if(TMath::Abs(photonY) < 1.0){
366 fhEPrimMC [mcFragmentation]->Fill(photonE ) ;
367 fhPtPrimMC [mcFragmentation]->Fill(photonPt) ;
368 fhPhiPrimMC[mcFragmentation]->Fill(photonE , photonPhi) ;
369 fhYPrimMC[mcFragmentation] ->Fill(photonE , photonEta) ;
370 }
371 if(inacceptance){
372 fhEPrimMCAcc[mcFragmentation] ->Fill(photonE ) ;
373 fhPtPrimMCAcc[mcFragmentation] ->Fill(photonPt) ;
374 fhPhiPrimMCAcc[mcFragmentation]->Fill(photonE , photonPhi) ;
375 fhYPrimMCAcc[mcFragmentation] ->Fill(photonE , photonY) ;
376 }//Accepted
377 }
378 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCISR))
379 {
380 fhYPrimMC[mcISR]->Fill(photonPt, photonY) ;
381 if(TMath::Abs(photonY) < 1.0){
382 fhEPrimMC [mcISR]->Fill(photonE ) ;
383 fhPtPrimMC [mcISR]->Fill(photonPt) ;
384 fhPhiPrimMC[mcISR]->Fill(photonE , photonPhi) ;
385 fhYPrimMC[mcISR]->Fill(photonE , photonEta) ;
386 }
387 if(inacceptance){
388 fhEPrimMCAcc[mcISR] ->Fill(photonE ) ;
389 fhPtPrimMCAcc[mcISR] ->Fill(photonPt) ;
390 fhPhiPrimMCAcc[mcISR]->Fill(photonE , photonPhi) ;
391 fhYPrimMCAcc[mcISR] ->Fill(photonE , photonY) ;
392 }//Accepted
393 }
394 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0Decay))
395 {
396 fhYPrimMC[mcPi0Decay]->Fill(photonPt, photonY) ;
397 if(TMath::Abs(photonY) < 1.0){
398 fhEPrimMC [mcPi0Decay]->Fill(photonE ) ;
399 fhPtPrimMC [mcPi0Decay]->Fill(photonPt) ;
400 fhPhiPrimMC[mcPi0Decay]->Fill(photonE , photonPhi) ;
401 fhYPrimMC[mcPi0Decay] ->Fill(photonE , photonEta) ;
402 }
403 if(inacceptance){
404 fhEPrimMCAcc[mcPi0Decay] ->Fill(photonE ) ;
405 fhPtPrimMCAcc[mcPi0Decay] ->Fill(photonPt) ;
406 fhPhiPrimMCAcc[mcPi0Decay]->Fill(photonE , photonPhi) ;
407 fhYPrimMCAcc[mcPi0Decay] ->Fill(photonE , photonY) ;
408 }//Accepted
409 }
410 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCEtaDecay) ||
411 GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCOtherDecay))
412 {
413 fhYPrimMC[mcOtherDecay]->Fill(photonPt, photonY) ;
414 if(TMath::Abs(photonY) < 1.0){
415 fhEPrimMC [mcOtherDecay]->Fill(photonE ) ;
416 fhPtPrimMC [mcOtherDecay]->Fill(photonPt) ;
417 fhPhiPrimMC[mcOtherDecay]->Fill(photonE , photonPhi) ;
418 fhYPrimMC[mcOtherDecay] ->Fill(photonE , photonEta) ;
419 }
420 if(inacceptance){
421 fhEPrimMCAcc[mcOtherDecay] ->Fill(photonE ) ;
422 fhPtPrimMCAcc[mcOtherDecay] ->Fill(photonPt) ;
423 fhPhiPrimMCAcc[mcOtherDecay]->Fill(photonE , photonPhi) ;
424 fhYPrimMCAcc[mcOtherDecay] ->Fill(photonE , photonY) ;
425 }//Accepted
426 }
427 else
428 {
429 fhYPrimMC[mcOther]->Fill(photonPt, photonY) ;
430 if(TMath::Abs(photonY) < 1.0){
431 fhEPrimMC [mcOther]->Fill(photonE ) ;
432 fhPtPrimMC [mcOther]->Fill(photonPt) ;
433 fhPhiPrimMC[mcOther]->Fill(photonE , photonPhi) ;
434 fhYPrimMC[mcOther] ->Fill(photonE , photonEta) ;
435 }
436 if(inacceptance){
437 fhEPrimMCAcc[mcOther] ->Fill(photonE ) ;
438 fhPtPrimMCAcc[mcOther] ->Fill(photonPt) ;
439 fhPhiPrimMCAcc[mcOther]->Fill(photonE , photonPhi) ;
440 fhYPrimMCAcc[mcOther] ->Fill(photonE , photonY) ;
441 }//Accepted
442 }//Other origin
443 }// Primary photon
444 }//loop on primaries
445 }//stack exists and data is MC
446 }//read stack
447 else if(GetReader()->ReadAODMCParticles()){
448 TClonesArray * mcparticles = GetReader()->GetAODMCParticles(0);
449 if(mcparticles){
450 Int_t nprim = mcparticles->GetEntriesFast();
451
452 for(Int_t i=0; i < nprim; i++)
453 {
454 AliAODMCParticle * prim = (AliAODMCParticle *) mcparticles->At(i);
455
456 Int_t pdg = prim->GetPdgCode();
457
458 if(pdg == 22){
459
460 // Get tag of this particle photon from fragmentation, decay, prompt ...
461 Int_t tag = GetMCAnalysisUtils()->CheckOrigin(i,GetReader(), 0);
462 if(!GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPhoton)){
463 //A conversion photon from a hadron, skip this kind of photon
464// 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,
465// GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCConversion),
466// GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0),
467// GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCOther),
468// GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCElectron),
469// GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCUnknown),
470// GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCMuon),
471// GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPion),
472// GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCProton),
473// GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCAntiNeutron),
474// GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCKaon),
475// GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCAntiProton),
476// GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCAntiNeutron));
477
478 return;
479 }
480
481 //Get photon kinematics
482 if(prim->E() == TMath::Abs(prim->Pz())) continue ; //Protection against floating point exception
483
484 Double_t photonY = 0.5*TMath::Log((prim->E()-prim->Pz())/(prim->E()+prim->Pz())) ;
485 Double_t photonE = prim->E() ;
486 Double_t photonPt = prim->Pt() ;
487 Double_t photonPhi = TMath::RadToDeg()*prim->Phi() ;
488 if(photonPhi < 0) photonPhi+=TMath::TwoPi();
489 Double_t photonEta = prim->Eta() ;
490
491 //Check if photons hit the Calorimeter
492 TLorentzVector lv;
493 lv.SetPxPyPzE(prim->Px(),prim->Py(),prim->Pz(),prim->E());
494 Bool_t inacceptance = kFALSE;
495 if (fCalorimeter == "PHOS"){
496 if(GetPHOSGeometry() && GetCaloUtils()->IsPHOSGeoMatrixSet()){
497 Int_t mod ;
498 Double_t x,z ;
499 Double_t vtx[]={prim->Xv(),prim->Yv(),prim->Zv()};
500 if(GetPHOSGeometry()->ImpactOnEmc(vtx, prim->Theta(),prim->Phi(),mod,z,x))
501 inacceptance = kTRUE;
502 if(GetDebug() > 2) printf("In %s Real acceptance? %d\n",fCalorimeter.Data(),inacceptance);
503 }
504 else{
505 if(GetFiducialCut()->IsInFiducialCut(lv,fCalorimeter))
506 inacceptance = kTRUE ;
507 if(GetDebug() > 2) printf("In %s fiducial cut acceptance? %d\n",fCalorimeter.Data(),inacceptance);
508 }
509 }
510 else if(fCalorimeter == "EMCAL" && GetCaloUtils()->IsEMCALGeoMatrixSet()){
511 if(GetEMCALGeometry()){
512
513 Int_t absID=0;
514
515 GetEMCALGeometry()->GetAbsCellIdFromEtaPhi(prim->Eta(),prim->Phi(),absID);
516
517 if( absID >= 0)
518 inacceptance = kTRUE;
519
520 if(GetDebug() > 2) printf("In %s Real acceptance? %d\n",fCalorimeter.Data(),inacceptance);
521 }
522 else{
523 if(GetFiducialCut()->IsInFiducialCut(lv,fCalorimeter))
524 inacceptance = kTRUE ;
525 if(GetDebug() > 2) printf("In %s fiducial cut acceptance? %d\n",fCalorimeter.Data(),inacceptance);
526 }
527 } //In EMCAL
528
529 //Fill histograms
530
531 fhYPrimMC[mcPhoton]->Fill(photonPt, photonY) ;
532 if(TMath::Abs(photonY) < 1.0)
533 {
534 fhEPrimMC [mcPhoton]->Fill(photonE ) ;
535 fhPtPrimMC [mcPhoton]->Fill(photonPt) ;
536 fhPhiPrimMC[mcPhoton]->Fill(photonE , photonPhi) ;
537 fhYPrimMC[mcPhoton]->Fill(photonE , photonEta) ;
538 }
539 if(inacceptance){
540 fhEPrimMCAcc[mcPhoton] ->Fill(photonE ) ;
541 fhPtPrimMCAcc[mcPhoton] ->Fill(photonPt) ;
542 fhPhiPrimMCAcc[mcPhoton]->Fill(photonE , photonPhi) ;
543 fhYPrimMCAcc[mcPhoton] ->Fill(photonE , photonY) ;
544 }//Accepted
545
546
547 //Origin of photon
548 if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPrompt))
549 {
550 fhYPrimMC[mcPrompt]->Fill(photonPt, photonY) ;
551 if(TMath::Abs(photonY) < 1.0){
552 fhEPrimMC [mcPrompt]->Fill(photonE ) ;
553 fhPtPrimMC [mcPrompt]->Fill(photonPt) ;
554 fhPhiPrimMC[mcPrompt]->Fill(photonE , photonPhi) ;
555 fhYPrimMC[mcPrompt]->Fill(photonE , photonEta) ;
556 }
557 if(inacceptance){
558 fhEPrimMCAcc[mcPrompt] ->Fill(photonE ) ;
559 fhPtPrimMCAcc[mcPrompt] ->Fill(photonPt) ;
560 fhPhiPrimMCAcc[mcPrompt]->Fill(photonE , photonPhi) ;
561 fhYPrimMCAcc[mcPrompt] ->Fill(photonE , photonY) ;
562 }//Accepted
563 }
564 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCFragmentation))
565 {
566 fhYPrimMC[mcFragmentation]->Fill(photonPt, photonY) ;
567 if(TMath::Abs(photonY) < 1.0){
568 fhEPrimMC [mcFragmentation]->Fill(photonE ) ;
569 fhPtPrimMC [mcFragmentation]->Fill(photonPt) ;
570 fhPhiPrimMC[mcFragmentation]->Fill(photonE , photonPhi) ;
571 fhYPrimMC[mcFragmentation]->Fill(photonE , photonEta) ;
572 }
573 if(inacceptance){
574 fhEPrimMCAcc[mcFragmentation] ->Fill(photonE ) ;
575 fhPtPrimMCAcc[mcFragmentation] ->Fill(photonPt) ;
576 fhPhiPrimMCAcc[mcFragmentation]->Fill(photonE , photonPhi) ;
577 fhYPrimMCAcc[mcFragmentation] ->Fill(photonE , photonY) ;
578 }//Accepted
579 }
580 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCISR))
581 {
582 fhYPrimMC[mcISR]->Fill(photonPt, photonY) ;
583 if(TMath::Abs(photonY) < 1.0){
584 fhEPrimMC [mcISR]->Fill(photonE ) ;
585 fhPtPrimMC [mcISR]->Fill(photonPt) ;
586 fhPhiPrimMC[mcISR]->Fill(photonE , photonPhi) ;
587 fhYPrimMC[mcISR]->Fill(photonE , photonEta) ;
588 }
589 if(inacceptance){
590 fhEPrimMCAcc[mcISR] ->Fill(photonE ) ;
591 fhPtPrimMCAcc[mcISR] ->Fill(photonPt) ;
592 fhPhiPrimMCAcc[mcISR]->Fill(photonE , photonPhi) ;
593 fhYPrimMCAcc[mcISR] ->Fill(photonE , photonY) ;
594 }//Accepted
595 }
596 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0Decay))
597 {
598 fhYPrimMC[mcPi0Decay]->Fill(photonPt, photonY) ;
599 if(TMath::Abs(photonY) < 1.0){
600 fhEPrimMC [mcPi0Decay]->Fill(photonE ) ;
601 fhPtPrimMC [mcPi0Decay]->Fill(photonPt) ;
602 fhPhiPrimMC[mcPi0Decay]->Fill(photonE , photonPhi) ;
603 fhYPrimMC[mcPi0Decay]->Fill(photonE , photonEta) ;
604 }
605 if(inacceptance){
606 fhEPrimMCAcc[mcPi0Decay] ->Fill(photonE ) ;
607 fhPtPrimMCAcc[mcPi0Decay] ->Fill(photonPt) ;
608 fhPhiPrimMCAcc[mcPi0Decay]->Fill(photonE , photonPhi) ;
609 fhYPrimMCAcc[mcPi0Decay] ->Fill(photonE , photonY) ;
610 }//Accepted
611 }
612 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCEtaDecay) ||
613 GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCOtherDecay))
614 {
615 fhYPrimMC[mcOtherDecay]->Fill(photonPt, photonY) ;
616 if(TMath::Abs(photonY) < 1.0){
617 fhEPrimMC [mcOtherDecay]->Fill(photonE ) ;
618 fhPtPrimMC [mcOtherDecay]->Fill(photonPt) ;
619 fhPhiPrimMC[mcOtherDecay]->Fill(photonE , photonPhi) ;
620 fhYPrimMC[mcOtherDecay]->Fill(photonE , photonEta) ;
621 }
622 if(inacceptance){
623 fhEPrimMCAcc[mcOtherDecay] ->Fill(photonE ) ;
624 fhPtPrimMCAcc[mcOtherDecay] ->Fill(photonPt) ;
625 fhPhiPrimMCAcc[mcOtherDecay]->Fill(photonE , photonPhi) ;
626 fhYPrimMCAcc[mcOtherDecay] ->Fill(photonE , photonY) ;
627 }//Accepted
628 }
629 else
630 {
631 fhYPrimMC[mcOther]->Fill(photonPt, photonY) ;
632 if(TMath::Abs(photonY) < 1.0){
633 fhEPrimMC [mcOther]->Fill(photonE ) ;
634 fhPtPrimMC [mcOther]->Fill(photonPt) ;
635 fhPhiPrimMC[mcOther]->Fill(photonE , photonPhi) ;
636 fhYPrimMC[mcOther]->Fill(photonE , photonEta) ;
637 }
638 if(inacceptance){
639 fhEPrimMCAcc[mcOther] ->Fill(photonE ) ;
640 fhPtPrimMCAcc[mcOther] ->Fill(photonPt) ;
641 fhPhiPrimMCAcc[mcOther]->Fill(photonE , photonPhi) ;
642 fhYPrimMCAcc[mcOther] ->Fill(photonE , photonY) ;
643 }//Accepted
644 }//Other origin
645 }// Primary photon
646 }//loop on primaries
647
648 }//mc array exists and data is MC
649 } // read AOD MC
650}
521636d2 651
652//__________________________________________________________________
653void AliAnaPhoton::FillShowerShapeHistograms(AliVCluster* cluster, const Int_t mcTag){
654
655 //Fill cluster Shower Shape histograms
656
657 if(!fFillSSHistograms || GetMixedEvent()) return;
658
659 Float_t energy = cluster->E();
660 Int_t ncells = cluster->GetNCells();
661 Int_t ncells2 = ncells*ncells;
662 Float_t lambda0 = cluster->GetM02();
663 Float_t lambda1 = cluster->GetM20();
664 Float_t disp = cluster->GetDispersion()*cluster->GetDispersion();
665
666 TLorentzVector mom;
667 if(GetReader()->GetDataType() != AliCaloTrackReader::kMC){
668 cluster->GetMomentum(mom,GetVertex(0)) ;}//Assume that come from vertex in straight line
669 else{
670 Double_t vertex[]={0,0,0};
671 cluster->GetMomentum(mom,vertex) ;
672 }
673
674 Float_t eta = mom.Eta();
675 Float_t phi = mom.Phi();
676 if(phi < 0) phi+=TMath::TwoPi();
677
678 fhLam0E ->Fill(energy,lambda0);
679 fhLam1E ->Fill(energy,lambda1);
680 fhDispE ->Fill(energy,disp);
681 fhdLam0E->Fill(energy,lambda0/ncells2);
682 fhdLam1E->Fill(energy,lambda1/ncells2);
683 fhdDispE->Fill(energy,disp/ncells2);
684
685 if(fCalorimeter == "EMCAL" && GetModuleNumber(cluster) > 5){
686 fhLam0ETRD->Fill(energy,lambda0);
687 fhLam1ETRD->Fill(energy,lambda1);
688 fhDispETRD->Fill(energy,disp);
689 fhdLam0ETRD->Fill(energy,lambda0/ncells2);
690 fhdLam1ETRD->Fill(energy,lambda1/ncells2);
691 fhdDispETRD->Fill(energy,disp/ncells2);
692 }
693
694 if(energy < 2){
695 fhNCellsLam0LowE ->Fill(ncells,lambda0);
696 fhNCellsLam1LowE ->Fill(ncells,lambda1);
697 fhNCellsDispLowE ->Fill(ncells,disp);
698 fhNCellsdLam0LowE->Fill(ncells,lambda0/ncells2);
699 fhNCellsdLam1LowE->Fill(ncells,lambda1/ncells2);
700 fhNCellsdDispLowE->Fill(ncells,disp/ncells2);
701
702 fhLam1Lam0LowE ->Fill(lambda1,lambda0);
703 fhLam0DispLowE ->Fill(lambda0,disp);
704 fhDispLam1LowE ->Fill(disp,lambda1);
705 fhEtaLam0LowE ->Fill(eta,lambda0);
706 fhPhiLam0LowE ->Fill(phi,lambda0);
707
708 fhdLam1dLam0LowE->Fill(lambda1/ncells2,lambda0/ncells2);
709 fhdLam0dDispLowE->Fill(lambda0/ncells2,disp/ncells2);
710 fhdDispdLam1LowE->Fill(disp/ncells2,lambda1/ncells2);
711 fhEtadLam0LowE ->Fill(eta,lambda0/ncells2);
712 fhPhidLam0LowE ->Fill(phi,lambda0/ncells2);
713 }
714 else {
715 fhNCellsLam0HighE ->Fill(ncells,lambda0);
716 fhNCellsLam1HighE ->Fill(ncells,lambda1);
717 fhNCellsDispHighE ->Fill(ncells,disp);
718 fhNCellsdLam0HighE->Fill(ncells,lambda0/ncells2);
719 fhNCellsdLam1HighE->Fill(ncells,lambda1/ncells2);
720 fhNCellsdDispHighE->Fill(ncells,disp/ncells2);
721
722 fhLam1Lam0HighE ->Fill(lambda1,lambda0);
723 fhLam0DispHighE ->Fill(lambda0,disp);
724 fhDispLam1HighE ->Fill(disp,lambda1);
725 fhEtaLam0HighE ->Fill(eta, lambda0);
726 fhPhiLam0HighE ->Fill(phi, lambda0);
727
728 fhdLam1dLam0HighE->Fill(lambda1/ncells2,lambda0/ncells2);
729 fhdLam0dDispHighE->Fill(lambda0/ncells2,disp/ncells2);
730 fhdDispdLam1HighE->Fill(disp/ncells2,lambda1/ncells2);
731 fhEtadLam0HighE ->Fill(eta, lambda0/ncells2);
732 fhPhidLam0HighE ->Fill(phi, lambda0/ncells2);
733 }
734
735 if(IsDataMC()){
3d5d5078 736
737
738 //Fill histograms to check shape of embedded clusters
739 Float_t fraction = 0;
740 if(GetReader()->IsEmbeddedClusterSelectionOn()){
741 AliVCaloCells * cells = 0;
742 if(fCalorimeter == "EMCAL") cells = GetReader()->GetEMCALCells();
743 else cells = GetReader()->GetPHOSCells(); //only working for EMCAL, but in case in future it works
744
745 Float_t clusterE = 0; // recalculate in case corrections applied.
746 Float_t cellE = 0;
747 for(Int_t icell = 0; icell < cluster->GetNCells(); icell++){
748 cellE = cells->GetCellAmplitude(cluster->GetCellAbsId(icell));
749 clusterE+=cellE;
750 fraction+=cellE*cluster->GetCellAmplitudeFraction(icell);
751 }
752
753 //Fraction of total energy due to the embedded signal
754 fraction/=clusterE;
755
756 if(GetDebug() > 1 ) printf("AliAnaPhoton::FillShowerShapeHistogram() - Energy fraction of embedded signal %2.3f, Energy %2.3f\n",fraction, clusterE);
757
758 fhEmbeddedSignalFractionEnergy->Fill(clusterE,fraction);
759
760 } // embedded fraction
761
521636d2 762 if( GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPhoton) &&
3d5d5078 763 !GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCConversion) &&
764 !GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPi0) &&
765 !GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCEta)){
766 fhMCELambda0[mcssPhoton] ->Fill(energy, lambda0);
767 fhMCEdLambda0[mcssPhoton] ->Fill(energy, lambda0/ncells2);
768 fhMCELambda1[mcssPhoton] ->Fill(energy, lambda1);
769 fhMCEdLambda1[mcssPhoton] ->Fill(energy, lambda1/ncells2);
770 fhMCEDispersion[mcssPhoton] ->Fill(energy, disp);
771 fhMCEdDispersion[mcssPhoton]->Fill(energy, disp/ncells2);
772
773
774 if(!GetReader()->IsEmbeddedClusterSelectionOn()){
775 //Check particle overlaps in cluster
776
777 //Compare the primary depositing more energy with the rest, if no photon/electron as comon ancestor (conversions), count as other particle
778 Int_t ancPDG = 0, ancStatus = -1;
779 TLorentzVector momentum; TVector3 prodVertex;
780 Int_t ancLabel = 0;
781 Int_t noverlaps = 1;
782 for (UInt_t ilab = 0; ilab < cluster->GetNLabels(); ilab++ ) {
783 ancLabel = GetMCAnalysisUtils()->CheckCommonAncestor(cluster->GetLabels()[0],cluster->GetLabels()[ilab], GetReader(),ancPDG,ancStatus,momentum,prodVertex);
784 if(ancPDG!=22 && TMath::Abs(ancPDG)!=11) noverlaps++;
785 }
786
787 if(noverlaps == 1){
788 fhMCPhotonELambda0NoOverlap ->Fill(energy, lambda0);
789 fhMCPhotonEdLambda0NoOverlap ->Fill(energy, lambda0/ncells2);
790 }
791 else if(noverlaps == 2){
792 fhMCPhotonELambda0TwoOverlap ->Fill(energy, lambda0);
793 fhMCPhotonEdLambda0TwoOverlap->Fill(energy, lambda0/ncells2);
794 }
795 else if(noverlaps > 2){
796 fhMCPhotonELambda0NOverlap ->Fill(energy, lambda0);
797 fhMCPhotonEdLambda0NOverlap ->Fill(energy, lambda0/ncells2);
798 }
799 else {
800 printf("AliAnaPhoton::FillShowerShapeHistogram() - n overlaps = %d!!", noverlaps);
801 }
802 }//No embedding
803
804 //Fill histograms to check shape of embedded clusters
805 if(GetReader()->IsEmbeddedClusterSelectionOn()){
806
807 if (fraction > 0.9)
808 {
809 fhEmbedPhotonELambda0FullSignal ->Fill(energy, lambda0);
810 fhEmbedPhotonEdLambda0FullSignal ->Fill(energy, lambda0/ncells2);
811 }
812 else if(fraction > 0.5)
813 {
814 fhEmbedPhotonELambda0MostlySignal ->Fill(energy, lambda0);
815 fhEmbedPhotonEdLambda0MostlySignal->Fill(energy, lambda0/ncells2);
816 }
817 else if(fraction > 0.1)
818 {
819 fhEmbedPhotonELambda0MostlyBkg ->Fill(energy, lambda0);
820 fhEmbedPhotonEdLambda0MostlyBkg ->Fill(energy, lambda0/ncells2);
821 }
822 else
823 {
824 fhEmbedPhotonELambda0FullBkg ->Fill(energy, lambda0);
825 fhEmbedPhotonEdLambda0FullBkg ->Fill(energy, lambda0/ncells2);
826 }
827 } // embedded
828
521636d2 829 }//photon no conversion
830 else if ( GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCElectron)){
3d5d5078 831 fhMCELambda0[mcssElectron] ->Fill(energy, lambda0);
832 fhMCEdLambda0[mcssElectron] ->Fill(energy, lambda0/ncells2);
833 fhMCELambda1[mcssElectron] ->Fill(energy, lambda1);
834 fhMCEdLambda1[mcssElectron] ->Fill(energy, lambda1/ncells2);
835 fhMCEDispersion[mcssElectron] ->Fill(energy, disp);
836 fhMCEdDispersion[mcssElectron]->Fill(energy, disp/ncells2);
521636d2 837 }//electron
3d5d5078 838 else if ( GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPhoton) &&
521636d2 839 GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCConversion) ){
3d5d5078 840 fhMCELambda0[mcssConversion] ->Fill(energy, lambda0);
841 fhMCEdLambda0[mcssConversion] ->Fill(energy, lambda0/ncells2);
842 fhMCELambda1[mcssConversion] ->Fill(energy, lambda1);
843 fhMCEdLambda1[mcssConversion] ->Fill(energy, lambda1/ncells2);
844 fhMCEDispersion[mcssConversion] ->Fill(energy, disp);
845 fhMCEdDispersion[mcssConversion]->Fill(energy, disp/ncells2);
846
521636d2 847 }//conversion photon
848 else if ( GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPi0) ){
3d5d5078 849 fhMCELambda0[mcssPi0] ->Fill(energy, lambda0);
850 fhMCEdLambda0[mcssPi0] ->Fill(energy, lambda0/ncells2);
851 fhMCELambda1[mcssPi0] ->Fill(energy, lambda1);
852 fhMCEdLambda1[mcssPi0] ->Fill(energy, lambda1/ncells2);
853 fhMCEDispersion[mcssPi0] ->Fill(energy, disp);
854 fhMCEdDispersion[mcssPi0]->Fill(energy, disp/ncells2);
855
856 //Fill histograms to check shape of embedded clusters
857 if(GetReader()->IsEmbeddedClusterSelectionOn()){
858
859 if (fraction > 0.9)
860 {
861 fhEmbedPi0ELambda0FullSignal ->Fill(energy, lambda0);
862 fhEmbedPi0EdLambda0FullSignal ->Fill(energy, lambda0/ncells2);
863 }
864 else if(fraction > 0.5)
865 {
866 fhEmbedPi0ELambda0MostlySignal ->Fill(energy, lambda0);
867 fhEmbedPi0EdLambda0MostlySignal->Fill(energy, lambda0/ncells2);
868 }
869 else if(fraction > 0.1)
870 {
871 fhEmbedPi0ELambda0MostlyBkg ->Fill(energy, lambda0);
872 fhEmbedPi0EdLambda0MostlyBkg ->Fill(energy, lambda0/ncells2);
873 }
874 else
875 {
876 fhEmbedPi0ELambda0FullBkg ->Fill(energy, lambda0);
877 fhEmbedPi0EdLambda0FullBkg ->Fill(energy, lambda0/ncells2);
878 }
879 } // embedded
880
521636d2 881 }//pi0
3d5d5078 882 else if ( GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCEta) ){
883 fhMCELambda0[mcssEta] ->Fill(energy, lambda0);
884 fhMCEdLambda0[mcssEta] ->Fill(energy, lambda0/ncells2);
885 fhMCELambda1[mcssEta] ->Fill(energy, lambda1);
886 fhMCEdLambda1[mcssEta] ->Fill(energy, lambda1/ncells2);
887 fhMCEDispersion[mcssEta] ->Fill(energy, disp);
888 fhMCEdDispersion[mcssEta]->Fill(energy, disp/ncells2);
889 }//eta
521636d2 890 else {
3d5d5078 891 fhMCELambda0[mcssOther] ->Fill(energy, lambda0);
892 fhMCEdLambda0[mcssOther] ->Fill(energy, lambda0/ncells2);
893 fhMCELambda1[mcssOther] ->Fill(energy, lambda1);
894 fhMCEdLambda1[mcssOther] ->Fill(energy, lambda1/ncells2);
895 fhMCEDispersion[mcssOther] ->Fill(energy, disp);
896 fhMCEdDispersion[mcssOther]->Fill(energy, disp/ncells2);
521636d2 897 }//other particles
898
899 }//MC data
900
901}
902
0c1383b5 903//________________________________________________________________________
904TObjString * AliAnaPhoton::GetAnalysisCuts()
905{
906 //Save parameters used for analysis
907 TString parList ; //this will be list of parameters used for this analysis.
5ae09196 908 const Int_t buffersize = 255;
909 char onePar[buffersize] ;
0c1383b5 910
5ae09196 911 snprintf(onePar,buffersize,"--- AliAnaPhoton ---\n") ;
0c1383b5 912 parList+=onePar ;
5ae09196 913 snprintf(onePar,buffersize,"Calorimeter: %s\n",fCalorimeter.Data()) ;
0c1383b5 914 parList+=onePar ;
5ae09196 915 snprintf(onePar,buffersize,"fMinDist =%2.2f (Minimal distance to bad channel to accept cluster) \n",fMinDist) ;
0c1383b5 916 parList+=onePar ;
5ae09196 917 snprintf(onePar,buffersize,"fMinDist2=%2.2f (Cuts on Minimal distance to study acceptance evaluation) \n",fMinDist2) ;
0c1383b5 918 parList+=onePar ;
5ae09196 919 snprintf(onePar,buffersize,"fMinDist3=%2.2f (One more cut on distance used for acceptance-efficiency study) \n",fMinDist3) ;
0c1383b5 920 parList+=onePar ;
5ae09196 921 snprintf(onePar,buffersize,"fRejectTrackMatch: %d\n",fRejectTrackMatch) ;
0c1383b5 922 parList+=onePar ;
41121cfe 923 snprintf(onePar,buffersize,"Conversion Selection: fConvAsymCut %1.2f, fConvDEtaCut %1.2f fConvDPhiCut (%1.2f,%1.2f)\n",
924 fConvAsymCut, fConvDEtaCut, fConvDPhiMinCut, fConvDPhiMaxCut) ;
6175da48 925 parList+=onePar ;
0c1383b5 926
927 //Get parameters set in base class.
928 parList += GetBaseParametersList() ;
929
930 //Get parameters set in PID class.
931 parList += GetCaloPID()->GetPIDParametersList() ;
932
933 //Get parameters set in FiducialCut class (not available yet)
934 //parlist += GetFidCut()->GetFidCutParametersList()
935
936 return new TObjString(parList) ;
937}
938
1c5acb87 939//________________________________________________________________________
940TList * AliAnaPhoton::GetCreateOutputObjects()
941{
477d6cee 942 // Create histograms to be saved in output file and
943 // store them in outputContainer
944 TList * outputContainer = new TList() ;
945 outputContainer->SetName("PhotonHistos") ;
4a745797 946
521636d2 947 Int_t nptbins = GetHistoPtBins(); Float_t ptmax = GetHistoPtMax(); Float_t ptmin = GetHistoPtMin();
948 Int_t nphibins = GetHistoPhiBins(); Float_t phimax = GetHistoPhiMax(); Float_t phimin = GetHistoPhiMin();
949 Int_t netabins = GetHistoEtaBins(); Float_t etamax = GetHistoEtaMax(); Float_t etamin = GetHistoEtaMin();
950 Int_t ssbins = GetHistoShowerShapeBins(); Float_t ssmax = GetHistoShowerShapeMax(); Float_t ssmin = GetHistoShowerShapeMin();
e1e62b89 951 Int_t nbins = GetHistoNClusterCellBins(); Int_t nmax = GetHistoNClusterCellMax(); Int_t nmin = GetHistoNClusterCellMin();
521636d2 952
e1e62b89 953 fhNCellsE = new TH2F ("hNCellsE","# of cells in cluster vs E of clusters", nptbins,ptmin,ptmax, nbins,nmin,nmax);
c4a7d28a 954 fhNCellsE->SetXTitle("E (GeV)");
955 fhNCellsE->SetYTitle("# of cells in cluster");
956 outputContainer->Add(fhNCellsE);
6175da48 957
20218aea 958 fhEPhoton = new TH1F("hEPhoton","Number of #gamma over calorimeter vs energy",nptbins,ptmin,ptmax);
959 fhEPhoton->SetYTitle("N");
960 fhEPhoton->SetXTitle("E_{#gamma}(GeV)");
961 outputContainer->Add(fhEPhoton) ;
962
963 fhPtPhoton = new TH1F("hPtPhoton","Number of #gamma over calorimeter vs p_{T}",nptbins,ptmin,ptmax);
477d6cee 964 fhPtPhoton->SetYTitle("N");
965 fhPtPhoton->SetXTitle("p_{T #gamma}(GeV/c)");
966 outputContainer->Add(fhPtPhoton) ;
967
968 fhPhiPhoton = new TH2F
20218aea 969 ("hPhiPhoton","#phi_{#gamma} vs p_{T}",nptbins,ptmin,ptmax,nphibins,phimin,phimax);
6175da48 970 fhPhiPhoton->SetYTitle("#phi (rad)");
477d6cee 971 fhPhiPhoton->SetXTitle("p_{T #gamma} (GeV/c)");
972 outputContainer->Add(fhPhiPhoton) ;
973
974 fhEtaPhoton = new TH2F
20218aea 975 ("hEtaPhoton","#eta_{#gamma} vs p_{T}",nptbins,ptmin,ptmax,netabins,etamin,etamax);
477d6cee 976 fhEtaPhoton->SetYTitle("#eta");
977 fhEtaPhoton->SetXTitle("p_{T #gamma} (GeV/c)");
978 outputContainer->Add(fhEtaPhoton) ;
979
6175da48 980 fhEtaPhiPhoton = new TH2F
981 ("hEtaPhiPhoton","#eta vs #phi",netabins,etamin,etamax,nphibins,phimin,phimax);
982 fhEtaPhiPhoton->SetYTitle("#phi (rad)");
983 fhEtaPhiPhoton->SetXTitle("#eta");
984 outputContainer->Add(fhEtaPhiPhoton) ;
20218aea 985 if(GetMinPt() < 0.5){
986 fhEtaPhi05Photon = new TH2F
987 ("hEtaPhi05Photon","#eta vs #phi, E > 0.5",netabins,etamin,etamax,nphibins,phimin,phimax);
988 fhEtaPhi05Photon->SetYTitle("#phi (rad)");
989 fhEtaPhi05Photon->SetXTitle("#eta");
990 outputContainer->Add(fhEtaPhi05Photon) ;
991 }
6175da48 992
993 //Conversion
20218aea 994 if(fCheckConversion){
521636d2 995
996 fhEtaPhiConversion = new TH2F
997 ("hEtaPhiConversion","#eta vs #phi for converted clusters",netabins,etamin,etamax,nphibins,phimin,phimax);
998 fhEtaPhiConversion->SetYTitle("#phi (rad)");
999 fhEtaPhiConversion->SetXTitle("#eta");
1000 outputContainer->Add(fhEtaPhiConversion) ;
1001 if(GetMinPt() < 0.5){
1002 fhEtaPhi05Conversion = new TH2F
1003 ("hEtaPhi05Conversion","#eta vs #phi, E > 0.5, for converted clusters",netabins,etamin,etamax,nphibins,phimin,phimax);
1004 fhEtaPhi05Conversion->SetYTitle("#phi (rad)");
1005 fhEtaPhi05Conversion->SetXTitle("#eta");
1006 outputContainer->Add(fhEtaPhi05Conversion) ;
1007 }
1008
20218aea 1009 fhPtPhotonConv = new TH1F("hPtPhotonConv","Number of #gamma over calorimeter, conversion",nptbins,ptmin,ptmax);
1010 fhPtPhotonConv->SetYTitle("N");
1011 fhPtPhotonConv->SetXTitle("p_{T #gamma}(GeV/c)");
1012 outputContainer->Add(fhPtPhotonConv) ;
1013
1014 fhEtaPhiPhotonConv = new TH2F
1015 ("hEtaPhiPhotonConv","#eta vs #phi",netabins,etamin,etamax,nphibins,phimin,phimax);
1016 fhEtaPhiPhotonConv->SetYTitle("#phi (rad)");
1017 fhEtaPhiPhotonConv->SetXTitle("#eta");
1018 outputContainer->Add(fhEtaPhiPhotonConv) ;
1019 if(GetMinPt() < 0.5){
1020 fhEtaPhi05PhotonConv = new TH2F
1021 ("hEtaPhi05PhotonConv","#eta vs #phi, E > 0.5",netabins,etamin,etamax,nphibins,phimin,phimax);
1022 fhEtaPhi05PhotonConv->SetYTitle("#phi (rad)");
1023 fhEtaPhi05PhotonConv->SetXTitle("#eta");
1024 outputContainer->Add(fhEtaPhi05PhotonConv) ;
1025 }
1026
1027 fhConvDeltaEta = new TH2F
1028 ("hConvDeltaEta","#Delta #eta of selected conversion pairs",100,0,fMassCut,netabins*2,-0.5,0.5);
1029 fhConvDeltaEta->SetYTitle("#Delta #eta");
1030 fhConvDeltaEta->SetXTitle("Pair Mass (GeV/c^2)");
1031 outputContainer->Add(fhConvDeltaEta) ;
1032
1033 fhConvDeltaPhi = new TH2F
1034 ("hConvDeltaPhi","#Delta #phi of selected conversion pairs",100,0,fMassCut,nphibins*2,-0.5,0.5);
1035 fhConvDeltaPhi->SetYTitle("#Delta #phi");
1036 fhConvDeltaPhi->SetXTitle("Pair Mass (GeV/c^2)");
1037 outputContainer->Add(fhConvDeltaPhi) ;
1038
1039 fhConvDeltaEtaPhi = new TH2F
1040 ("hConvDeltaEtaPhi","#Delta #eta vs #Delta #phi of selected conversion pairs",netabins,-0.5,0.5,nphibins,-0.5,0.5);
1041 fhConvDeltaEtaPhi->SetYTitle("#Delta #phi");
1042 fhConvDeltaEtaPhi->SetXTitle("#Delta #eta");
1043 outputContainer->Add(fhConvDeltaEtaPhi) ;
1044
1045 fhConvAsym = new TH2F
1046 ("hConvAsym","Asymmetry of selected conversion pairs",100,0,fMassCut,100,0,1);
1047 fhConvAsym->SetYTitle("Asymmetry");
1048 fhConvAsym->SetXTitle("Pair Mass (GeV/c^2)");
1049 outputContainer->Add(fhConvAsym) ;
1050
1051 fhConvPt = new TH2F
1052 ("hConvPt","p_{T} of selected conversion pairs",100,0,fMassCut,100,0.,10.);
1053 fhConvPt->SetYTitle("Pair p_{T} (GeV/c)");
1054 fhConvPt->SetXTitle("Pair Mass (GeV/c^2)");
1055 outputContainer->Add(fhConvPt) ;
c4a7d28a 1056
1057 fhConvDistEta = new TH2F
1058 ("hConvDistEta","distance to conversion vertex",100,-0.7,0.7,100,0.,5.);
1059 fhConvDistEta->SetXTitle("#eta");
1060 fhConvDistEta->SetYTitle(" distance (m)");
1061 outputContainer->Add(fhConvDistEta) ;
1062
1063 fhConvDistEn = new TH2F
1064 ("hConvDistEn","distance to conversion vertex",nptbins,ptmin,ptmax,100,0.,5.);
1065 fhConvDistEn->SetXTitle("E (GeV)");
1066 fhConvDistEn->SetYTitle(" distance (m)");
1067 outputContainer->Add(fhConvDistEn) ;
1068
1069 fhConvDistMass = new TH2F
1070 ("hConvDistMass","distance to conversion vertex",100,0,fMassCut,100,0.,5.);
1071 fhConvDistMass->SetXTitle("m (GeV/c^2)");
1072 fhConvDistMass->SetYTitle(" distance (m)");
1073 outputContainer->Add(fhConvDistMass) ;
1074
1075 fhConvDistEtaCutEta = new TH2F
1076 ("hConvDistEtaCutEta","distance to conversion vertex, dEta < 0.05",100,-0.7,0.7,100,0.,5.);
1077 fhConvDistEtaCutEta->SetXTitle("#eta");
1078 fhConvDistEtaCutEta->SetYTitle(" distance (m)");
1079 outputContainer->Add(fhConvDistEtaCutEta) ;
1080
1081 fhConvDistEnCutEta = new TH2F
1082 ("hConvDistEnCutEta","distance to conversion vertex, dEta < 0.05",nptbins,ptmin,ptmax,100,0.,5.);
1083 fhConvDistEnCutEta->SetXTitle("E (GeV)");
1084 fhConvDistEnCutEta->SetYTitle(" distance (m)");
1085 outputContainer->Add(fhConvDistEnCutEta) ;
1086
1087 fhConvDistMassCutEta = new TH2F
1088 ("hConvDistMassCutEta","distance to conversion vertex, dEta < 0.05",100,0,fMassCut,100,0.,5.);
1089 fhConvDistMassCutEta->SetXTitle("m (GeV/c^2)");
1090 fhConvDistMassCutEta->SetYTitle(" distance (m)");
1091 outputContainer->Add(fhConvDistMassCutEta) ;
1092
1093 fhConvDistEtaCutMass = new TH2F
1094 ("hConvDistEtaCutMass","distance to conversion vertex, dEta < 0.05, m < 10 MeV",100,-0.7,0.7,100,0.,5.);
1095 fhConvDistEtaCutMass->SetXTitle("#eta");
1096 fhConvDistEtaCutMass->SetYTitle(" distance (m)");
1097 outputContainer->Add(fhConvDistEtaCutMass) ;
1098
1099 fhConvDistEnCutMass = new TH2F
1100 ("hConvDistEnCutMass","distance to conversion vertex, dEta < 0.05, m < 10 MeV",nptbins,ptmin,ptmax,100,0.,5.);
1101 fhConvDistEnCutMass->SetXTitle("E (GeV)");
1102 fhConvDistEnCutMass->SetYTitle(" distance (m)");
1103 outputContainer->Add(fhConvDistEnCutMass) ;
1104
1105 fhConvDistEtaCutAsy = new TH2F
1106 ("hConvDistEtaCutAsy","distance to conversion vertex, dEta < 0.05, m < 10 MeV, A < 0.1",100,-0.7,0.7,100,0.,5.);
1107 fhConvDistEtaCutAsy->SetXTitle("#eta");
1108 fhConvDistEtaCutAsy->SetYTitle(" distance (m)");
1109 outputContainer->Add(fhConvDistEtaCutAsy) ;
1110
1111 fhConvDistEnCutAsy = new TH2F
1112 ("hConvDistEnCutAsy","distance to conversion vertex, dEta < 0.05, m < 10 MeV, A < 0.1",nptbins,ptmin,ptmax,100,0.,5.);
1113 fhConvDistEnCutAsy->SetXTitle("E (GeV)");
1114 fhConvDistEnCutAsy->SetYTitle(" distance (m)");
1115 outputContainer->Add(fhConvDistEnCutAsy) ;
1116
521636d2 1117 } // check conversion
1118
1119
1120 //Shower shape
1121 if(fFillSSHistograms){
1122
1123 fhLam0E = new TH2F ("hLam0E","#lambda_{0}^{2} vs E", nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
1124 fhLam0E->SetYTitle("#lambda_{0}^{2}");
1125 fhLam0E->SetXTitle("E (GeV)");
1126 outputContainer->Add(fhLam0E);
1127
1128 fhLam1E = new TH2F ("hLam1E","#lambda_{1}^{2} vs E", nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
1129 fhLam1E->SetYTitle("#lambda_{1}^{2}");
1130 fhLam1E->SetXTitle("E (GeV)");
1131 outputContainer->Add(fhLam1E);
1132
1133 fhDispE = new TH2F ("hDispE"," dispersion^{2} vs E", nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
1134 fhDispE->SetYTitle("D^{2}");
1135 fhDispE->SetXTitle("E (GeV) ");
1136 outputContainer->Add(fhDispE);
1137
1138 fhdLam0E = new TH2F ("hdLam0E","#lambda_{0}^{2}/N_{cells}^{2} vs E", nptbins,ptmin,ptmax,ssbins,ssmin,ssmax/50);
1139 fhdLam0E->SetYTitle("d#lambda_{0}^{2}");
1140 fhdLam0E->SetXTitle("E (GeV)");
1141 outputContainer->Add(fhdLam0E);
1142
1143 fhdLam1E = new TH2F ("hdLam1E","#lambda_{1}^{2}/N_{cells}^{2} vs E", nptbins,ptmin,ptmax,ssbins,ssmin,ssmax/50);
1144 fhdLam1E->SetYTitle("d#lambda_{1}^{2}");
1145 fhdLam1E->SetXTitle("E (GeV)");
1146 outputContainer->Add(fhdLam1E);
1147
1148 fhdDispE = new TH2F ("hdDispE"," dispersion^{2}/N_{cells}^{2} vs E", nptbins,ptmin,ptmax,ssbins,ssmin,ssmax/50);
1149 fhdDispE->SetYTitle("dD^{2}");
1150 fhdDispE->SetXTitle("E (GeV) ");
1151 outputContainer->Add(fhdDispE);
1152
1153
1154 if(fCalorimeter == "EMCAL"){
1155 fhLam0ETRD = new TH2F ("hLam0ETRD","#lambda_{0}^{2} vs E, EMCAL SM covered by TRD", nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
1156 fhLam0ETRD->SetYTitle("#lambda_{0}^{2}");
1157 fhLam0ETRD->SetXTitle("E (GeV)");
1158 outputContainer->Add(fhLam0ETRD);
1159
1160 fhLam1ETRD = new TH2F ("hLam1ETRD","#lambda_{1}^{2} vs E, EMCAL SM covered by TRD", nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
1161 fhLam1ETRD->SetYTitle("#lambda_{1}^{2}");
1162 fhLam1ETRD->SetXTitle("E (GeV)");
1163 outputContainer->Add(fhLam1ETRD);
1164
1165 fhDispETRD = new TH2F ("hDispETRD"," dispersion^{2} vs E, EMCAL SM covered by TRD", nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
1166 fhDispETRD->SetYTitle("Dispersion^{2}");
1167 fhDispETRD->SetXTitle("E (GeV) ");
1168 outputContainer->Add(fhDispETRD);
1169
1170 fhdLam0ETRD = new TH2F ("hdLam0ETRD","#lambda_{0}^{2}/N_{cells}^{2} vs E, EMCAL SM covered by TRD", nptbins,ptmin,ptmax,ssbins,ssmin,ssmax/50);
1171 fhdLam0ETRD->SetYTitle("d#lambda_{0}^{2}");
1172 fhdLam0ETRD->SetXTitle("E (GeV)");
1173 outputContainer->Add(fhdLam0ETRD);
1174
1175 fhdLam1ETRD = new TH2F ("hdLam1ETRD","#lambda_{1}^{2}/N_{cells}^{2} vs E, EMCAL SM covered by TRD", nptbins,ptmin,ptmax,ssbins,ssmin,ssmax/50);
1176 fhdLam1ETRD->SetYTitle("d#lambda_{1}^{2}");
1177 fhdLam1ETRD->SetXTitle("E (GeV)");
1178 outputContainer->Add(fhdLam1ETRD);
1179
1180 fhdDispETRD = new TH2F ("hdDispETRD"," dispersion^{2}/N_{cells}^{2} vs E, EMCAL SM covered by TRD", nptbins,ptmin,ptmax,ssbins,ssmin,ssmax/50);
1181 fhdDispETRD->SetYTitle("dD^{2}");
1182 fhdDispETRD->SetXTitle("E (GeV) ");
1183 outputContainer->Add(fhdDispETRD);
1184
1185 }
1186
1187 fhNCellsLam0LowE = new TH2F ("hNCellsLam0LowE","N_{cells} in cluster vs #lambda_{0}^{2}, E < 2 GeV", 20,0, 20, ssbins,ssmin,ssmax);
1188 fhNCellsLam0LowE->SetXTitle("N_{cells}");
1189 fhNCellsLam0LowE->SetYTitle("#lambda_{0}^{2}");
1190 outputContainer->Add(fhNCellsLam0LowE);
1191
1192 fhNCellsLam0HighE = new TH2F ("hNCellsLam0HighE","N_{cells} in cluster vs #lambda_{0}^{2}, E > 2 GeV", 20,0, 20, ssbins,ssmin,ssmax);
1193 fhNCellsLam0HighE->SetXTitle("N_{cells}");
1194 fhNCellsLam0HighE->SetYTitle("#lambda_{0}^{2}");
1195 outputContainer->Add(fhNCellsLam0HighE);
1196
1197 fhNCellsLam1LowE = new TH2F ("hNCellsLam1LowE","N_{cells} in cluster vs #lambda_{1}^{2}, E < 2 GeV", 20,0, 20, ssbins,ssmin,ssmax);
1198 fhNCellsLam1LowE->SetXTitle("N_{cells}");
1199 fhNCellsLam1LowE->SetYTitle("#lambda_{0}^{2}");
1200 outputContainer->Add(fhNCellsLam1LowE);
1201
1202 fhNCellsLam1HighE = new TH2F ("hNCellsLam1HighE","N_{cells} in cluster vs #lambda_{1}^{2}, E > 2 GeV", 20,0, 20, ssbins,ssmin,ssmax);
1203 fhNCellsLam1HighE->SetXTitle("N_{cells}");
1204 fhNCellsLam1HighE->SetYTitle("#lambda_{0}^{2}");
1205 outputContainer->Add(fhNCellsLam1HighE);
1206
1207 fhNCellsDispLowE = new TH2F ("hNCellsDispLowE","N_{cells} in cluster vs dispersion^{2}, E < 2 GeV", 20,0, 20, ssbins,ssmin,ssmax);
1208 fhNCellsDispLowE->SetXTitle("N_{cells}");
1209 fhNCellsDispLowE->SetYTitle("D^{2}");
1210 outputContainer->Add(fhNCellsDispLowE);
1211
1212 fhNCellsDispHighE = new TH2F ("hNCellsDispHighE","N_{cells} in cluster vs dispersion^{2}, E < 2 GeV", 20,0, 20, ssbins,ssmin,ssmax);
1213 fhNCellsDispHighE->SetXTitle("N_{cells}");
1214 fhNCellsDispHighE->SetYTitle("D^{2}");
1215 outputContainer->Add(fhNCellsDispHighE);
1216
1217
1218 fhNCellsdLam0LowE = new TH2F ("hNCellsdLam0LowE","N_{cells} in cluster vs #lambda_{0}^{2}/N_{cells}^{2}, E < 2 GeV", 20,0, 20, ssbins,ssmin,ssmax/50);
1219 fhNCellsdLam0LowE->SetXTitle("N_{cells}");
1220 fhNCellsdLam0LowE->SetYTitle("#lambda_{0}^{2}");
1221 outputContainer->Add(fhNCellsdLam0LowE);
1222
1223 fhNCellsdLam0HighE = new TH2F ("hNCellsdLam0HighE","N_{cells} in cluster vs #lambda_{0}^{2}/N_{cells}^{2}, E > 2 GeV", 20,0, 20, ssbins,ssmin,ssmax/50);
1224 fhNCellsdLam0HighE->SetXTitle("N_{cells}");
1225 fhNCellsdLam0HighE->SetYTitle("#lambda_{0}^{2}");
1226 outputContainer->Add(fhNCellsdLam0HighE);
1227
1228 fhNCellsdLam1LowE = new TH2F ("hNCellsdLam1LowE","N_{cells} in cluster vs #lambda_{1}^{2}/N_{cells}^{2}, E < 2 GeV", 20,0, 20, ssbins,ssmin,ssmax/50);
1229 fhNCellsdLam1LowE->SetXTitle("N_{cells}");
1230 fhNCellsdLam1LowE->SetYTitle("#lambda_{0}^{2}");
1231 outputContainer->Add(fhNCellsdLam1LowE);
1232
1233 fhNCellsdLam1HighE = new TH2F ("hNCellsdLam1HighE","N_{cells} in cluster vs #lambda_{1}^{2}/N_{cells}^{2}, E > 2 GeV", 20,0, 20, ssbins,ssmin,ssmax/50);
1234 fhNCellsdLam1HighE->SetXTitle("N_{cells}");
1235 fhNCellsdLam1HighE->SetYTitle("#lambda_{0}^{2}");
1236 outputContainer->Add(fhNCellsdLam1HighE);
1237
1238 fhNCellsdDispLowE = new TH2F ("hNCellsdDispLowE","N_{cells} in cluster vs dispersion^{2}/N_{cells}^{2}, E < 2 GeV", 20,0, 20, ssbins,ssmin,ssmax/50);
1239 fhNCellsdDispLowE->SetXTitle("N_{cells}");
1240 fhNCellsdDispLowE->SetYTitle("D^{2}");
1241 outputContainer->Add(fhNCellsdDispLowE);
1242
1243 fhNCellsdDispHighE = new TH2F ("hNCellsdDispHighE","N_{cells} in cluster vs dispersion^{2}/N_{cells}^{2}, E < 2 GeV", 20,0, 20, ssbins,ssmin,ssmax/50);
1244 fhNCellsdDispHighE->SetXTitle("N_{cells}");
1245 fhNCellsdDispHighE->SetYTitle("D^{2}");
1246 outputContainer->Add(fhNCellsdDispHighE);
1247
1248
1249 fhEtaLam0LowE = new TH2F ("hEtaLam0LowE","#eta vs #lambda_{0}^{2}, E < 2 GeV", netabins,etamin,etamax, ssbins,ssmin,ssmax);
1250 fhEtaLam0LowE->SetYTitle("#lambda_{0}^{2}");
1251 fhEtaLam0LowE->SetXTitle("#eta");
1252 outputContainer->Add(fhEtaLam0LowE);
1253
1254 fhPhiLam0LowE = new TH2F ("hPhiLam0LowE","#phi vs #lambda_{0}^{2}, E < 2 GeV", nphibins,phimin,phimax, ssbins,ssmin,ssmax);
1255 fhPhiLam0LowE->SetYTitle("#lambda_{0}^{2}");
1256 fhPhiLam0LowE->SetXTitle("#phi");
1257 outputContainer->Add(fhPhiLam0LowE);
1258
1259 fhEtaLam0HighE = new TH2F ("hEtaLam0HighE","#eta vs #lambda_{0}^{2}, E > 2 GeV", netabins,etamin,etamax, ssbins,ssmin,ssmax);
1260 fhEtaLam0HighE->SetYTitle("#lambda_{0}^{2}");
1261 fhEtaLam0HighE->SetXTitle("#eta");
1262 outputContainer->Add(fhEtaLam0HighE);
1263
1264 fhPhiLam0HighE = new TH2F ("hPhiLam0HighE","#phi vs #lambda_{0}^{2}, E > 2 GeV", nphibins,phimin,phimax, ssbins,ssmin,ssmax);
1265 fhPhiLam0HighE->SetYTitle("#lambda_{0}^{2}");
1266 fhPhiLam0HighE->SetXTitle("#phi");
1267 outputContainer->Add(fhPhiLam0HighE);
1268
1269 fhLam1Lam0LowE = new TH2F ("hLam1Lam0LowE","#lambda_{0}^{2} vs #lambda_{1}^{2} in cluster of E < 2 GeV", ssbins,ssmin,ssmax, ssbins,ssmin,ssmax);
1270 fhLam1Lam0LowE->SetYTitle("#lambda_{0}^{2}");
1271 fhLam1Lam0LowE->SetXTitle("#lambda_{1}^{2}");
1272 outputContainer->Add(fhLam1Lam0LowE);
1273
1274 fhLam1Lam0HighE = new TH2F ("hLam1Lam0HighE","#lambda_{0}^{2} vs #lambda_{1}^{2} in cluster of E > 2 GeV", ssbins,ssmin,ssmax, ssbins,ssmin,ssmax);
1275 fhLam1Lam0HighE->SetYTitle("#lambda_{0}^{2}");
1276 fhLam1Lam0HighE->SetXTitle("#lambda_{1}^{2}");
1277 outputContainer->Add(fhLam1Lam0HighE);
1278
1279 fhLam0DispLowE = new TH2F ("hLam0DispLowE","#lambda_{0}^{2} vs dispersion^{2} in cluster of E < 2 GeV", ssbins,ssmin,ssmax, ssbins,ssmin,ssmax);
1280 fhLam0DispLowE->SetXTitle("#lambda_{0}^{2}");
1281 fhLam0DispLowE->SetYTitle("D^{2}");
1282 outputContainer->Add(fhLam0DispLowE);
1283
1284 fhLam0DispHighE = new TH2F ("hLam0DispHighE","#lambda_{0}^{2} vs dispersion^{2} in cluster of E > 2 GeV", ssbins,ssmin,ssmax, ssbins,ssmin,ssmax);
1285 fhLam0DispHighE->SetXTitle("#lambda_{0}^{2}");
1286 fhLam0DispHighE->SetYTitle("D^{2}");
1287 outputContainer->Add(fhLam0DispHighE);
1288
1289 fhDispLam1LowE = new TH2F ("hDispLam1LowE","Dispersion^{2} vs #lambda_{1}^{2} in cluster of E < 2 GeV", ssbins,ssmin,ssmax, ssbins,ssmin,ssmax);
1290 fhDispLam1LowE->SetXTitle("D^{2}");
1291 fhDispLam1LowE->SetYTitle("#lambda_{1}^{2}");
1292 outputContainer->Add(fhDispLam1LowE);
1293
1294 fhDispLam1HighE = new TH2F ("hDispLam1HighE","Dispersion^{2} vs #lambda_{1^{2}} in cluster of E > 2 GeV", ssbins,ssmin,ssmax, ssbins,ssmin,ssmax);
1295 fhDispLam1HighE->SetXTitle("D^{2}");
1296 fhDispLam1HighE->SetYTitle("#lambda_{1}^{2}");
1297 outputContainer->Add(fhDispLam1HighE);
1298
1299
1300
1301 fhEtadLam0LowE = new TH2F ("hEtadLam0LowE","#eta vs #lambda_{0}^{2}/N_{cells}^{2}, E < 2 GeV", netabins,etamin,etamax, ssbins,ssmin,ssmax/50);
1302 fhEtadLam0LowE->SetYTitle("d#lambda_{0}^{2}");
1303 fhEtadLam0LowE->SetXTitle("#eta");
1304 outputContainer->Add(fhEtadLam0LowE);
1305
1306 fhPhidLam0LowE = new TH2F ("hPhidLam0LowE","#phi vs #lambda_{0}^{2}/N_{cells}^{2}, E < 2 GeV", nphibins,phimin,phimax, ssbins,ssmin,ssmax/50);
1307 fhPhidLam0LowE->SetYTitle("d#lambda_{0}^{2}");
1308 fhPhidLam0LowE->SetXTitle("#phi");
1309 outputContainer->Add(fhPhidLam0LowE);
1310
1311 fhEtadLam0HighE = new TH2F ("hEtadLam0HighE","#eta vs #lambda_{0}^{2}/N_{cells}^{2}", netabins,etamin,etamax, ssbins,ssmin,ssmax/50);
1312 fhEtadLam0HighE->SetYTitle("d#lambda_{0}^{2}");
1313 fhEtadLam0HighE->SetXTitle("#eta");
1314 outputContainer->Add(fhEtadLam0HighE);
1315
1316 fhPhidLam0HighE = new TH2F ("hPhidLam0HighE","#phi vs #lambda_{0}^{2}/N_{cells}^{2}, E > 2 GeV", nphibins,phimin,phimax, ssbins,ssmin,ssmax/50);
1317 fhPhidLam0HighE->SetYTitle("d#lambda_{0}^{2}");
1318 fhPhidLam0HighE->SetXTitle("#phi");
1319 outputContainer->Add(fhPhidLam0HighE);
1320
1321 fhdLam1dLam0LowE = new TH2F ("hdLam1dLam0LowE","#lambda_{0}^{2}/N_{cells}^{2} vs #lambda_{1}^{2}/N_{cells}^{2} in cluster of E < 2 GeV", ssbins,ssmin,ssmax/50, ssbins,ssmin,ssmax/50);
1322 fhdLam1dLam0LowE->SetYTitle("d#lambda_{0}^{2}");
1323 fhdLam1dLam0LowE->SetXTitle("d#lambda_{1}^{2}");
1324 outputContainer->Add(fhdLam1dLam0LowE);
1325
1326 fhdLam1dLam0HighE = new TH2F ("hdLam1dLam0HighE","#lambda_{0}^{2}/N_{cells}^{2}vs #lambda_{1}^{2}/N_{cells}^{2} in cluster of E > 2 GeV", ssbins,ssmin,ssmax/50, ssbins,ssmin,ssmax/50);
1327 fhdLam1dLam0HighE->SetYTitle("d#lambda_{0}^{2}");
1328 fhdLam1dLam0HighE->SetXTitle("d#lambda_{1}^{2}");
1329 outputContainer->Add(fhdLam1dLam0HighE);
1330
1331 fhdLam0dDispLowE = new TH2F ("hdLam0dDispLowE","#lambda_{0}^{2}/N_{cells}^{2} vs dispersion^{2}/N_{cells}^{2} in cluster of E < 2 GeV", ssbins,ssmin,ssmax/50, ssbins,ssmin,ssmax/50);
1332 fhdLam0dDispLowE->SetXTitle("d#lambda_{0}^{2}");
1333 fhdLam0dDispLowE->SetYTitle("dD^{2}");
1334 outputContainer->Add(fhdLam0dDispLowE);
1335
1336 fhdLam0dDispHighE = new TH2F ("hdLam0dDispHighE","#lambda_{0}^{2}/N_{cells}^{2} vs dispersion^{2}/N_{cells}^{2} in cluster of E > 2 GeV", ssbins,ssmin,ssmax/50, ssbins,ssmin,ssmax/50);
1337 fhdLam0dDispHighE->SetXTitle("d#lambda_{0}^{2}");
1338 fhdLam0dDispHighE->SetYTitle("dD^{2}");
1339 outputContainer->Add(fhdLam0dDispHighE);
1340
1341 fhdDispdLam1LowE = new TH2F ("hdDispddLam1LowE","Dispersion^{2}/N_{cells}^{2} vs #lambda_{1}^{2}/N_{cells}^{2} in cluster of E < 2 GeV", ssbins,ssmin,ssmax/50, ssbins,ssmin,ssmax/50);
1342 fhdDispdLam1LowE->SetXTitle("dD^{2}");
1343 fhdDispdLam1LowE->SetYTitle("d#lambda_{1}^{2}");
1344 outputContainer->Add(fhdDispdLam1LowE);
1345
1346 fhdDispdLam1HighE = new TH2F ("hdDispdLam1HighE","Dispersion^{2}/N_{cells}^{2} vs #lambda_{1^{2}}/N_{cells}^{2} in cluster of E > 2 GeV", ssbins,ssmin,ssmax/50, ssbins,ssmin,ssmax/50);
1347 fhdDispdLam1HighE->SetXTitle("dD^{2}");
1348 fhdDispdLam1HighE->SetYTitle("d#lambda_{1}^{2}");
1349 outputContainer->Add(fhdDispdLam1HighE);
3d5d5078 1350
521636d2 1351
1352 } // Shower shape
1353
6175da48 1354
477d6cee 1355 if(IsDataMC()){
123fc3bd 1356 fhDeltaE = new TH1F ("hDeltaE","MC - Reco E ", 200,-50,50);
1357 fhDeltaE->SetXTitle("#Delta E (GeV)");
1358 outputContainer->Add(fhDeltaE);
1359
1360 fhDeltaPt = new TH1F ("hDeltaPt","MC - Reco p_{T} ", 200,-50,50);
1361 fhDeltaPt->SetXTitle("#Delta p_{T} (GeV/c)");
1362 outputContainer->Add(fhDeltaPt);
1363
1364 fhRatioE = new TH1F ("hRatioE","Reco/MC E ", 200,0,2);
1365 fhRatioE->SetXTitle("E_{reco}/E_{gen}");
1366 outputContainer->Add(fhRatioE);
477d6cee 1367
123fc3bd 1368 fhRatioPt = new TH1F ("hRatioPt","Reco/MC p_{T} ", 200,0,2);
1369 fhRatioPt->SetXTitle("p_{T, reco}/p_{T, gen}");
1370 outputContainer->Add(fhRatioPt);
1371
1372 fh2E = new TH2F ("h2E","E distribution, reconstructed vs generated", nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
41e886c8 1373 fh2E->SetXTitle("E_{rec} (GeV)");
1374 fh2E->SetYTitle("E_{gen} (GeV)");
123fc3bd 1375 outputContainer->Add(fh2E);
1376
1377 fh2Pt = new TH2F ("h2Pt","p_T distribution, reconstructed vs generated", nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
41e886c8 1378 fh2Pt->SetXTitle("p_{T,rec} (GeV/c)");
1379 fh2Pt->SetYTitle("p_{T,gen} (GeV/c)");
123fc3bd 1380 outputContainer->Add(fh2Pt);
1381
3d5d5078 1382 TString ptype[] = { "#gamma","#gamma_{prompt}","#gamma_{fragmentation}","#gamma_{ISR}",
1383 "#gamma_{#pi decay}","#gamma_{decay}","hadron?",
1384 "#pi^{0}","#gamma->e^{#pm}","e^{#pm}",
1385 "Anti-N","Anti-P","String" } ;
1386
1387 TString pname[] = { "Photon","PhotonPrompt","PhotonFragmentation","PhotonISR",
1388 "PhotonPi0Decay", "PhotonOtherDecay", "Hadron",
1389 "Pi0","Conversion","Electron",
1390 "AntiNeutron","AntiProton","String" } ;
521636d2 1391
3d5d5078 1392 for(Int_t i = 0; i < 13; i++){
521636d2 1393
3d5d5078 1394 fhMCE[i] = new TH1F(Form("hE_MC%s",pname[i].Data()),
521636d2 1395 Form("cluster from %s : E ",ptype[i].Data()),
1396 nptbins,ptmin,ptmax);
3d5d5078 1397 fhMCE[i]->SetXTitle("E (GeV)");
1398 outputContainer->Add(fhMCE[i]) ;
521636d2 1399
1400 fhPtMC[i] = new TH1F(Form("hPt_MC%s",pname[i].Data()),
1401 Form("cluster from %s : p_{T} ",ptype[i].Data()),
1402 nptbins,ptmin,ptmax);
1403 fhPtMC[i]->SetXTitle("p_{T} (GeV/c)");
1404 outputContainer->Add(fhPtMC[i]) ;
1405
1406 fhEtaMC[i] = new TH2F(Form("hEta_MC%s",pname[i].Data()),
1407 Form("cluster from %s : #eta ",ptype[i].Data()),
1408 nptbins,ptmin,ptmax,netabins,etamin,etamax);
1409 fhEtaMC[i]->SetYTitle("#eta");
1410 fhEtaMC[i]->SetXTitle("E (GeV)");
1411 outputContainer->Add(fhEtaMC[i]) ;
1412
1413 fhPhiMC[i] = new TH2F(Form("hPhi_MC%s",pname[i].Data()),
1414 Form("cluster from %s : #phi ",ptype[i].Data()),
1415 nptbins,ptmin,ptmax,nphibins,phimin,phimax);
1416 fhPhiMC[i]->SetYTitle("#phi (rad)");
1417 fhPhiMC[i]->SetXTitle("E (GeV)");
1418 outputContainer->Add(fhPhiMC[i]) ;
1419
1420 }
3d5d5078 1421
1422 for(Int_t i = 0; i < 7; i++){
1423 fhEPrimMC[i] = new TH1F(Form("hEPrim_MC%s",pname[i].Data()),
1424 Form("primary photon %s : E ",ptype[i].Data()),
1425 nptbins,ptmin,ptmax);
1426 fhEPrimMC[i]->SetXTitle("E (GeV)");
1427 outputContainer->Add(fhEPrimMC[i]) ;
1428
1429 fhPtPrimMC[i] = new TH1F(Form("hPtPrim_MC%s",pname[i].Data()),
1430 Form("primary photon %s : p_{T} ",ptype[i].Data()),
1431 nptbins,ptmin,ptmax);
1432 fhPtPrimMC[i]->SetXTitle("p_{T} (GeV/c)");
1433 outputContainer->Add(fhPtPrimMC[i]) ;
1434
1435 fhYPrimMC[i] = new TH2F(Form("hYPrim_MC%s",pname[i].Data()),
1436 Form("primary photon %s : Rapidity ",ptype[i].Data()),
1437 nptbins,ptmin,ptmax,800,-8,8);
1438 fhYPrimMC[i]->SetYTitle("Rapidity");
1439 fhYPrimMC[i]->SetXTitle("E (GeV)");
1440 outputContainer->Add(fhYPrimMC[i]) ;
1441
1442 fhPhiPrimMC[i] = new TH2F(Form("hPhiPrim_MC%s",pname[i].Data()),
1443 Form("primary photon %s : #phi ",ptype[i].Data()),
1444 nptbins,ptmin,ptmax,nphibins,phimin,phimax);
1445 fhPhiPrimMC[i]->SetYTitle("#phi (rad)");
1446 fhPhiPrimMC[i]->SetXTitle("E (GeV)");
1447 outputContainer->Add(fhPhiPrimMC[i]) ;
1448
1449
1450 fhEPrimMCAcc[i] = new TH1F(Form("hEPrimAcc_MC%s",pname[i].Data()),
1451 Form("primary photon %s in acceptance: E ",ptype[i].Data()),
1452 nptbins,ptmin,ptmax);
1453 fhEPrimMCAcc[i]->SetXTitle("E (GeV)");
1454 outputContainer->Add(fhEPrimMCAcc[i]) ;
1455
1456 fhPtPrimMCAcc[i] = new TH1F(Form("hPtPrimAcc_MC%s",pname[i].Data()),
1457 Form("primary photon %s in acceptance: p_{T} ",ptype[i].Data()),
1458 nptbins,ptmin,ptmax);
1459 fhPtPrimMCAcc[i]->SetXTitle("p_{T} (GeV/c)");
1460 outputContainer->Add(fhPtPrimMCAcc[i]) ;
1461
1462 fhYPrimMCAcc[i] = new TH2F(Form("hYPrimAcc_MC%s",pname[i].Data()),
1463 Form("primary photon %s in acceptance: Rapidity ",ptype[i].Data()),
1464 nptbins,ptmin,ptmax,100,-1,1);
1465 fhYPrimMCAcc[i]->SetYTitle("Rapidity");
1466 fhYPrimMCAcc[i]->SetXTitle("E (GeV)");
1467 outputContainer->Add(fhYPrimMCAcc[i]) ;
1468
1469 fhPhiPrimMCAcc[i] = new TH2F(Form("hPhiPrimAcc_MC%s",pname[i].Data()),
1470 Form("primary photon %s in acceptance: #phi ",ptype[i].Data()),
1471 nptbins,ptmin,ptmax,nphibins,phimin,phimax);
1472 fhPhiPrimMCAcc[i]->SetYTitle("#phi (rad)");
1473 fhPhiPrimMCAcc[i]->SetXTitle("E (GeV)");
1474 outputContainer->Add(fhPhiPrimMCAcc[i]) ;
1475
1476 }
1477
20218aea 1478 if(fCheckConversion){
1479 fhPtConversionTagged = new TH1F("hPtMCConversionTagged","Number of converted #gamma over calorimeter, tagged as converted",nptbins,ptmin,ptmax);
1480 fhPtConversionTagged->SetYTitle("N");
1481 fhPtConversionTagged->SetXTitle("p_{T #gamma}(GeV/c)");
1482 outputContainer->Add(fhPtConversionTagged) ;
1483
1484
1485 fhPtAntiNeutronTagged = new TH1F("hPtMCAntiNeutronTagged","Number of AntiNeutron id as Photon over calorimeter, tagged as converted",nptbins,ptmin,ptmax);
1486 fhPtAntiNeutronTagged->SetYTitle("N");
1487 fhPtAntiNeutronTagged->SetXTitle("p_{T #gamma}(GeV/c)");
1488 outputContainer->Add(fhPtAntiNeutronTagged) ;
1489
1490 fhPtAntiProtonTagged = new TH1F("hPtMCAntiProtonTagged","Number of AntiProton id as Photon over calorimeter, tagged as converted",nptbins,ptmin,ptmax);
1491 fhPtAntiProtonTagged->SetYTitle("N");
1492 fhPtAntiProtonTagged->SetXTitle("p_{T #gamma}(GeV/c)");
1493 outputContainer->Add(fhPtAntiProtonTagged) ;
1494
1495 fhPtUnknownTagged = new TH1F("hPtMCUnknownTagged","Number of Unknown id as Photon over calorimeter, tagged as converted",nptbins,ptmin,ptmax);
1496 fhPtUnknownTagged->SetYTitle("N");
1497 fhPtUnknownTagged->SetXTitle("p_{T #gamma}(GeV/c)");
1498 outputContainer->Add(fhPtUnknownTagged) ;
1499
1500 fhConvDeltaEtaMCConversion = new TH2F
1501 ("hConvDeltaEtaMCConversion","#Delta #eta of selected conversion pairs from real conversions",100,0,fMassCut,netabins,-0.5,0.5);
1502 fhConvDeltaEtaMCConversion->SetYTitle("#Delta #eta");
1503 fhConvDeltaEtaMCConversion->SetXTitle("Pair Mass (GeV/c^2)");
1504 outputContainer->Add(fhConvDeltaEtaMCConversion) ;
1505
1506 fhConvDeltaPhiMCConversion = new TH2F
1507 ("hConvDeltaPhiMCConversion","#Delta #phi of selected conversion pairs from real conversions",100,0,fMassCut,nphibins,-0.5,0.5);
1508 fhConvDeltaPhiMCConversion->SetYTitle("#Delta #phi");
1509 fhConvDeltaPhiMCConversion->SetXTitle("Pair Mass (GeV/c^2)");
1510 outputContainer->Add(fhConvDeltaPhiMCConversion) ;
1511
1512 fhConvDeltaEtaPhiMCConversion = new TH2F
1513 ("hConvDeltaEtaPhiMCConversion","#Delta #eta vs #Delta #phi of selected conversion pairs, from real conversions",netabins,-0.5,0.5,nphibins,-0.5,0.5);
1514 fhConvDeltaEtaPhiMCConversion->SetYTitle("#Delta #phi");
1515 fhConvDeltaEtaPhiMCConversion->SetXTitle("#Delta #eta");
1516 outputContainer->Add(fhConvDeltaEtaPhiMCConversion) ;
1517
1518 fhConvAsymMCConversion = new TH2F
1519 ("hConvAsymMCConversion","Asymmetry of selected conversion pairs from real conversions",100,0,fMassCut,100,0,1);
1520 fhConvAsymMCConversion->SetYTitle("Asymmetry");
1521 fhConvAsymMCConversion->SetXTitle("Pair Mass (GeV/c^2)");
1522 outputContainer->Add(fhConvAsymMCConversion) ;
1523
1524 fhConvPtMCConversion = new TH2F
1525 ("hConvPtMCConversion","p_{T} of selected conversion pairs from real conversions",100,0,fMassCut,100,0.,10.);
1526 fhConvPtMCConversion->SetYTitle("Pair p_{T} (GeV/c)");
1527 fhConvPtMCConversion->SetXTitle("Pair Mass (GeV/c^2)");
1528 outputContainer->Add(fhConvPtMCConversion) ;
1529
1530 fhConvDispersionMCConversion = new TH2F
1531 ("hConvDispersionMCConversion","p_{T} of selected conversion pairs from real conversions",100,0.,1.,100,0.,1.);
1532 fhConvDispersionMCConversion->SetYTitle("Dispersion cluster 1");
1533 fhConvDispersionMCConversion->SetXTitle("Dispersion cluster 2");
1534 outputContainer->Add(fhConvDispersionMCConversion) ;
1535
1536 fhConvM02MCConversion = new TH2F
1537 ("hConvM02MCConversion","p_{T} of selected conversion pairs from string",100,0.,1.,100,0.,1.);
1538 fhConvM02MCConversion->SetYTitle("M02 cluster 1");
1539 fhConvM02MCConversion->SetXTitle("M02 cluster 2");
1540 outputContainer->Add(fhConvM02MCConversion) ;
1541
1542 fhConvDeltaEtaMCAntiNeutron = new TH2F
1543 ("hConvDeltaEtaMCAntiNeutron","#Delta #eta of selected conversion pairs from anti-neutrons",100,0,fMassCut,netabins,-0.5,0.5);
1544 fhConvDeltaEtaMCAntiNeutron->SetYTitle("#Delta #eta");
1545 fhConvDeltaEtaMCAntiNeutron->SetXTitle("Pair Mass (GeV/c^2)");
1546 outputContainer->Add(fhConvDeltaEtaMCAntiNeutron) ;
1547
1548 fhConvDeltaPhiMCAntiNeutron = new TH2F
1549 ("hConvDeltaPhiMCAntiNeutron","#Delta #phi of selected conversion pairs from anti-neutrons",100,0,fMassCut,nphibins,-0.5,0.5);
1550 fhConvDeltaPhiMCAntiNeutron->SetYTitle("#Delta #phi");
1551 fhConvDeltaPhiMCAntiNeutron->SetXTitle("Pair Mass (GeV/c^2)");
1552 outputContainer->Add(fhConvDeltaPhiMCAntiNeutron) ;
1553
1554 fhConvDeltaEtaPhiMCAntiNeutron = new TH2F
1555 ("hConvDeltaEtaPhiMCAntiNeutron","#Delta #eta vs #Delta #phi of selected conversion pairs from anti-neutrons",netabins,-0.5,0.5,nphibins,-0.5,0.5);
1556 fhConvDeltaEtaPhiMCAntiNeutron->SetYTitle("#Delta #phi");
1557 fhConvDeltaEtaPhiMCAntiNeutron->SetXTitle("#Delta #eta");
1558 outputContainer->Add(fhConvDeltaEtaPhiMCAntiNeutron) ;
1559
1560 fhConvAsymMCAntiNeutron = new TH2F
1561 ("hConvAsymMCAntiNeutron","Asymmetry of selected conversion pairs from anti-neutrons",100,0,fMassCut,100,0,1);
1562 fhConvAsymMCAntiNeutron->SetYTitle("Asymmetry");
1563 fhConvAsymMCAntiNeutron->SetXTitle("Pair Mass (GeV/c^2)");
1564 outputContainer->Add(fhConvAsymMCAntiNeutron) ;
1565
1566 fhConvPtMCAntiNeutron = new TH2F
1567 ("hConvPtMCAntiNeutron","p_{T} of selected conversion pairs from anti-neutrons",100,0,fMassCut,100,0.,10.);
1568 fhConvPtMCAntiNeutron->SetYTitle("Pair p_{T} (GeV/c)");
1569 fhConvPtMCAntiNeutron->SetXTitle("Pair Mass (GeV/c^2)");
1570 outputContainer->Add(fhConvPtMCAntiNeutron) ;
1571
1572 fhConvDispersionMCAntiNeutron = new TH2F
1573 ("hConvDispersionMCAntiNeutron","p_{T} of selected conversion pairs from anti-neutrons",100,0.,1.,100,0.,1.);
1574 fhConvDispersionMCAntiNeutron->SetYTitle("Dispersion cluster 1");
1575 fhConvDispersionMCAntiNeutron->SetXTitle("Dispersion cluster 2");
1576 outputContainer->Add(fhConvDispersionMCAntiNeutron) ;
1577
1578 fhConvM02MCAntiNeutron = new TH2F
1579 ("hConvM02MCAntiNeutron","p_{T} of selected conversion pairs from string",100,0.,1.,100,0.,1.);
1580 fhConvM02MCAntiNeutron->SetYTitle("M02 cluster 1");
1581 fhConvM02MCAntiNeutron->SetXTitle("M02 cluster 2");
1582 outputContainer->Add(fhConvM02MCAntiNeutron) ;
1583
1584 fhConvDeltaEtaMCAntiProton = new TH2F
1585 ("hConvDeltaEtaMCAntiProton","#Delta #eta of selected conversion pairs from anti-protons",100,0,fMassCut,netabins,-0.5,0.5);
1586 fhConvDeltaEtaMCAntiProton->SetYTitle("#Delta #eta");
1587 fhConvDeltaEtaMCAntiProton->SetXTitle("Pair Mass (GeV/c^2)");
1588 outputContainer->Add(fhConvDeltaEtaMCAntiProton) ;
1589
1590 fhConvDeltaPhiMCAntiProton = new TH2F
1591 ("hConvDeltaPhiMCAntiProton","#Delta #phi of selected conversion pairs from anti-protons",100,0,fMassCut,nphibins,-0.5,0.5);
1592 fhConvDeltaPhiMCAntiProton->SetYTitle("#Delta #phi");
1593 fhConvDeltaPhiMCAntiProton->SetXTitle("Pair Mass (GeV/c^2)");
1594 outputContainer->Add(fhConvDeltaPhiMCAntiProton) ;
1595
1596 fhConvDeltaEtaPhiMCAntiProton = new TH2F
1597 ("hConvDeltaEtaPhiMCAntiProton","#Delta #eta vs #Delta #phi of selected conversion pairs from anti-protons",netabins,-0.5,0.5,nphibins,-0.5,0.5);
1598 fhConvDeltaEtaPhiMCAntiProton->SetYTitle("#Delta #phi");
1599 fhConvDeltaEtaPhiMCAntiProton->SetXTitle("#Delta #eta");
1600 outputContainer->Add(fhConvDeltaEtaPhiMCAntiProton) ;
1601
1602 fhConvAsymMCAntiProton = new TH2F
1603 ("hConvAsymMCAntiProton","Asymmetry of selected conversion pairs from anti-protons",100,0,fMassCut,100,0,1);
1604 fhConvAsymMCAntiProton->SetYTitle("Asymmetry");
1605 fhConvAsymMCAntiProton->SetXTitle("Pair Mass (GeV/c^2)");
1606 outputContainer->Add(fhConvAsymMCAntiProton) ;
1607
1608 fhConvPtMCAntiProton = new TH2F
1609 ("hConvPtMCAntiProton","p_{T} of selected conversion pairs from anti-protons",100,0,fMassCut,100,0.,10.);
1610 fhConvPtMCAntiProton->SetYTitle("Pair p_{T} (GeV/c)");
1611 fhConvPtMCAntiProton->SetXTitle("Pair Mass (GeV/c^2)");
1612 outputContainer->Add(fhConvPtMCAntiProton) ;
1613
1614 fhConvDispersionMCAntiProton = new TH2F
1615 ("hConvDispersionMCAntiProton","p_{T} of selected conversion pairs from anti-protons",100,0.,1.,100,0.,1.);
1616 fhConvDispersionMCAntiProton->SetYTitle("Dispersion cluster 1");
1617 fhConvDispersionMCAntiProton->SetXTitle("Dispersion cluster 2");
1618 outputContainer->Add(fhConvDispersionMCAntiProton) ;
1619
1620 fhConvM02MCAntiProton = new TH2F
1621 ("hConvM02MCAntiProton","p_{T} of selected conversion pairs from string",100,0.,1.,100,0.,1.);
1622 fhConvM02MCAntiProton->SetYTitle("M02 cluster 1");
1623 fhConvM02MCAntiProton->SetXTitle("M02 cluster 2");
1624 outputContainer->Add(fhConvM02MCAntiProton) ;
1625
1626 fhConvDeltaEtaMCString = new TH2F
1627 ("hConvDeltaEtaMCString","#Delta #eta of selected conversion pairs from string",100,0,fMassCut,netabins,-0.5,0.5);
1628 fhConvDeltaEtaMCString->SetYTitle("#Delta #eta");
1629 fhConvDeltaEtaMCString->SetXTitle("Pair Mass (GeV/c^2)");
1630 outputContainer->Add(fhConvDeltaEtaMCString) ;
1631
1632 fhConvDeltaPhiMCString = new TH2F
1633 ("hConvDeltaPhiMCString","#Delta #phi of selected conversion pairs from string",100,0,fMassCut,nphibins,-0.5,0.5);
1634 fhConvDeltaPhiMCString->SetYTitle("#Delta #phi");
1635 fhConvDeltaPhiMCString->SetXTitle("Pair Mass (GeV/c^2)");
1636 outputContainer->Add(fhConvDeltaPhiMCString) ;
1637
1638 fhConvDeltaEtaPhiMCString = new TH2F
1639 ("hConvDeltaEtaPhiMCString","#Delta #eta vs #Delta #phi of selected conversion pairs from string",netabins,-0.5,0.5,nphibins,-0.5,0.5);
1640 fhConvDeltaEtaPhiMCString->SetYTitle("#Delta #phi");
1641 fhConvDeltaEtaPhiMCString->SetXTitle("#Delta #eta");
1642 outputContainer->Add(fhConvDeltaEtaPhiMCString) ;
1643
1644 fhConvAsymMCString = new TH2F
1645 ("hConvAsymMCString","Asymmetry of selected conversion pairs from string",100,0,fMassCut,100,0,1);
1646 fhConvAsymMCString->SetYTitle("Asymmetry");
1647 fhConvAsymMCString->SetXTitle("Pair Mass (GeV/c^2)");
1648 outputContainer->Add(fhConvAsymMCString) ;
1649
1650 fhConvPtMCString = new TH2F
1651 ("hConvPtMCString","p_{T} of selected conversion pairs from string",100,0,fMassCut,100,0.,10.);
1652 fhConvPtMCString->SetYTitle("Pair p_{T} (GeV/c)");
1653 fhConvPtMCString->SetXTitle("Pair Mass (GeV/c^2)");
1654 outputContainer->Add(fhConvPtMCString) ;
1655
1656 fhConvDispersionMCString = new TH2F
1657 ("hConvDispersionMCString","p_{T} of selected conversion pairs from string",100,0.,1.,100,0.,1.);
1658 fhConvDispersionMCString->SetYTitle("Dispersion cluster 1");
1659 fhConvDispersionMCString->SetXTitle("Dispersion cluster 2");
1660 outputContainer->Add(fhConvDispersionMCString) ;
1661
1662 fhConvM02MCString = new TH2F
1663 ("hConvM02MCString","p_{T} of selected conversion pairs from string",100,0.,1.,100,0.,1.);
1664 fhConvM02MCString->SetYTitle("M02 cluster 1");
1665 fhConvM02MCString->SetXTitle("M02 cluster 2");
c4a7d28a 1666 outputContainer->Add(fhConvM02MCString) ;
1667
1668 fhConvDistMCConversion = new TH2F
1669 ("hConvDistMCConversion","calculated conversion distance vs real vertes for MC conversion",100,0.,5.,100,0.,5.);
1670 fhConvDistMCConversion->SetYTitle("distance");
1671 fhConvDistMCConversion->SetXTitle("vertex R");
1672 outputContainer->Add(fhConvDistMCConversion) ;
1673
1674 fhConvDistMCConversionCuts = new TH2F
1675 ("hConvDistMCConversionCuts","calculated conversion distance vs real vertes for MC conversion, deta < 0.05, m < 10 MeV, asym < 0.1",100,0.,5.,100,0.,5.);
1676 fhConvDistMCConversionCuts->SetYTitle("distance");
1677 fhConvDistMCConversionCuts->SetXTitle("vertex R");
1678 outputContainer->Add(fhConvDistMCConversionCuts) ;
20218aea 1679 }
6175da48 1680
521636d2 1681 if(fFillSSHistograms){
1682
3d5d5078 1683 TString ptypess[] = { "#gamma","hadron?","#pi^{0}","#eta","#gamma->e^{#pm}","e^{#pm}"} ;
1684
1685 TString pnamess[] = { "Photon","Hadron","Pi0","Eta","Conversion","Electron"} ;
1686
1687 for(Int_t i = 0; i < 6; i++){
521636d2 1688
3d5d5078 1689 fhMCELambda0[i] = new TH2F(Form("hELambda0_MC%s",pnamess[i].Data()),
1690 Form("cluster from %s : E vs #lambda_{0}^{2}",ptypess[i].Data()),
521636d2 1691 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
3d5d5078 1692 fhMCELambda0[i]->SetYTitle("#lambda_{0}^{2}");
1693 fhMCELambda0[i]->SetXTitle("E (GeV)");
1694 outputContainer->Add(fhMCELambda0[i]) ;
521636d2 1695
3d5d5078 1696 fhMCEdLambda0[i] = new TH2F(Form("hEdLambda0_MC%s",pnamess[i].Data()),
1697 Form("cluster from %s : E vs #lambda_{0}^{2}/N_{cells}^{2}",ptypess[i].Data()),
521636d2 1698 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax/50);
3d5d5078 1699 fhMCEdLambda0[i]->SetYTitle("d#lambda_{0}^{2}");
1700 fhMCEdLambda0[i]->SetXTitle("E (GeV)");
1701 outputContainer->Add(fhMCEdLambda0[i]) ;
521636d2 1702
3d5d5078 1703 fhMCELambda1[i] = new TH2F(Form("hELambda1_MC%s",pnamess[i].Data()),
1704 Form("cluster from %s : E vs #lambda_{1}^{2}",ptypess[i].Data()),
521636d2 1705 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
3d5d5078 1706 fhMCELambda1[i]->SetYTitle("#lambda_{1}^{2}");
1707 fhMCELambda1[i]->SetXTitle("E (GeV)");
1708 outputContainer->Add(fhMCELambda1[i]) ;
521636d2 1709
3d5d5078 1710 fhMCEdLambda1[i] = new TH2F(Form("hEdLambda1_MC%s",pnamess[i].Data()),
1711 Form("cluster from %s : E vs d#lambda_{1}^{2}/N_{cells}^{2}",ptypess[i].Data()),
521636d2 1712 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax/50);
3d5d5078 1713 fhMCEdLambda1[i]->SetYTitle("d#lambda_{1}^{2}");
1714 fhMCEdLambda1[i]->SetXTitle("E (GeV)");
1715 outputContainer->Add(fhMCEdLambda1[i]) ;
521636d2 1716
3d5d5078 1717 fhMCEDispersion[i] = new TH2F(Form("hEDispersion_MC%s",pnamess[i].Data()),
1718 Form("cluster from %s : E vs dispersion^{2}",ptypess[i].Data()),
521636d2 1719 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
3d5d5078 1720 fhMCEDispersion[i]->SetYTitle("D^{2}");
1721 fhMCEDispersion[i]->SetXTitle("E (GeV)");
1722 outputContainer->Add(fhMCEDispersion[i]) ;
521636d2 1723
3d5d5078 1724 fhMCEdDispersion[i] = new TH2F(Form("hEdDispersion_MC%s",pnamess[i].Data()),
1725 Form("cluster from %s : E vs dispersion^{2}/N_{cells}^{2}",ptypess[i].Data()),
521636d2 1726 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax/50);
3d5d5078 1727 fhMCEdDispersion[i]->SetYTitle("dD^{2}");
1728 fhMCEdDispersion[i]->SetXTitle("E (GeV)");
1729 outputContainer->Add(fhMCEdDispersion[i]) ;
1730
1731 }// loop
1732
1733 if(!GetReader()->IsEmbeddedClusterSelectionOn())
1734 {
1735 fhMCPhotonELambda0NoOverlap = new TH2F("hELambda0_MCPhoton_NoOverlap",
1736 "cluster from Photon : E vs #lambda_{0}^{2}",
1737 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
1738 fhMCPhotonELambda0NoOverlap->SetYTitle("#lambda_{0}^{2}");
1739 fhMCPhotonELambda0NoOverlap->SetXTitle("E (GeV)");
1740 outputContainer->Add(fhMCPhotonELambda0NoOverlap) ;
1741
1742 fhMCPhotonEdLambda0NoOverlap = new TH2F("hEdLambda0_MCPhoton_NoOverlap",
1743 "cluster from Photon : E vs #lambda_{0}^{2}/N_{cells}^{2}",
1744 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax/50);
1745 fhMCPhotonEdLambda0NoOverlap->SetYTitle("d#lambda_{0}^{2}");
1746 fhMCPhotonEdLambda0NoOverlap->SetXTitle("E (GeV)");
1747 outputContainer->Add(fhMCPhotonEdLambda0NoOverlap) ;
1748
1749 fhMCPhotonELambda0TwoOverlap = new TH2F("hELambda0_MCPhoton_TwoOverlap",
1750 "cluster from Photon : E vs #lambda_{0}^{2}",
1751 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
1752 fhMCPhotonELambda0TwoOverlap->SetYTitle("#lambda_{0}^{2}");
1753 fhMCPhotonELambda0TwoOverlap->SetXTitle("E (GeV)");
1754 outputContainer->Add(fhMCPhotonELambda0TwoOverlap) ;
1755
1756 fhMCPhotonEdLambda0TwoOverlap = new TH2F("hEdLambda0_MCPhoton_TwoOverlap",
1757 "cluster from Photon : E vs #lambda_{0}^{2}/N_{cells}^{2}",
1758 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax/50);
1759 fhMCPhotonEdLambda0TwoOverlap->SetYTitle("d#lambda_{0}^{2}");
1760 fhMCPhotonEdLambda0TwoOverlap->SetXTitle("E (GeV)");
1761 outputContainer->Add(fhMCPhotonEdLambda0TwoOverlap) ;
1762
1763 fhMCPhotonELambda0NOverlap = new TH2F("hELambda0_MCPhoton_NOverlap",
1764 "cluster from Photon : E vs #lambda_{0}^{2}",
1765 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
1766 fhMCPhotonELambda0NOverlap->SetYTitle("#lambda_{0}^{2}");
1767 fhMCPhotonELambda0NOverlap->SetXTitle("E (GeV)");
1768 outputContainer->Add(fhMCPhotonELambda0NOverlap) ;
521636d2 1769
3d5d5078 1770 fhMCPhotonEdLambda0NOverlap = new TH2F("hEdLambda0_MCPhoton_NOverlap",
1771 "cluster from Photon : E vs #lambda_{0}^{2}/N_{cells}^{2}",
1772 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax/50);
1773 fhMCPhotonEdLambda0NOverlap->SetYTitle("d#lambda_{0}^{2}");
1774 fhMCPhotonEdLambda0NOverlap->SetXTitle("E (GeV)");
1775 outputContainer->Add(fhMCPhotonEdLambda0NOverlap) ;
1776 } //No embedding
1777
1778 //Fill histograms to check shape of embedded clusters
1779 if(GetReader()->IsEmbeddedClusterSelectionOn())
1780 {
1781
1782 fhEmbeddedSignalFractionEnergy = new TH2F("hEmbeddedSignal_FractionEnergy",
1783 "Energy Fraction of embedded signal versus cluster energy",
1784 nptbins,ptmin,ptmax,100,0.,1.);
1785 fhEmbeddedSignalFractionEnergy->SetYTitle("Fraction");
1786 fhEmbeddedSignalFractionEnergy->SetXTitle("E (GeV)");
1787 outputContainer->Add(fhEmbeddedSignalFractionEnergy) ;
1788
1789 fhEmbedPhotonELambda0FullSignal = new TH2F("hELambda0_EmbedPhoton_FullSignal",
1790 "cluster from Photon embedded with more than 90% energy in cluster : E vs #lambda_{0}^{2}",
1791 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
1792 fhEmbedPhotonELambda0FullSignal->SetYTitle("#lambda_{0}^{2}");
1793 fhEmbedPhotonELambda0FullSignal->SetXTitle("E (GeV)");
1794 outputContainer->Add(fhEmbedPhotonELambda0FullSignal) ;
1795
1796 fhEmbedPhotonEdLambda0FullSignal = new TH2F("hEdLambda0_EmbedPhoton_FullSignal",
1797 "cluster from Photon with more than 90% energy in cluster: E vs d#lambda_{0}^{2}",
1798 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
1799 fhEmbedPhotonEdLambda0FullSignal->SetYTitle("#lambda_{0}^{2}");
1800 fhEmbedPhotonEdLambda0FullSignal->SetXTitle("E (GeV)");
1801 outputContainer->Add(fhEmbedPhotonEdLambda0FullSignal) ;
1802
1803
1804 fhEmbedPhotonELambda0MostlySignal = new TH2F("hELambda0_EmbedPhoton_MostlySignal",
1805 "cluster from Photon embedded with 50% to 90% energy in cluster : E vs #lambda_{0}^{2}",
1806 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
1807 fhEmbedPhotonELambda0MostlySignal->SetYTitle("#lambda_{0}^{2}");
1808 fhEmbedPhotonELambda0MostlySignal->SetXTitle("E (GeV)");
1809 outputContainer->Add(fhEmbedPhotonELambda0MostlySignal) ;
1810
1811 fhEmbedPhotonEdLambda0MostlySignal = new TH2F("hEdLambda0_EmbedPhoton_MostlySignal",
1812 "cluster from Photon with 50% to 90% energy in cluster: E vs d#lambda_{0}^{2}",
1813 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
1814 fhEmbedPhotonEdLambda0MostlySignal->SetYTitle("#lambda_{0}^{2}");
1815 fhEmbedPhotonEdLambda0MostlySignal->SetXTitle("E (GeV)");
1816 outputContainer->Add(fhEmbedPhotonEdLambda0MostlySignal) ;
1817
1818
1819 fhEmbedPhotonELambda0MostlyBkg = new TH2F("hELambda0_EmbedPhoton_MostlyBkg",
1820 "cluster from Photon embedded with 10% to 50% energy in cluster : E vs #lambda_{0}^{2}",
1821 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
1822 fhEmbedPhotonELambda0MostlyBkg->SetYTitle("#lambda_{0}^{2}");
1823 fhEmbedPhotonELambda0MostlyBkg->SetXTitle("E (GeV)");
1824 outputContainer->Add(fhEmbedPhotonELambda0MostlyBkg) ;
1825
1826 fhEmbedPhotonEdLambda0MostlyBkg = new TH2F("hEdLambda0_EmbedPhoton_MostlyBkg",
1827 "cluster from Photon with 10% to 50% energy in cluster: E vs d#lambda_{0}^{2}",
1828 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
1829 fhEmbedPhotonEdLambda0MostlyBkg->SetYTitle("#lambda_{0}^{2}");
1830 fhEmbedPhotonEdLambda0MostlyBkg->SetXTitle("E (GeV)");
1831 outputContainer->Add(fhEmbedPhotonEdLambda0MostlyBkg) ;
1832
1833
1834 fhEmbedPhotonELambda0FullBkg = new TH2F("hELambda0_EmbedPhoton_FullBkg",
1835 "cluster from Photonm embedded with 0% to 10% energy in cluster : E vs #lambda_{0}^{2}",
1836 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
1837 fhEmbedPhotonELambda0FullBkg->SetYTitle("#lambda_{0}^{2}");
1838 fhEmbedPhotonELambda0FullBkg->SetXTitle("E (GeV)");
1839 outputContainer->Add(fhEmbedPhotonELambda0FullBkg) ;
1840
1841 fhEmbedPhotonEdLambda0FullBkg = new TH2F("hEdLambda0_EmbedPhoton_FullBkg",
1842 "cluster from Photon with 0% to 10% energy in cluster: E vs d#lambda_{0}^{2}",
1843 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
1844 fhEmbedPhotonEdLambda0FullBkg->SetYTitle("#lambda_{0}^{2}");
1845 fhEmbedPhotonEdLambda0FullBkg->SetXTitle("E (GeV)");
1846 outputContainer->Add(fhEmbedPhotonEdLambda0FullBkg) ;
1847
1848
1849 fhEmbedPi0ELambda0FullSignal = new TH2F("hELambda0_EmbedPi0_FullSignal",
1850 "cluster from Pi0 embedded with more than 90% energy in cluster : E vs #lambda_{0}^{2}",
1851 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
1852 fhEmbedPi0ELambda0FullSignal->SetYTitle("#lambda_{0}^{2}");
1853 fhEmbedPi0ELambda0FullSignal->SetXTitle("E (GeV)");
1854 outputContainer->Add(fhEmbedPi0ELambda0FullSignal) ;
1855
1856 fhEmbedPi0EdLambda0FullSignal = new TH2F("hEdLambda0_EmbedPi0_FullSignal",
1857 "cluster from Pi0 with more than 90% energy in cluster: E vs d#lambda_{0}^{2}",
1858 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
1859 fhEmbedPi0EdLambda0FullSignal->SetYTitle("#lambda_{0}^{2}");
1860 fhEmbedPi0EdLambda0FullSignal->SetXTitle("E (GeV)");
1861 outputContainer->Add(fhEmbedPi0EdLambda0FullSignal) ;
1862
1863
1864 fhEmbedPi0ELambda0MostlySignal = new TH2F("hELambda0_EmbedPi0_MostlySignal",
1865 "cluster from Pi0 embedded with 50% to 90% energy in cluster : E vs #lambda_{0}^{2}",
1866 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
1867 fhEmbedPi0ELambda0MostlySignal->SetYTitle("#lambda_{0}^{2}");
1868 fhEmbedPi0ELambda0MostlySignal->SetXTitle("E (GeV)");
1869 outputContainer->Add(fhEmbedPi0ELambda0MostlySignal) ;
1870
1871 fhEmbedPi0EdLambda0MostlySignal = new TH2F("hEdLambda0_EmbedPi0_MostlySignal",
1872 "cluster from Pi0 with 50% to 90% energy in cluster: E vs d#lambda_{0}^{2}",
1873 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
1874 fhEmbedPi0EdLambda0MostlySignal->SetYTitle("#lambda_{0}^{2}");
1875 fhEmbedPi0EdLambda0MostlySignal->SetXTitle("E (GeV)");
1876 outputContainer->Add(fhEmbedPi0EdLambda0MostlySignal) ;
1877
1878
1879 fhEmbedPi0ELambda0MostlyBkg = new TH2F("hELambda0_EmbedPi0_MostlyBkg",
1880 "cluster from Pi0 embedded with 10% to 50% energy in cluster : E vs #lambda_{0}^{2}",
1881 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
1882 fhEmbedPi0ELambda0MostlyBkg->SetYTitle("#lambda_{0}^{2}");
1883 fhEmbedPi0ELambda0MostlyBkg->SetXTitle("E (GeV)");
1884 outputContainer->Add(fhEmbedPi0ELambda0MostlyBkg) ;
1885
1886 fhEmbedPi0EdLambda0MostlyBkg = new TH2F("hEdLambda0_EmbedPi0_MostlyBkg",
1887 "cluster from Pi0 with 10% to 50% energy in cluster: E vs d#lambda_{0}^{2}",
1888 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
1889 fhEmbedPi0EdLambda0MostlyBkg->SetYTitle("#lambda_{0}^{2}");
1890 fhEmbedPi0EdLambda0MostlyBkg->SetXTitle("E (GeV)");
1891 outputContainer->Add(fhEmbedPi0EdLambda0MostlyBkg) ;
1892
1893
1894 fhEmbedPi0ELambda0FullBkg = new TH2F("hELambda0_EmbedPi0_FullBkg",
1895 "cluster from Pi0 embedded with 0% to 10% energy in cluster : E vs #lambda_{0}^{2}",
1896 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
1897 fhEmbedPi0ELambda0FullBkg->SetYTitle("#lambda_{0}^{2}");
1898 fhEmbedPi0ELambda0FullBkg->SetXTitle("E (GeV)");
1899 outputContainer->Add(fhEmbedPi0ELambda0FullBkg) ;
1900
1901 fhEmbedPi0EdLambda0FullBkg = new TH2F("hEdLambda0_EmbedPi0_FullBkg",
1902 "cluster from Pi0 with 0% to 10% energy in cluster: E vs d#lambda_{0}^{2}",
1903 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
1904 fhEmbedPi0EdLambda0FullBkg->SetYTitle("#lambda_{0}^{2}");
1905 fhEmbedPi0EdLambda0FullBkg->SetXTitle("E (GeV)");
1906 outputContainer->Add(fhEmbedPi0EdLambda0FullBkg) ;
1907
1908 }// embedded histograms
1909
521636d2 1910
1911 }// Fill SS MC histograms
1912
477d6cee 1913 }//Histos with MC
0c1383b5 1914
d39cba7e 1915 //Store calo PID histograms
1916 if(fRejectTrackMatch){
1917 TList * caloPIDHistos = GetCaloPID()->GetCreateOutputObjects() ;
1918 for(Int_t i = 0; i < caloPIDHistos->GetEntries(); i++) {
1919 outputContainer->Add(caloPIDHistos->At(i)) ;
1920 }
1921 delete caloPIDHistos;
1922 }
1923
477d6cee 1924 return outputContainer ;
1925
1c5acb87 1926}
1927
6639984f 1928//____________________________________________________________________________
1929void AliAnaPhoton::Init()
1930{
1931
1932 //Init
1933 //Do some checks
1e86c71e 1934 if(fCalorimeter == "PHOS" && !GetReader()->IsPHOSSwitchedOn() && NewOutputAOD()){
591cc579 1935 printf("AliAnaPhoton::Init() - !!STOP: You want to use PHOS in analysis but it is not read!! \n!!Check the configuration file!!\n");
6639984f 1936 abort();
1937 }
1e86c71e 1938 else if(fCalorimeter == "EMCAL" && !GetReader()->IsEMCALSwitchedOn() && NewOutputAOD()){
591cc579 1939 printf("AliAnaPhoton::Init() - !!STOP: You want to use EMCAL in analysis but it is not read!! \n!!Check the configuration file!!\n");
6639984f 1940 abort();
1941 }
1942
1943}
1944
1c5acb87 1945//____________________________________________________________________________
1946void AliAnaPhoton::InitParameters()
1947{
1948
1949 //Initialize the parameters of the analysis.
a3aebfff 1950 AddToHistogramsName("AnaPhoton_");
521636d2 1951
6175da48 1952 fCalorimeter = "EMCAL" ;
1953 fMinDist = 2.;
1954 fMinDist2 = 4.;
1955 fMinDist3 = 5.;
1956 fMassCut = 0.03; //30 MeV
1e86c71e 1957
4cf55759 1958 fTimeCutMin = -1;
1959 fTimeCutMax = 9999999;
6175da48 1960 fNCellsCut = 0;
2ac125bf 1961
1e86c71e 1962 fRejectTrackMatch = kTRUE ;
1963 fCheckConversion = kFALSE;
20218aea 1964 fRemoveConvertedPair = kFALSE;
1e86c71e 1965 fAddConvertedPairsToAOD = kFALSE;
1966
1c5acb87 1967}
1968
1969//__________________________________________________________________
1970void AliAnaPhoton::MakeAnalysisFillAOD()
1971{
f8006433 1972 //Do photon analysis and fill aods
f37fa8d2 1973
6175da48 1974 //Get the vertex
5025c139 1975 Double_t v[3] = {0,0,0}; //vertex ;
1976 GetReader()->GetVertex(v);
f8006433 1977
f37fa8d2 1978 //Select the Calorimeter of the photon
c8fe2783 1979 TObjArray * pl = 0x0;
477d6cee 1980 if(fCalorimeter == "PHOS")
be518ab0 1981 pl = GetPHOSClusters();
477d6cee 1982 else if (fCalorimeter == "EMCAL")
be518ab0 1983 pl = GetEMCALClusters();
5ae09196 1984
1985 if(!pl) {
1986 Info("MakeAnalysisFillAOD","TObjArray with %s clusters is NULL!\n",fCalorimeter.Data());
1987 return;
1988 }
521636d2 1989
6175da48 1990 //Init arrays, variables, get number of clusters
1e86c71e 1991 TLorentzVector mom, mom2 ;
1992 Int_t nCaloClusters = pl->GetEntriesFast();
6175da48 1993 //List to be used in conversion analysis, to tag the cluster as candidate for conversion
20218aea 1994 Bool_t * indexConverted = 0x0;
1995 if(fCheckConversion){
1996 indexConverted = new Bool_t[nCaloClusters];
1997 for (Int_t i = 0; i < nCaloClusters; i++)
1998 indexConverted[i] = kFALSE;
1999 }
2000
6175da48 2001 if(GetDebug() > 0) printf("AliAnaPhoton::MakeAnalysisFillAOD() - input %s cluster entries %d\n", fCalorimeter.Data(), nCaloClusters);
521636d2 2002
6175da48 2003 //----------------------------------------------------
2004 // Fill AOD with PHOS/EMCAL AliAODPWG4Particle objects
2005 //----------------------------------------------------
2006 // Loop on clusters
1e86c71e 2007 for(Int_t icalo = 0; icalo < nCaloClusters; icalo++){
2008
0ae57829 2009 AliVCluster * calo = (AliVCluster*) (pl->At(icalo));
2010 //printf("calo %d, %f\n",icalo,calo->E());
521636d2 2011
f8006433 2012 //Get the index where the cluster comes, to retrieve the corresponding vertex
c8fe2783 2013 Int_t evtIndex = 0 ;
2014 if (GetMixedEvent()) {
2015 evtIndex=GetMixedEvent()->EventIndexForCaloCluster(calo->GetID()) ;
5025c139 2016 //Get the vertex and check it is not too large in z
96539743 2017 if(TMath::Abs(GetVertex(evtIndex)[2])> GetZvertexCut()) continue;
c8fe2783 2018 }
521636d2 2019
2020 //Cluster selection, not charged, with photon id and in fiducial cut
f8006433 2021 if(GetReader()->GetDataType() != AliCaloTrackReader::kMC){
2022 calo->GetMomentum(mom,GetVertex(evtIndex)) ;}//Assume that come from vertex in straight line
2023 else{
2024 Double_t vertex[]={0,0,0};
2025 calo->GetMomentum(mom,vertex) ;
2026 }
c8fe2783 2027
6175da48 2028 //--------------------------------------
2029 // Cluster selection
2030 //--------------------------------------
c4a7d28a 2031 if(!ClusterSelected(calo,mom)) continue;
6175da48 2032
2033 //----------------------------
2034 //Create AOD for analysis
2035 //----------------------------
2036 AliAODPWG4Particle aodph = AliAODPWG4Particle(mom);
2037
2038 //...............................................
2039 //Set the indeces of the original caloclusters (MC, ID), and calorimeter
2040 Int_t label = calo->GetLabel();
2041 aodph.SetLabel(label);
6175da48 2042 aodph.SetCaloLabel(calo->GetID(),-1);
2043 aodph.SetDetector(fCalorimeter);
c4a7d28a 2044 //printf("Index %d, Id %d, iaod %d\n",icalo, calo->GetID(),GetOutputAODBranch()->GetEntriesFast());
521636d2 2045
6175da48 2046 //...............................................
2047 //Set bad channel distance bit
c4a7d28a 2048 Double_t distBad=calo->GetDistanceToBadChannel() ; //Distance to bad channel
f37fa8d2 2049 if (distBad > fMinDist3) aodph.SetDistToBad(2) ;
477d6cee 2050 else if(distBad > fMinDist2) aodph.SetDistToBad(1) ;
f37fa8d2 2051 else aodph.SetDistToBad(0) ;
af7b3903 2052 //printf("DistBad %f Bit %d\n",distBad, aodph.DistToBad());
c8fe2783 2053
521636d2 2054 //--------------------------------------------------------------------------------------
2055 //Play with the MC stack if available
2056 //--------------------------------------------------------------------------------------
2057
2058 //Check origin of the candidates
2059 if(IsDataMC()){
2060 aodph.SetTag(GetMCAnalysisUtils()->CheckOrigin(calo->GetLabels(),calo->GetNLabels(),GetReader(), aodph.GetInputFileIndex()));
3d5d5078 2061
521636d2 2062 if(GetDebug() > 0)
2063 printf("AliAnaPhoton::MakeAnalysisFillAOD() - Origin of candidate, bit map %d\n",aodph.GetTag());
2064 }//Work with stack also
2065
2066 //--------------------------------------------------------------------------------------
2067 //Fill some shower shape histograms before PID is applied
2068 //--------------------------------------------------------------------------------------
2069
2070 FillShowerShapeHistograms(calo,aodph.GetTag());
6175da48 2071
2072 //-------------------------------------
f37fa8d2 2073 //PID selection or bit setting
6175da48 2074 //-------------------------------------
2075 // MC
477d6cee 2076 if(GetReader()->GetDataType() == AliCaloTrackReader::kMC){
f37fa8d2 2077 //Get most probable PID, check PID weights (in MC this option is mandatory)
21a4b1c0 2078 aodph.SetIdentifiedParticleType(GetCaloPID()->GetIdentifiedParticleType(fCalorimeter,calo->GetPID(),mom.E()));//PID with weights
2079 if(GetDebug() > 1) printf("AliAnaPhoton::MakeAnalysisFillAOD() - PDG of identified particle %d\n",aodph.GetIdentifiedParticleType());
f37fa8d2 2080 //If primary is not photon, skip it.
21a4b1c0 2081 if(aodph.GetIdentifiedParticleType() != AliCaloPID::kPhoton) continue ;
6175da48 2082 }
2083 //...............................................
2084 // Data, PID check on
477d6cee 2085 else if(IsCaloPIDOn()){
f37fa8d2 2086 //Get most probable PID, 2 options check PID weights
3d5d5078 2087 //or redo PID, recommended option for MCEal.
477d6cee 2088 if(!IsCaloPIDRecalculationOn())
21a4b1c0 2089 aodph.SetIdentifiedParticleType(GetCaloPID()->GetIdentifiedParticleType(fCalorimeter,calo->GetPID(),mom.E()));//PID with weights
477d6cee 2090 else
21a4b1c0 2091 aodph.SetIdentifiedParticleType(GetCaloPID()->GetIdentifiedParticleType(fCalorimeter,mom,calo));//PID recalculated
477d6cee 2092
21a4b1c0 2093 if(GetDebug() > 1) printf("AliAnaPhoton::MakeAnalysisFillAOD() - PDG of identified particle %d\n",aodph.GetIdentifiedParticleType());
477d6cee 2094
f37fa8d2 2095 //If cluster does not pass pid, not photon, skip it.
21a4b1c0 2096 if(aodph.GetIdentifiedParticleType() != AliCaloPID::kPhoton) continue ;
477d6cee 2097
2098 }
6175da48 2099 //...............................................
2100 // Data, PID check off
477d6cee 2101 else{
f37fa8d2 2102 //Set PID bits for later selection (AliAnaPi0 for example)
2103 //GetPDG already called in SetPIDBits.
f2ccb5b8 2104 GetCaloPID()->SetPIDBits(fCalorimeter,calo,&aodph, GetCaloUtils());
a3aebfff 2105 if(GetDebug() > 1) printf("AliAnaPhoton::MakeAnalysisFillAOD() - PID Bits set \n");
477d6cee 2106 }
2107
21a4b1c0 2108 if(GetDebug() > 1) printf("AliAnaPhoton::MakeAnalysisFillAOD() - Photon selection cuts passed: pT %3.2f, pdg %d\n",aodph.Pt(), aodph.GetIdentifiedParticleType());
477d6cee 2109
477d6cee 2110
6175da48 2111 //--------------------------------------------------------------------------------------
2112 // Conversions pairs analysis
f37fa8d2 2113 // Check if cluster comes from a conversion in the material in front of the calorimeter
2114 // Do invariant mass of all pairs, if mass is close to 0, then it is conversion.
6175da48 2115 //--------------------------------------------------------------------------------------
521636d2 2116
6175da48 2117 // Do analysis only if there are more than one cluster
20218aea 2118 if( nCaloClusters > 1 && fCheckConversion){
c8fe2783 2119 Bool_t bConverted = kFALSE;
2120 Int_t id2 = -1;
1e86c71e 2121
f37fa8d2 2122 //Check if set previously as converted couple, if so skip its use.
20218aea 2123 if (indexConverted[icalo]) continue;
1e86c71e 2124
6175da48 2125 // Second cluster loop
c8fe2783 2126 for(Int_t jcalo = icalo + 1 ; jcalo < nCaloClusters ; jcalo++) {
f37fa8d2 2127 //Check if set previously as converted couple, if so skip its use.
c8fe2783 2128 if (indexConverted[jcalo]) continue;
f37fa8d2 2129 //printf("Check Conversion indeces %d and %d\n",icalo,jcalo);
521636d2 2130 AliVCluster * calo2 = (AliVCluster*) (pl->At(jcalo)); //Get cluster kinematics
6175da48 2131
c4a7d28a 2132
6175da48 2133 //Mixed event, get index of event
c8fe2783 2134 Int_t evtIndex2 = 0 ;
2135 if (GetMixedEvent()) {
2136 evtIndex2=GetMixedEvent()->EventIndexForCaloCluster(calo2->GetID()) ;
f8006433 2137
6175da48 2138 }
2139
2140 //Get kinematics of second cluster
f8006433 2141 if(GetReader()->GetDataType() != AliCaloTrackReader::kMC){
8cdc266d 2142 calo2->GetMomentum(mom2,GetVertex(evtIndex2)) ;}//Assume that come from vertex in straight line
f8006433 2143 else{
2144 Double_t vertex[]={0,0,0};
8cdc266d 2145 calo2->GetMomentum(mom2,vertex) ;
f8006433 2146 }
2147
c4a7d28a 2148 //--------------------------------------
2149 // Cluster selection
2150 //--------------------------------------
2151
2152 if(!ClusterSelected(calo2,mom2)) continue;
2153
2154 //................................................
2155 // Get TOF of each cluster in pair, calculate difference if small,
2156 // take this pair. Way to reject clusters from hadrons (or pileup?)
2157 Double_t t12diff = calo2->GetTOF()-calo->GetTOF()*1e9;
2158 if(TMath::Abs(t12diff) > GetPairTimeCut()) continue;
c8fe2783 2159
6175da48 2160 //................................................
f37fa8d2 2161 //Get mass of pair, if small, take this pair.
41121cfe 2162 Float_t pairM = (mom+mom2).M();
2163 //printf("\t both in calo, mass %f, cut %f\n",pairM,fMassCut);
2164 if(pairM < fMassCut){
2165 aodph.SetTagged(kFALSE);
c8fe2783 2166 id2 = calo2->GetID();
6175da48 2167 indexConverted[icalo]=kTRUE;
c4a7d28a 2168 indexConverted[jcalo]=kTRUE;
41121cfe 2169 Float_t asymmetry = TMath::Abs(mom.E()-mom2.E())/(mom.E()+mom2.E());
2170 Float_t dPhi = mom.Phi()-mom2.Phi();
d39cba7e 2171 Float_t dEta = mom.Eta()-mom2.Eta();
41121cfe 2172
c4a7d28a 2173 //...............................................
2174 //Fill few histograms with kinematics of the pair
2175 //FIXME, move all this to MakeAnalysisFillHistograms ...
2176
2177 fhConvDeltaEta ->Fill( pairM, dPhi );
2178 fhConvDeltaPhi ->Fill( pairM, dEta );
2179 fhConvAsym ->Fill( pairM, asymmetry );
2180 fhConvDeltaEtaPhi->Fill( dEta , dPhi );
2181 fhConvPt ->Fill( pairM, (mom+mom2).Pt());
2182
d39cba7e 2183 //Estimate conversion distance, T. Awes, M. Ivanov
2184 //Under the assumption that the pair has zero mass, and that each electron
2185 //of the pair has the same momentum, they will each have the same bend radius
2186 //given by R=p/(qB) = p / (300 B) with p in [MeV/c], B in [Tesla] and R in [m].
2187 //With nominal ALICE magnet current of 30kA B=0.5T, and so with E_cluster=p,
2188 //R = E/1.5 [cm]. Under these assumptions, the distance from the conversion
3d5d5078 2189 //point to the MCEal can be related to the separation distance, L=2y, on the MCEal
d39cba7e 2190 //as d = sqrt(R^2 -(R-y)^2) = sqrt(2Ry - y^2). And since R>>y we can write as
2191 //d = sqrt(E*L/1.5) where E is the cluster energy and L is the distance in cm between
2192 //the clusters.
c4a7d28a 2193 Float_t pos1[3];
2194 calo->GetPosition(pos1);
2195 Float_t pos2[3];
2196 calo2->GetPosition(pos2);
2197 Float_t clustDist = TMath::Sqrt((pos1[0]-pos2[0])*(pos1[0]-pos2[0])+
521636d2 2198 (pos1[1]-pos2[1])*(pos1[1]-pos2[1])+
2199 (pos1[2]-pos2[2])*(pos1[2]-pos2[2]));
c4a7d28a 2200
2201 Float_t convDist = TMath::Sqrt(mom.E() *clustDist*0.01/0.15);
2202 Float_t convDist2 = TMath::Sqrt(mom2.E()*clustDist*0.01/0.15);
2203 //printf("l = %f, e1 = %f, d1=%f, e2 = %f, d2=%f\n",clustDist,mom.E(),convDist,mom2.E(),convDist2);
8cdc266d 2204 if(GetDebug() > 2)
41121cfe 2205 printf("AliAnaPhoton::MakeAnalysisFillAOD(): Pair with mass %2.3f < %2.3f, %1.2f < dPhi %2.2f < %2.2f, dEta %f < %2.2f, asymmetry %2.2f< %2.2f; \n cluster1 id %d, e %2.3f SM %d, eta %2.3f, phi %2.3f ; \n cluster2 id %d, e %2.3f, SM %d,eta %2.3f, phi %2.3f\n",
2206 pairM,fMassCut,fConvDPhiMinCut, dPhi, fConvDPhiMaxCut, dEta, fConvDEtaCut, asymmetry, fConvAsymCut,
8cdc266d 2207 calo->GetID(),calo->E(),GetCaloUtils()->GetModuleNumber(calo), mom.Eta(), mom.Phi(),
2208 id2, calo2->E(), GetCaloUtils()->GetModuleNumber(calo2),mom2.Eta(), mom2.Phi());
6175da48 2209
c4a7d28a 2210 fhConvDistEta ->Fill(mom .Eta(), convDist );
2211 fhConvDistEta ->Fill(mom2.Eta(), convDist2);
521636d2 2212 fhConvDistEn ->Fill(mom .E(), convDist );
2213 fhConvDistEn ->Fill(mom2.E(), convDist2);
c4a7d28a 2214 fhConvDistMass->Fill((mom+mom2).M(), convDist );
2215 //dEta cut
2216 if(dEta<0.05){
2217 fhConvDistEtaCutEta ->Fill(mom .Eta(), convDist );
2218 fhConvDistEtaCutEta ->Fill(mom2.Eta(), convDist2);
2219 fhConvDistEnCutEta ->Fill(mom .E(), convDist );
2220 fhConvDistEnCutEta ->Fill(mom2.E(), convDist2);
2221 fhConvDistMassCutEta->Fill((mom+mom2).M(), convDist );
2222 //mass cut
2223 if(pairM<0.01){//10 MeV
2224 fhConvDistEtaCutMass ->Fill(mom .Eta(), convDist );
2225 fhConvDistEtaCutMass ->Fill(mom2.Eta(), convDist2);
2226 fhConvDistEnCutMass ->Fill(mom .E(), convDist );
2227 fhConvDistEnCutMass ->Fill(mom2.E(), convDist2);
2228 // asymmetry cut
2229 if(asymmetry<0.1){
2230 fhConvDistEtaCutAsy ->Fill(mom .Eta(), convDist );
2231 fhConvDistEtaCutAsy ->Fill(mom2.Eta(), convDist2);
2232 fhConvDistEnCutAsy ->Fill(mom .E(), convDist );
2233 fhConvDistEnCutAsy ->Fill(mom2.E(), convDist2);
2234 }//asymmetry cut
2235 }//mass cut
2236 }//dEta cut
6175da48 2237
2238 //...............................................
2239 //Select pairs in a eta-phi window
41121cfe 2240 if(TMath::Abs(dEta) < fConvDEtaCut &&
2241 TMath::Abs(dPhi) < fConvDPhiMaxCut &&
2242 TMath::Abs(dPhi) > fConvDPhiMinCut &&
2243 asymmetry < fConvAsymCut ){
521636d2 2244 bConverted = kTRUE;
41121cfe 2245 }
2246 //printf("Accepted? %d\n",bConverted);
6175da48 2247 //...........................................
2248 //Fill more histograms, simulated data
2249 //FIXME, move all this to MakeAnalysisFillHistograms ...
2250 if(IsDataMC()){
2251
2252 //Check the origin of the pair, look for conversion, antinucleons or jet correlations (strings)
2253 Int_t ancPDG = 0;
2254 Int_t ancStatus = 0;
2255 TLorentzVector momentum;
c4a7d28a 2256 TVector3 prodVertex;
6175da48 2257 Int_t ancLabel = GetMCAnalysisUtils()->CheckCommonAncestor(calo->GetLabel(), calo2->GetLabel(),
c4a7d28a 2258 GetReader(), ancPDG, ancStatus, momentum, prodVertex);
6175da48 2259
2260 // printf("AliAnaPhoton::MakeAnalysisFillHistograms() - Common ancestor label %d, pdg %d, name %s, status %d; \n",
2261 // ancLabel,ancPDG,TDatabasePDG::Instance()->GetParticle(ancPDG)->GetName(),ancStatus);
2262
2263 Int_t tag2 = GetMCAnalysisUtils()->CheckOrigin(calo2->GetLabels(),calo2->GetNLabels(),GetReader(), 0);
2264 if(GetMCAnalysisUtils()->CheckTagBit(aodph.GetTag(),AliMCAnalysisUtils::kMCConversion)){
2265 if(GetMCAnalysisUtils()->CheckTagBit(tag2,AliMCAnalysisUtils::kMCConversion) && (ancPDG==22 || TMath::Abs(ancPDG)==11) && ancLabel > -1){
41121cfe 2266 fhConvDeltaEtaMCConversion ->Fill( pairM, dEta );
2267 fhConvDeltaPhiMCConversion ->Fill( pairM, dPhi );
2268 fhConvAsymMCConversion ->Fill( pairM, asymmetry );
2269 fhConvDeltaEtaPhiMCConversion->Fill( dEta , dPhi );
2270 fhConvPtMCConversion ->Fill( pairM, (mom+mom2).Pt());
6175da48 2271 fhConvDispersionMCConversion ->Fill( calo->GetDispersion(), calo2->GetDispersion());
2272 fhConvM02MCConversion ->Fill( calo->GetM02(), calo2->GetM02());
c4a7d28a 2273 fhConvDistMCConversion ->Fill( convDist , prodVertex.Mag() );
2274 fhConvDistMCConversion ->Fill( convDist2, prodVertex.Mag() );
2275
2276 if(dEta<0.05 && pairM<0.01 && asymmetry<0.1){
2277 fhConvDistMCConversionCuts->Fill( convDist , prodVertex.Mag() );
2278 fhConvDistMCConversionCuts->Fill( convDist2, prodVertex.Mag() );
2279 }
521636d2 2280
6175da48 2281 }
2282 }
2283 else if(GetMCAnalysisUtils()->CheckTagBit(aodph.GetTag(),AliMCAnalysisUtils::kMCAntiNeutron)){
2284 if(GetMCAnalysisUtils()->CheckTagBit(tag2,AliMCAnalysisUtils::kMCAntiNeutron) && ancPDG==-2112 && ancLabel > -1){
41121cfe 2285 fhConvDeltaEtaMCAntiNeutron ->Fill( pairM, dEta );
2286 fhConvDeltaPhiMCAntiNeutron ->Fill( pairM, dPhi );
2287 fhConvAsymMCAntiNeutron ->Fill( pairM, asymmetry );
2288 fhConvDeltaEtaPhiMCAntiNeutron ->Fill( dEta , dPhi );
2289 fhConvPtMCAntiNeutron ->Fill( pairM, (mom+mom2).Pt());
6175da48 2290 fhConvDispersionMCAntiNeutron ->Fill( calo->GetDispersion(), calo2->GetDispersion());
2291 fhConvM02MCAntiNeutron ->Fill( calo->GetM02(), calo2->GetM02());
2292 }
2293 }
2294 else if(GetMCAnalysisUtils()->CheckTagBit(aodph.GetTag(),AliMCAnalysisUtils::kMCAntiProton)){
2295 if(GetMCAnalysisUtils()->CheckTagBit(tag2,AliMCAnalysisUtils::kMCAntiProton) && ancPDG==-2212 && ancLabel > -1){
41121cfe 2296 fhConvDeltaEtaMCAntiProton ->Fill( pairM, dEta );
2297 fhConvDeltaPhiMCAntiProton ->Fill( pairM, dPhi );
2298 fhConvAsymMCAntiProton ->Fill( pairM, asymmetry );
2299 fhConvDeltaEtaPhiMCAntiProton ->Fill( dEta , dPhi );
2300 fhConvPtMCAntiProton ->Fill( pairM, (mom+mom2).Pt());
6175da48 2301 fhConvDispersionMCAntiProton ->Fill( calo->GetDispersion(), calo2->GetDispersion());
2302 fhConvM02MCAntiProton ->Fill( calo->GetM02(), calo2->GetM02());
2303 }
2304 }
2305
2306 //Pairs coming from fragmenting pairs.
2307 if(ancPDG < 22 && ancLabel > 7 && (ancStatus == 11 || ancStatus == 12) ){
41121cfe 2308 fhConvDeltaEtaMCString ->Fill( pairM, dPhi);
2309 fhConvDeltaPhiMCString ->Fill( pairM, dPhi);
2310 fhConvAsymMCString ->Fill( pairM, TMath::Abs(mom.E()-mom2.E())/(mom.E()+mom2.E()) );
2311 fhConvDeltaEtaPhiMCString ->Fill( dPhi, dPhi );
2312 fhConvPtMCString ->Fill( pairM, (mom+mom2).Pt());
6175da48 2313 fhConvDispersionMCString ->Fill( calo->GetDispersion(), calo2->GetDispersion());
2314 fhConvM02MCString ->Fill( calo->GetM02(), calo2->GetM02());
2315 }
2316
2317 }// Data MC
521636d2 2318
c8fe2783 2319 break;
2320 }
1e86c71e 2321
c8fe2783 2322 }//Mass loop
1e86c71e 2323
6175da48 2324 //..........................................................................................................
2325 //Pair selected as converted, remove both clusters or recombine them into a photon and put them in the AOD
c8fe2783 2326 if(bConverted){
6175da48 2327 //Add to AOD
c8fe2783 2328 if(fAddConvertedPairsToAOD){
f37fa8d2 2329 //Create AOD of pair analysis
c8fe2783 2330 TLorentzVector mpair = mom+mom2;
2331 AliAODPWG4Particle aodpair = AliAODPWG4Particle(mpair);
2332 aodpair.SetLabel(aodph.GetLabel());
f8006433 2333 //aodpair.SetInputFileIndex(input);
c8fe2783 2334
f37fa8d2 2335 //printf("Index %d, Id %d\n",icalo, calo->GetID());
2336 //Set the indeces of the original caloclusters
c8fe2783 2337 aodpair.SetCaloLabel(calo->GetID(),id2);
2338 aodpair.SetDetector(fCalorimeter);
21a4b1c0 2339 aodpair.SetIdentifiedParticleType(aodph.GetIdentifiedParticleType());
c8fe2783 2340 aodpair.SetTag(aodph.GetTag());
6175da48 2341 aodpair.SetTagged(kTRUE);
f37fa8d2 2342 //Add AOD with pair object to aod branch
c8fe2783 2343 AddAODParticle(aodpair);
f37fa8d2 2344 //printf("\t \t both added pair\n");
c8fe2783 2345 }
2346
f37fa8d2 2347 //Do not add the current calocluster
20218aea 2348 if(fRemoveConvertedPair) continue;
6175da48 2349 else {
2350 //printf("TAGGED\n");
2351 //Tag this cluster as likely conversion
2352 aodph.SetTagged(kTRUE);
2353 }
c8fe2783 2354 }//converted pair
2355 }//check conversion
f37fa8d2 2356 //printf("\t \t added single cluster %d\n",icalo);
1e86c71e 2357
6175da48 2358 //FIXME, this to MakeAnalysisFillHistograms ...
c4a7d28a 2359 fhNCellsE->Fill(aodph.E(),calo->GetNCells());
6175da48 2360
f37fa8d2 2361 //Add AOD with photon object to aod branch
477d6cee 2362 AddAODParticle(aodph);
2363
2364 }//loop
2365
4a745797 2366 delete [] indexConverted;
2367
f37fa8d2 2368 if(GetDebug() > 1) printf("AliAnaPhoton::MakeAnalysisFillAOD() End fill AODs, with %d entries \n",GetOutputAODBranch()->GetEntriesFast());
477d6cee 2369
1c5acb87 2370}
2371
2372//__________________________________________________________________
2373void AliAnaPhoton::MakeAnalysisFillHistograms()
2374{
6175da48 2375 //Fill histograms
f8006433 2376
6175da48 2377 //-------------------------------------------------------------------
577d9801 2378 // Access MC information in stack if requested, check that it exists.
521636d2 2379 AliStack * stack = 0x0;
2380 TParticle * primary = 0x0;
2381 TClonesArray * mcparticles = 0x0;
2382 AliAODMCParticle * aodprimary = 0x0;
2383
577d9801 2384 if(IsDataMC()){
521636d2 2385
577d9801 2386 if(GetReader()->ReadStack()){
2387 stack = GetMCStack() ;
2388 if(!stack) {
3d5d5078 2389 printf("AliAnaPhoton::MakeAnalysisFillHistograms() - Stack not available, is the MC handler called? STOP\n");
2390 abort();
577d9801 2391 }
f8006433 2392
577d9801 2393 }
2394 else if(GetReader()->ReadAODMCParticles()){
f8006433 2395
577d9801 2396 //Get the list of MC particles
521636d2 2397 mcparticles = GetReader()->GetAODMCParticles(0);
2398 if(!mcparticles && GetDebug() > 0) {
3d5d5078 2399 printf("AliAnaPhoton::MakeAnalysisFillHistograms() - Standard MCParticles not available!\n");
577d9801 2400 }
577d9801 2401 }
2402 }// is data and MC
521636d2 2403
6175da48 2404
2405 // Get vertex
2244659d 2406 Double_t v[3] = {0,0,0}; //vertex ;
2407 GetReader()->GetVertex(v);
6175da48 2408 //fhVertex->Fill(v[0],v[1],v[2]);
2409 if(TMath::Abs(v[2]) > GetZvertexCut()) return ; // done elsewhere for Single Event analysis, but there for mixed event
2410
2411 //----------------------------------
577d9801 2412 //Loop on stored AOD photons
2413 Int_t naod = GetOutputAODBranch()->GetEntriesFast();
577d9801 2414 if(GetDebug() > 0) printf("AliAnaPhoton::MakeAnalysisFillHistograms() - aod branch entries %d\n", naod);
521636d2 2415
577d9801 2416 for(Int_t iaod = 0; iaod < naod ; iaod++){
2417 AliAODPWG4Particle* ph = (AliAODPWG4Particle*) (GetOutputAODBranch()->At(iaod));
2418 Int_t pdg = ph->GetIdentifiedParticleType();
521636d2 2419
577d9801 2420 if(GetDebug() > 3)
2421 printf("AliAnaPhoton::MakeAnalysisFillHistograms() - PDG %d, MC TAG %d, Calorimeter %s\n", ph->GetIdentifiedParticleType(),ph->GetTag(), (ph->GetDetector()).Data()) ;
521636d2 2422
577d9801 2423 //If PID used, fill histos with photons in Calorimeter fCalorimeter
2424 if(IsCaloPIDOn() && pdg != AliCaloPID::kPhoton) continue;
2425 if(ph->GetDetector() != fCalorimeter) continue;
521636d2 2426
577d9801 2427 if(GetDebug() > 2)
2428 printf("AliAnaPhoton::MakeAnalysisFillHistograms() - ID Photon: pt %f, phi %f, eta %f\n", ph->Pt(),ph->Phi(),ph->Eta()) ;
521636d2 2429
6175da48 2430 //................................
577d9801 2431 //Fill photon histograms
2432 Float_t ptcluster = ph->Pt();
2433 Float_t phicluster = ph->Phi();
2434 Float_t etacluster = ph->Eta();
2435 Float_t ecluster = ph->E();
521636d2 2436
20218aea 2437 fhEPhoton ->Fill(ecluster);
577d9801 2438 fhPtPhoton ->Fill(ptcluster);
2439 fhPhiPhoton ->Fill(ptcluster,phicluster);
2440 fhEtaPhoton ->Fill(ptcluster,etacluster);
521636d2 2441 if (ecluster > 0.5) fhEtaPhiPhoton ->Fill(etacluster, phicluster);
20218aea 2442 else if(GetMinPt() < 0.5) fhEtaPhi05Photon->Fill(etacluster, phicluster);
2443
2444 if(fCheckConversion &&ph->IsTagged()){
2445 fhPtPhotonConv->Fill(ptcluster);
2446 if(ecluster > 0.5) fhEtaPhiPhotonConv ->Fill(etacluster, phicluster);
2447 else if(GetMinPt() < 0.5) fhEtaPhi05PhotonConv->Fill(etacluster, phicluster);
6175da48 2448 }
20218aea 2449
6175da48 2450 //.......................................
577d9801 2451 //Play with the MC data if available
2452 if(IsDataMC()){
521636d2 2453
3d5d5078 2454 FillAcceptanceHistograms();
2455
577d9801 2456 Int_t tag =ph->GetTag();
521636d2 2457
577d9801 2458 if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPhoton))
3d5d5078 2459 {
2460 fhMCE [mcPhoton] ->Fill(ecluster);
2461 fhPtMC [mcPhoton] ->Fill(ptcluster);
2462 fhPhiMC[mcPhoton] ->Fill(ecluster,phicluster);
2463 fhEtaMC[mcPhoton] ->Fill(ecluster,etacluster);
2464
2465 if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCConversion))
2466 {
2467 fhMCE [mcConversion] ->Fill(ecluster);
2468 fhPtMC [mcConversion] ->Fill(ptcluster);
2469 fhPhiMC[mcConversion] ->Fill(ecluster,phicluster);
2470 fhEtaMC[mcConversion] ->Fill(ecluster,etacluster);
2471
2472 if(fCheckConversion){
2473 if(ph->IsTagged()) fhPtConversionTagged ->Fill(ptcluster);
2474 if(ptcluster > 0.5)fhEtaPhiConversion ->Fill(etacluster,phicluster);
2475 else fhEtaPhi05Conversion ->Fill(etacluster,phicluster);
2476 }
2477 }
2478
2479 if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPrompt)){
2480 fhMCE [mcPrompt] ->Fill(ecluster);
2481 fhPtMC [mcPrompt] ->Fill(ptcluster);
2482 fhPhiMC[mcPrompt] ->Fill(ecluster,phicluster);
2483 fhEtaMC[mcPrompt] ->Fill(ecluster,etacluster);
2484 }
2485 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCFragmentation))
2486 {
2487 fhMCE [mcFragmentation] ->Fill(ecluster);
2488 fhPtMC [mcFragmentation] ->Fill(ptcluster);
2489 fhPhiMC[mcFragmentation] ->Fill(ecluster,phicluster);
2490 fhEtaMC[mcFragmentation] ->Fill(ecluster,etacluster);
2491
2492 }
2493 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCISR))
2494 {
2495 fhMCE [mcISR] ->Fill(ecluster);
2496 fhPtMC [mcISR] ->Fill(ptcluster);
2497 fhPhiMC[mcISR] ->Fill(ecluster,phicluster);
2498 fhEtaMC[mcISR] ->Fill(ecluster,etacluster);
2499 }
2500 else if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0Decay) &&
2501 !GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0))
2502 {
2503 fhMCE [mcPi0Decay] ->Fill(ecluster);
2504 fhPtMC [mcPi0Decay] ->Fill(ptcluster);
2505 fhPhiMC[mcPi0Decay] ->Fill(ecluster,phicluster);
2506 fhEtaMC[mcPi0Decay] ->Fill(ecluster,etacluster);
2507 }
2508 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCEtaDecay) ||
2509 GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCOtherDecay))
2510 {
2511 fhMCE [mcOtherDecay] ->Fill(ecluster);
2512 fhPtMC [mcOtherDecay] ->Fill(ptcluster);
2513 fhPhiMC[mcOtherDecay] ->Fill(ecluster,phicluster);
2514 fhEtaMC[mcOtherDecay] ->Fill(ecluster,etacluster);
2515 }
2516 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0))
2517 {
2518 fhMCE [mcPi0] ->Fill(ecluster);
2519 fhPtMC [mcPi0] ->Fill(ptcluster);
2520 fhPhiMC[mcPi0] ->Fill(ecluster,phicluster);
2521 fhEtaMC[mcPi0] ->Fill(ecluster,etacluster);
2522 }
2523 }
6175da48 2524 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCAntiNeutron))
3d5d5078 2525 {
2526 fhMCE [mcAntiNeutron] ->Fill(ecluster);
2527 fhPtMC [mcAntiNeutron] ->Fill(ptcluster);
2528 fhPhiMC[mcAntiNeutron] ->Fill(ecluster,phicluster);
2529 fhEtaMC[mcAntiNeutron] ->Fill(ecluster,etacluster);
2530 if(ph->IsTagged() && fCheckConversion) fhPtAntiNeutronTagged ->Fill(ptcluster);
2531
2532 }
6175da48 2533 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCAntiProton))
3d5d5078 2534 {
2535 fhMCE [mcAntiProton] ->Fill(ecluster);
2536 fhPtMC [mcAntiProton] ->Fill(ptcluster);
2537 fhPhiMC[mcAntiProton] ->Fill(ecluster,phicluster);
2538 fhEtaMC[mcAntiProton] ->Fill(ecluster,etacluster);
2539 if(ph->IsTagged() && fCheckConversion) fhPtAntiProtonTagged ->Fill(ptcluster);
2540 }
521636d2 2541 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCElectron))
3d5d5078 2542 {
2543 fhMCE [mcElectron] ->Fill(ecluster);
2544 fhPtMC [mcElectron] ->Fill(ptcluster);
2545 fhPhiMC[mcElectron] ->Fill(ecluster,phicluster);
2546 fhEtaMC[mcElectron] ->Fill(ecluster,etacluster);
2547 }
577d9801 2548 else{
3d5d5078 2549 fhMCE [mcOther] ->Fill(ecluster);
521636d2 2550 fhPtMC [mcOther] ->Fill(ptcluster);
2551 fhPhiMC[mcOther] ->Fill(ecluster,phicluster);
2552 fhEtaMC[mcOther] ->Fill(ecluster,etacluster);
20218aea 2553 if(ph->IsTagged() && fCheckConversion) fhPtUnknownTagged ->Fill(ptcluster);
521636d2 2554
3d5d5078 2555
f8006433 2556 // printf(" AliAnaPhoton::MakeAnalysisFillHistograms() - Label %d, pT %2.3f Unknown, bits set: ",
2557 // ph->GetLabel(),ph->Pt());
2558 // for(Int_t i = 0; i < 20; i++) {
2559 // if(GetMCAnalysisUtils()->CheckTagBit(tag,i)) printf(" %d, ",i);
2560 // }
2561 // printf("\n");
2562
577d9801 2563 }
521636d2 2564
577d9801 2565 //....................................................................
2566 // Access MC information in stack if requested, check that it exists.
2567 Int_t label =ph->GetLabel();
2568 if(label < 0) {
3d5d5078 2569 if(GetDebug() > 1) printf("AliAnaPhoton::MakeAnalysisFillHistograms() *** bad label ***: label %d \n", label);
2570 continue;
577d9801 2571 }
521636d2 2572
577d9801 2573 Float_t eprim = 0;
2574 Float_t ptprim = 0;
2575 if(GetReader()->ReadStack()){
3d5d5078 2576
2577 if(label >= stack->GetNtrack()) {
f8006433 2578 if(GetDebug() > 2) printf("AliAnaPhoton::MakeAnalysisFillHistograms() *** large label ***: label %d, n tracks %d \n", label, stack->GetNtrack());
2579 continue ;
3d5d5078 2580 }
2581
2582 primary = stack->Particle(label);
2583 if(!primary){
2584 printf("AliAnaPhoton::MakeAnalysisFillHistograms() *** no primary ***: label %d \n", label);
2585 continue;
2586 }
2587 eprim = primary->Energy();
2588 ptprim = primary->Pt();
2589
577d9801 2590 }
2591 else if(GetReader()->ReadAODMCParticles()){
3d5d5078 2592 //Check which is the input
2593 if(ph->GetInputFileIndex() == 0){
2594 if(!mcparticles) continue;
2595 if(label >= mcparticles->GetEntriesFast()) {
2596 if(GetDebug() > 2) printf("AliAnaPhoton::MakeAnalysisFillHistograms() *** large label ***: label %d, n tracks %d \n",
521636d2 2597 label, mcparticles->GetEntriesFast());
3d5d5078 2598 continue ;
f8006433 2599 }
2600 //Get the particle
521636d2 2601 aodprimary = (AliAODMCParticle*) mcparticles->At(label);
f8006433 2602
3d5d5078 2603 }
2604
2605 if(!aodprimary){
2606 printf("AliAnaPhoton::MakeAnalysisFillHistograms() *** no primary ***: label %d \n", label);
2607 continue;
2608 }
2609
2610 eprim = aodprimary->E();
2611 ptprim = aodprimary->Pt();
2612
577d9801 2613 }
521636d2 2614
577d9801 2615 fh2E ->Fill(ecluster, eprim);
2616 fh2Pt ->Fill(ptcluster, ptprim);
2617 fhDeltaE ->Fill(eprim-ecluster);
2618 fhDeltaPt->Fill(ptprim-ptcluster);
2619 if(eprim > 0) fhRatioE ->Fill(ecluster/eprim);
2620 if(ptprim > 0) fhRatioPt ->Fill(ptcluster/ptprim);
521636d2 2621
577d9801 2622 }//Histograms with MC
521636d2 2623
577d9801 2624 }// aod loop
521636d2 2625
1c5acb87 2626}
2627
2628
2629//__________________________________________________________________
2630void AliAnaPhoton::Print(const Option_t * opt) const
2631{
477d6cee 2632 //Print some relevant parameters set for the analysis
2633
2634 if(! opt)
2635 return;
2636
2637 printf("**** Print %s %s ****\n", GetName(), GetTitle() ) ;
2638 AliAnaPartCorrBaseClass::Print(" ");
a3aebfff 2639
477d6cee 2640 printf("Calorimeter = %s\n", fCalorimeter.Data()) ;
2641 printf("Min Distance to Bad Channel = %2.1f\n",fMinDist);
2642 printf("Min Distance to Bad Channel 2 = %2.1f\n",fMinDist2);
2643 printf("Min Distance to Bad Channel 3 = %2.1f\n",fMinDist3);
a3aebfff 2644 printf("Reject clusters with a track matched = %d\n",fRejectTrackMatch);
1e86c71e 2645 printf("Check Pair Conversion = %d\n",fCheckConversion);
2646 printf("Add conversion pair to AOD = %d\n",fAddConvertedPairsToAOD);
2647 printf("Conversion pair mass cut = %f\n",fMassCut);
41121cfe 2648 printf("Conversion selection cut : A < %1.2f; %1.3f < Dphi < %1.3f; Deta < %1.3f\n",
2649 fConvAsymCut,fConvDPhiMinCut, fConvDPhiMaxCut, fConvDEtaCut);
4cf55759 2650 printf("Time Cut: %3.1f < TOF < %3.1f\n", fTimeCutMin, fTimeCutMax);
2ac125bf 2651 printf("Number of cells in cluster is > %d \n", fNCellsCut);
477d6cee 2652 printf(" \n") ;
1c5acb87 2653
2654}