]> git.uio.no Git - u/mrichter/AliRoot.git/blame_incremental - PWGGA/GammaConv/AliAnalysisTaskGammaConvCalo.cxx
MC changes to evsel for debugging
[u/mrichter/AliRoot.git] / PWGGA / GammaConv / AliAnalysisTaskGammaConvCalo.cxx
... / ...
CommitLineData
1/**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: Baldo Sahlmueller, Friederike Bock *
5 * Version 1.0 *
6 * *
7 * *
8 * Permission to use, copy, modify and distribute this software and its *
9 * documentation strictly for non-commercial purposes is hereby granted *
10 * without fee, provided that the above copyright notice appears in all *
11 * copies and that both the copyright notice and this permission notice *
12 * appear in the supporting documentation. The authors make no claims *
13 * about the suitability of this software for any purpose. It is *
14 * provided "as is" without express or implied warranty. *
15 **************************************************************************/
16
17//////////////////////////////////////////////////////////////////
18//----------------------------------------------------------------
19// Class used to do analysis on conversion photons + calo photons
20//----------------------------------------------------------------
21//////////////////////////////////////////////////////////////////
22#include "TChain.h"
23#include "TTree.h"
24#include "TBranch.h"
25#include "TFile.h"
26#include "TH1F.h"
27#include "TH2F.h"
28#include "TH3F.h"
29#include "THnSparse.h"
30#include "TCanvas.h"
31#include "TNtuple.h"
32#include "AliAnalysisTask.h"
33#include "AliAnalysisManager.h"
34#include "AliESDEvent.h"
35#include "AliESDInputHandler.h"
36#include "AliMCEventHandler.h"
37#include "AliMCEvent.h"
38#include "AliMCParticle.h"
39#include "AliCentrality.h"
40#include "AliESDVZERO.h"
41#include "AliESDpid.h"
42#include "AliAnalysisTaskGammaConvCalo.h"
43#include "AliVParticle.h"
44#include "AliESDtrack.h"
45#include "AliESDtrackCuts.h"
46#include "AliKFVertex.h"
47#include "AliV0ReaderV1.h"
48#include "AliGenCocktailEventHeader.h"
49#include "AliConversionAODBGHandlerRP.h"
50#include "AliAODMCParticle.h"
51#include "AliAODMCHeader.h"
52#include "AliEventplane.h"
53#include "AliAnalysisTaskEMCALClusterizeFast.h"
54#include "AliAODEvent.h"
55#include "AliAODInputHandler.h"
56#include "AliESDEvent.h"
57#include "AliESDInputHandler.h"
58#include "AliInputEventHandler.h"
59
60ClassImp(AliAnalysisTaskGammaConvCalo)
61
62//________________________________________________________________________
63AliAnalysisTaskGammaConvCalo::AliAnalysisTaskGammaConvCalo(): AliAnalysisTaskSE(),
64 fV0Reader(NULL),
65 fBGHandler(NULL),
66 fBGHandlerRP(NULL),
67 fBGClusHandler(NULL),
68 fBGClusHandlerRP(NULL),
69 fInputEvent(NULL),
70 fMCEvent(NULL),
71 fMCStack(NULL),
72 fCutFolder(NULL),
73 fESDList(NULL),
74 fBackList(NULL),
75 fMotherList(NULL),
76 fPhotonDCAList(NULL),
77 fMesonDCAList(NULL),
78 fTrueList(NULL),
79 fMCList(NULL),
80 fHeaderNameList(NULL),
81 fTagOutputList(NULL),
82 fOutputContainer(NULL),
83 fReaderGammas(NULL),
84 fGammaCandidates(NULL),
85 fClusterCandidates(NULL),
86 fCutArray(NULL),
87 fConversionCuts(NULL),
88 fClusterCutArray(NULL),
89 fCaloPhotonCuts(NULL),
90 fMesonCutArray(NULL),
91 fMesonCuts(NULL),
92 fHistoConvGammaPt(NULL),
93 fHistoConvGammaR(NULL),
94 fHistoConvGammaEta(NULL),
95 fTreeConvGammaPtDcazCat(NULL),
96 fPtGamma(0),
97 fDCAzPhoton(0),
98 fRConvPhoton(0),
99 fEtaPhoton(0),
100 fCharCatPhoton(0),
101 fCharPhotonMCInfo(0),
102 fHistoMotherInvMassPt(NULL),
103 fSparseMotherInvMassPtZM(NULL),
104 fHistoMotherBackInvMassPt(NULL),
105 fSparseMotherBackInvMassPtZM(NULL),
106 fHistoMotherInvMassEalpha(NULL),
107 fHistoMotherPi0PtY(NULL),
108 fHistoMotherEtaPtY(NULL),
109 fHistoMotherPi0PtAlpha(NULL),
110 fHistoMotherEtaPtAlpha(NULL),
111 fHistoMotherPi0PtOpenAngle(NULL),
112 fHistoMotherEtaPtOpenAngle(NULL),
113 fTreeMesonsInvMassPtDcazMinDcazMaxFlag(NULL),
114 fInvMass(0),
115 fPt(0),
116 fDCAzGammaMin(0),
117 fDCAzGammaMax(0),
118 fCharFlag(0),
119 fCharMesonMCInfo(0),
120 fHistoConvGammaUntagged(NULL),
121 fHistoConvGammaTagged(NULL),
122 fHistoConvGammaPi0Tagged(NULL),
123 fHistoConvGammaEtaTagged(NULL),
124 fHistoPhotonPairAll(NULL),
125 fHistoPhotonPairAllGam(NULL),
126 fHistoClusGammaPt(NULL),
127 fHistoMCHeaders(NULL),
128 fHistoMCAllGammaPt(NULL),
129 fHistoMCDecayGammaPi0Pt(NULL),
130 fHistoMCDecayGammaRhoPt(NULL),
131 fHistoMCDecayGammaEtaPt(NULL),
132 fHistoMCDecayGammaOmegaPt(NULL),
133 fHistoMCDecayGammaEtapPt(NULL),
134 fHistoMCDecayGammaPhiPt(NULL),
135 fHistoMCDecayGammaSigmaPt(NULL),
136 fHistoMCConvGammaPt(NULL),
137 fHistoMCConvGammaR(NULL),
138 fHistoMCConvGammaEta(NULL),
139 fHistoMCPi0Pt(NULL),
140 fHistoMCPi0WOWeightPt(NULL),
141 fHistoMCEtaPt(NULL),
142 fHistoMCEtaWOWeightPt(NULL),
143 fHistoMCPi0InAccPt(NULL),
144 fHistoMCEtaInAccPt(NULL),
145 fHistoMCPi0PtY(NULL),
146 fHistoMCEtaPtY(NULL),
147 fHistoMCK0sPt(NULL),
148 fHistoMCK0sWOWeightPt(NULL),
149 fHistoMCK0sPtY(NULL),
150 fHistoMCSecPi0PtvsSource(NULL),
151 fHistoMCSecPi0Source(NULL),
152 fHistoMCSecEtaPt(NULL),
153 fHistoMCSecEtaSource(NULL),
154 fHistoTrueMotherInvMassPt(NULL),
155 fHistoTruePrimaryMotherInvMassPt(NULL),
156 fHistoTruePrimaryMotherW0WeightingInvMassPt(NULL),
157 fProfileTruePrimaryMotherWeightsInvMassPt(NULL),
158 fHistoTruePrimaryPi0MCPtResolPt(NULL),
159 fHistoTruePrimaryEtaMCPtResolPt(NULL),
160 fHistoTrueSecondaryMotherInvMassPt(NULL),
161 fHistoTrueSecondaryMotherFromK0sInvMassPt(NULL),
162 fHistoTrueK0sWithPi0DaughterMCPt(NULL),
163 fHistoTrueSecondaryMotherFromEtaInvMassPt(NULL),
164 fHistoTrueEtaWithPi0DaughterMCPt(NULL),
165 fHistoTrueSecondaryMotherFromLambdaInvMassPt(NULL),
166 fHistoTrueLambdaWithPi0DaughterMCPt(NULL),
167 fHistoTrueBckGGInvMassPt(NULL),
168 fHistoTrueBckContInvMassPt(NULL),
169 fHistoTruePi0PtY(NULL),
170 fHistoTrueEtaPtY(NULL),
171 fHistoTruePi0PtAlpha(NULL),
172 fHistoTrueEtaPtAlpha(NULL),
173 fHistoTruePi0PtOpenAngle(NULL),
174 fHistoTrueEtaPtOpenAngle(NULL),
175 fHistoTrueMotherDalitzInvMassPt(NULL),
176 fHistoTrueConvGammaPt(NULL),
177 fHistoTrueConvPi0GammaPt(NULL),
178 fHistoTrueConvGammaEta(NULL),
179 fHistoCombinatorialPt(NULL),
180 fHistoTruePrimaryConvGammaPt(NULL),
181 fHistoTruePrimaryConvGammaESDPtMCPt(NULL),
182 fHistoTrueSecondaryConvGammaPt(NULL),
183 fHistoTrueSecondaryConvGammaFromXFromK0sPt(NULL),
184 fHistoTrueSecondaryConvGammaFromXFromLambdaPt(NULL),
185 fHistoTrueDalitzPsiPairDeltaPhi(NULL),
186 fHistoTrueGammaPsiPairDeltaPhi(NULL),
187 fHistoTrueClusGammaPt(NULL),
188 fHistoTruePrimaryClusGammaPt(NULL),
189 fHistoTruePrimaryClusGammaESDPtMCPt(NULL),
190 fHistoNEvents(NULL),
191 fHistoNGoodESDTracks(NULL),
192 fHistoNGammaCandidates(NULL),
193 fHistoNGoodESDTracksVsNGammaCanditates(NULL),
194 fHistoNV0Tracks(NULL),
195 fProfileEtaShift(NULL),
196 fEventPlaneAngle(-100),
197 fRandom(0),
198 fNGammaCandidates(0),
199 fUnsmearedPx(NULL),
200 fUnsmearedPy(NULL),
201 fUnsmearedPz(NULL),
202 fUnsmearedE(NULL),
203 fMCStackPos(NULL),
204 fMCStackNeg(NULL),
205 fESDArrayPos(NULL),
206 fESDArrayNeg(NULL),
207 fnCuts(0),
208 fiCut(0),
209 fMoveParticleAccordingToVertex(kTRUE),
210 fIsHeavyIon(0),
211 fDoMesonAnalysis(kTRUE),
212 fDoMesonQA(0),
213 fDoPhotonQA(0),
214 fIsFromMBHeader(kTRUE),
215 fIsMC(kFALSE),
216 fMinE(0.1),
217 fNminCells(2),
218 fEMCm02cut(0.5)
219{
220
221}
222
223//________________________________________________________________________
224AliAnalysisTaskGammaConvCalo::AliAnalysisTaskGammaConvCalo(const char *name):
225 AliAnalysisTaskSE(name),
226 fV0Reader(NULL),
227 fBGHandler(NULL),
228 fBGHandlerRP(NULL),
229 fBGClusHandler(NULL),
230 fBGClusHandlerRP(NULL),
231 fInputEvent(NULL),
232 fMCEvent(NULL),
233 fMCStack(NULL),
234 fCutFolder(NULL),
235 fESDList(NULL),
236 fBackList(NULL),
237 fMotherList(NULL),
238 fPhotonDCAList(NULL),
239 fMesonDCAList(NULL),
240 fTrueList(NULL),
241 fMCList(NULL),
242 fHeaderNameList(NULL),
243 fTagOutputList(NULL),
244 fOutputContainer(0),
245 fReaderGammas(NULL),
246 fGammaCandidates(NULL),
247 fClusterCandidates(NULL),
248 fCutArray(NULL),
249 fConversionCuts(NULL),
250 fClusterCutArray(NULL),
251 fCaloPhotonCuts(NULL),
252 fMesonCutArray(NULL),
253 fMesonCuts(NULL),
254 fHistoConvGammaPt(NULL),
255 fHistoConvGammaR(NULL),
256 fHistoConvGammaEta(NULL),
257 fTreeConvGammaPtDcazCat(NULL),
258 fPtGamma(0),
259 fDCAzPhoton(0),
260 fRConvPhoton(0),
261 fEtaPhoton(0),
262 fCharCatPhoton(0),
263 fCharPhotonMCInfo(0),
264 fHistoMotherInvMassPt(NULL),
265 fSparseMotherInvMassPtZM(NULL),
266 fHistoMotherBackInvMassPt(NULL),
267 fSparseMotherBackInvMassPtZM(NULL),
268 fHistoMotherInvMassEalpha(NULL),
269 fHistoMotherPi0PtY(NULL),
270 fHistoMotherEtaPtY(NULL),
271 fHistoMotherPi0PtAlpha(NULL),
272 fHistoMotherEtaPtAlpha(NULL),
273 fHistoMotherPi0PtOpenAngle(NULL),
274 fHistoMotherEtaPtOpenAngle(NULL),
275 fTreeMesonsInvMassPtDcazMinDcazMaxFlag(NULL),
276 fInvMass(0),
277 fPt(0),
278 fDCAzGammaMin(0),
279 fDCAzGammaMax(0),
280 fCharFlag(0),
281 fCharMesonMCInfo(0),
282 fHistoConvGammaUntagged(NULL),
283 fHistoConvGammaTagged(NULL),
284 fHistoConvGammaPi0Tagged(NULL),
285 fHistoConvGammaEtaTagged(NULL),
286 fHistoPhotonPairAll(NULL),
287 fHistoPhotonPairAllGam(NULL),
288 fHistoClusGammaPt(NULL),
289 fHistoMCHeaders(NULL),
290 fHistoMCAllGammaPt(NULL),
291 fHistoMCDecayGammaPi0Pt(NULL),
292 fHistoMCDecayGammaRhoPt(NULL),
293 fHistoMCDecayGammaEtaPt(NULL),
294 fHistoMCDecayGammaOmegaPt(NULL),
295 fHistoMCDecayGammaEtapPt(NULL),
296 fHistoMCDecayGammaPhiPt(NULL),
297 fHistoMCDecayGammaSigmaPt(NULL),
298 fHistoMCConvGammaPt(NULL),
299 fHistoMCConvGammaR(NULL),
300 fHistoMCConvGammaEta(NULL),
301 fHistoMCPi0Pt(NULL),
302 fHistoMCPi0WOWeightPt(NULL),
303 fHistoMCEtaPt(NULL),
304 fHistoMCEtaWOWeightPt(NULL),
305 fHistoMCPi0InAccPt(NULL),
306 fHistoMCEtaInAccPt(NULL),
307 fHistoMCPi0PtY(NULL),
308 fHistoMCEtaPtY(NULL),
309 fHistoMCK0sPt(NULL),
310 fHistoMCK0sWOWeightPt(NULL),
311 fHistoMCK0sPtY(NULL),
312 fHistoMCSecPi0PtvsSource(NULL),
313 fHistoMCSecPi0Source(NULL),
314 fHistoMCSecEtaPt(NULL),
315 fHistoMCSecEtaSource(NULL),
316 fHistoTrueMotherInvMassPt(NULL),
317 fHistoTruePrimaryMotherInvMassPt(NULL),
318 fHistoTruePrimaryMotherW0WeightingInvMassPt(NULL),
319 fProfileTruePrimaryMotherWeightsInvMassPt(NULL),
320 fHistoTruePrimaryPi0MCPtResolPt(NULL),
321 fHistoTruePrimaryEtaMCPtResolPt(NULL),
322 fHistoTrueSecondaryMotherInvMassPt(NULL),
323 fHistoTrueSecondaryMotherFromK0sInvMassPt(NULL),
324 fHistoTrueK0sWithPi0DaughterMCPt(NULL),
325 fHistoTrueSecondaryMotherFromEtaInvMassPt(NULL),
326 fHistoTrueEtaWithPi0DaughterMCPt(NULL),
327 fHistoTrueSecondaryMotherFromLambdaInvMassPt(NULL),
328 fHistoTrueLambdaWithPi0DaughterMCPt(NULL),
329 fHistoTrueBckGGInvMassPt(NULL),
330 fHistoTrueBckContInvMassPt(NULL),
331 fHistoTruePi0PtY(NULL),
332 fHistoTrueEtaPtY(NULL),
333 fHistoTruePi0PtAlpha(NULL),
334 fHistoTrueEtaPtAlpha(NULL),
335 fHistoTruePi0PtOpenAngle(NULL),
336 fHistoTrueEtaPtOpenAngle(NULL),
337 fHistoTrueMotherDalitzInvMassPt(NULL),
338 fHistoTrueConvGammaPt(NULL),
339 fHistoTrueConvPi0GammaPt(NULL),
340 fHistoTrueConvGammaEta(NULL),
341 fHistoCombinatorialPt(NULL),
342 fHistoTruePrimaryConvGammaPt(NULL),
343 fHistoTruePrimaryConvGammaESDPtMCPt(NULL),
344 fHistoTrueSecondaryConvGammaPt(NULL),
345 fHistoTrueSecondaryConvGammaFromXFromK0sPt(NULL),
346 fHistoTrueSecondaryConvGammaFromXFromLambdaPt(NULL),
347 fHistoTrueDalitzPsiPairDeltaPhi(NULL),
348 fHistoTrueGammaPsiPairDeltaPhi(NULL),
349 fHistoTrueClusGammaPt(NULL),
350 fHistoTruePrimaryClusGammaPt(NULL),
351 fHistoTruePrimaryClusGammaESDPtMCPt(NULL),
352 fHistoNEvents(NULL),
353 fHistoNGoodESDTracks(NULL),
354 fHistoNGammaCandidates(NULL),
355 fHistoNGoodESDTracksVsNGammaCanditates(NULL),
356 fHistoNV0Tracks(NULL),
357 fProfileEtaShift(NULL),
358 fEventPlaneAngle(-100),
359 fRandom(0),
360 fNGammaCandidates(0),
361 fUnsmearedPx(NULL),
362 fUnsmearedPy(NULL),
363 fUnsmearedPz(NULL),
364 fUnsmearedE(NULL),
365 fMCStackPos(NULL),
366 fMCStackNeg(NULL),
367 fESDArrayPos(NULL),
368 fESDArrayNeg(NULL),
369 fnCuts(0),
370 fiCut(0),
371 fMoveParticleAccordingToVertex(kTRUE),
372 fIsHeavyIon(0),
373 fDoMesonAnalysis(kTRUE),
374 fDoMesonQA(0),
375 fDoPhotonQA(0),
376 fIsFromMBHeader(kTRUE),
377 fIsMC(kFALSE),
378 fMinE(0.1),
379 fNminCells(2),
380 fEMCm02cut(0.5)
381{
382 // Define output slots here
383 DefineOutput(1, TList::Class());
384}
385
386AliAnalysisTaskGammaConvCalo::~AliAnalysisTaskGammaConvCalo()
387{
388 if(fGammaCandidates){
389 delete fGammaCandidates;
390 fGammaCandidates = 0x0;
391 }
392 if(fClusterCandidates){
393 delete fClusterCandidates;
394 fClusterCandidates = 0x0;
395 }
396 if(fBGHandler){
397 delete[] fBGHandler;
398 fBGHandler = 0x0;
399 }
400 if(fBGHandlerRP){
401 delete[] fBGHandlerRP;
402 fBGHandlerRP = 0x0;
403 }
404 if(fBGClusHandler){
405 delete[] fBGClusHandler;
406 fBGClusHandler = 0x0;
407 }
408 if(fBGClusHandlerRP){
409 delete[] fBGClusHandlerRP;
410 fBGClusHandlerRP = 0x0;
411 }
412}
413//___________________________________________________________
414void AliAnalysisTaskGammaConvCalo::InitBack(){
415
416 const Int_t nDim = 4;
417 Int_t nBins[nDim] = {800,250,7,4};
418 Double_t xMin[nDim] = {0,0, 0,0};
419 Double_t xMax[nDim] = {0.8,25,7,4};
420
421 fSparseMotherInvMassPtZM = new THnSparseF*[fnCuts];
422 fSparseMotherBackInvMassPtZM = new THnSparseF*[fnCuts];
423
424 fBGHandler = new AliGammaConversionAODBGHandler*[fnCuts];
425 fBGHandlerRP = new AliConversionAODBGHandlerRP*[fnCuts];
426
427 fBGClusHandler = new AliGammaConversionAODBGHandler*[fnCuts];
428 fBGClusHandlerRP = new AliConversionAODBGHandlerRP*[fnCuts];
429
430 for(Int_t iCut = 0; iCut<fnCuts;iCut++){
431 if (((AliConversionMesonCuts*)fMesonCutArray->At(iCut))->DoBGCalculation()){
432 TString cutstring = ((AliConversionCuts*)fCutArray->At(iCut))->GetCutNumber();
433 TString cutstringCalo = ((AliCaloPhotonCuts*)fClusterCutArray->At(iCut))->GetCutNumber();
434 TString cutstringMeson = ((AliConversionMesonCuts*)fMesonCutArray->At(iCut))->GetCutNumber();
435
436 Int_t collisionSystem = atoi((TString)(((AliConversionCuts*)fCutArray->At(iCut))->GetCutNumber())(0,1));
437 Int_t centMin = atoi((TString)(((AliConversionCuts*)fCutArray->At(iCut))->GetCutNumber())(1,1));
438 Int_t centMax = atoi((TString)(((AliConversionCuts*)fCutArray->At(iCut))->GetCutNumber())(2,1));
439
440 if(collisionSystem == 1 || collisionSystem == 2 ||
441 collisionSystem == 5 || collisionSystem == 8 ||
442 collisionSystem == 9){
443 centMin = centMin*10;
444 centMax = centMax*10;
445 if(centMax ==0 && centMax!=centMin) centMax=100;
446 } else if(collisionSystem == 3 || collisionSystem == 6){
447 centMin = centMin*5;
448 centMax = centMax*5;
449 } else if(collisionSystem == 4 || collisionSystem == 7){
450 centMin = ((centMin*5)+45);
451 centMax = ((centMax*5)+45);
452 }
453
454 fBackList[iCut] = new TList();
455 fBackList[iCut]->SetName(Form("%s_%s_%s Back histograms",cutstring.Data(),cutstringCalo.Data(), cutstringMeson.Data()));
456 fBackList[iCut]->SetOwner(kTRUE);
457 fCutFolder[iCut]->Add(fBackList[iCut]);
458
459 fSparseMotherBackInvMassPtZM[iCut] = new THnSparseF("Back_Back_InvMass_Pt_z_m","Back_Back_InvMass_Pt_z_m",nDim,nBins,xMin,xMax);
460 fBackList[iCut]->Add(fSparseMotherBackInvMassPtZM[iCut]);
461
462 fMotherList[iCut] = new TList();
463 fMotherList[iCut]->SetName(Form("%s_%s_%s Mother histograms",cutstring.Data(),cutstringCalo.Data(),cutstringMeson.Data()));
464 fMotherList[iCut]->SetOwner(kTRUE);
465 fCutFolder[iCut]->Add(fMotherList[iCut]);
466
467 fSparseMotherInvMassPtZM[iCut] = new THnSparseF("Back_Mother_InvMass_Pt_z_m","Back_Mother_InvMass_Pt_z_m",nDim,nBins,xMin,xMax);
468 fMotherList[iCut]->Add(fSparseMotherInvMassPtZM[iCut]);
469
470 if(((AliConversionMesonCuts*)fMesonCutArray->At(iCut))->BackgroundHandlerType() == 0){
471 fBGHandler[iCut] = new AliGammaConversionAODBGHandler(
472 collisionSystem,centMin,centMax,
473 ((AliConversionMesonCuts*)fMesonCutArray->At(iCut))->GetNumberOfBGEvents(),
474 ((AliConversionMesonCuts*)fMesonCutArray->At(iCut))->UseTrackMultiplicity());
475 fBGClusHandler[iCut] = new AliGammaConversionAODBGHandler(
476 collisionSystem,centMin,centMax,
477 ((AliConversionMesonCuts*)fMesonCutArray->At(iCut))->GetNumberOfBGEvents(),
478 ((AliConversionMesonCuts*)fMesonCutArray->At(iCut))->UseTrackMultiplicity());
479 fBGHandlerRP[iCut] = NULL;
480 } else{
481 fBGHandlerRP[iCut] = new AliConversionAODBGHandlerRP(
482 ((AliConversionCuts*)fCutArray->At(fiCut))->IsHeavyIon(),
483 ((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->UseTrackMultiplicity(),
484 ((AliConversionMesonCuts*)fMesonCutArray->At(iCut))->GetNumberOfBGEvents());
485 fBGClusHandlerRP[iCut] = new AliConversionAODBGHandlerRP(
486 ((AliConversionCuts*)fCutArray->At(fiCut))->IsHeavyIon(),
487 ((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->UseTrackMultiplicity(),
488 ((AliConversionMesonCuts*)fMesonCutArray->At(iCut))->GetNumberOfBGEvents());
489 fBGHandler[iCut] = NULL;
490 }
491 }
492 }
493}
494//________________________________________________________________________
495void AliAnalysisTaskGammaConvCalo::UserCreateOutputObjects(){
496
497 // Create histograms
498 if(fOutputContainer != NULL){
499 delete fOutputContainer;
500 fOutputContainer = NULL;
501 }
502 if(fOutputContainer == NULL){
503 fOutputContainer = new TList();
504 fOutputContainer->SetOwner(kTRUE);
505 }
506
507 // Array of current cut's gammas
508 fGammaCandidates = new TList();
509 fClusterCandidates = new TList();
510
511 fCutFolder = new TList*[fnCuts];
512 fESDList = new TList*[fnCuts];
513 fBackList = new TList*[fnCuts];
514 fMotherList = new TList*[fnCuts];
515 fHistoNEvents = new TH1I*[fnCuts];
516 fHistoNGoodESDTracks = new TH1I*[fnCuts];
517 fHistoNGammaCandidates = new TH1I*[fnCuts];
518 fHistoNGoodESDTracksVsNGammaCanditates = new TH2F*[fnCuts];
519 fHistoNV0Tracks = new TH1I*[fnCuts];
520 fProfileEtaShift = new TProfile*[fnCuts];
521 fHistoConvGammaPt = new TH1F*[fnCuts];
522
523 if (fDoPhotonQA == 2){
524 fPhotonDCAList = new TList*[fnCuts];
525 fTreeConvGammaPtDcazCat = new TTree*[fnCuts];
526 }
527 if (fDoPhotonQA > 0){
528 fHistoConvGammaR = new TH1F*[fnCuts];
529 fHistoConvGammaEta = new TH1F*[fnCuts];
530 }
531
532 if(fDoMesonAnalysis){
533 fHistoMotherInvMassPt = new TH2F*[fnCuts];
534 fHistoMotherBackInvMassPt = new TH2F*[fnCuts];
535 fHistoMotherInvMassEalpha = new TH2F*[fnCuts];
536 if (fDoMesonQA == 2){
537 fMesonDCAList = new TList*[fnCuts];
538 fTreeMesonsInvMassPtDcazMinDcazMaxFlag = new TTree*[fnCuts];
539 }
540 if (fDoMesonQA > 0){
541 fHistoMotherPi0PtY = new TH2F*[fnCuts];
542 fHistoMotherEtaPtY = new TH2F*[fnCuts];
543 fHistoMotherPi0PtAlpha = new TH2F*[fnCuts];
544 fHistoMotherEtaPtAlpha = new TH2F*[fnCuts];
545 fHistoMotherPi0PtOpenAngle = new TH2F*[fnCuts];
546 fHistoMotherEtaPtOpenAngle = new TH2F*[fnCuts];
547 }
548 }
549 fTagOutputList = new TList*[fnCuts];
550
551 fHistoConvGammaUntagged = new TH1F*[fnCuts];
552 fHistoConvGammaTagged = new TH1F*[fnCuts];
553 fHistoConvGammaPi0Tagged = new TH1F*[fnCuts];
554 fHistoConvGammaEtaTagged = new TH1F*[fnCuts];
555 fHistoPhotonPairAll = new TH2F*[fnCuts];
556 fHistoPhotonPairAllGam = new TH2F*[fnCuts];
557
558 fHistoClusGammaPt = new TH1F*[fnCuts];
559
560 for(Int_t iCut = 0; iCut<fnCuts;iCut++){
561
562 TString cutstring = ((AliConversionCuts*)fCutArray->At(iCut))->GetCutNumber();
563 TString cutstringCalo = ((AliCaloPhotonCuts*)fClusterCutArray->At(iCut))->GetCutNumber();
564 TString cutstringMeson = "NoMesonCut";
565 if(fDoMesonAnalysis)cutstringMeson = ((AliConversionMesonCuts*)fMesonCutArray->At(iCut))->GetCutNumber();
566
567 fCutFolder[iCut] = new TList();
568 fCutFolder[iCut]->SetName(Form("Cut Number %s_%s_%s",cutstring.Data(),cutstringCalo.Data(),cutstringMeson.Data()));
569 fCutFolder[iCut]->SetOwner(kTRUE);
570 fOutputContainer->Add(fCutFolder[iCut]);
571 fESDList[iCut] = new TList();
572 fESDList[iCut]->SetName(Form("%s_%s_%s ESD histograms",cutstring.Data(),cutstringCalo.Data(),cutstringMeson.Data()));
573 fESDList[iCut]->SetOwner(kTRUE);
574 fCutFolder[iCut]->Add(fESDList[iCut]);
575
576 fHistoNEvents[iCut] = new TH1I("NEvents","NEvents",9,-0.5,8.5);
577 fHistoNEvents[iCut]->GetXaxis()->SetBinLabel(1,"Accepted");
578 fHistoNEvents[iCut]->GetXaxis()->SetBinLabel(2,"Centrality");
579 fHistoNEvents[iCut]->GetXaxis()->SetBinLabel(3,"Missing MC");
580 if (((AliConversionCuts*)fCutArray->At(iCut))->IsSpecialTrigger() == 4 ){
581 TString TriggerNames = "Not Trigger: ";
582 TriggerNames = TriggerNames+ ( (AliConversionCuts*)fCutArray->At(iCut))->GetSpecialTriggerName();
583 fHistoNEvents[iCut]->GetXaxis()->SetBinLabel(4,TriggerNames.Data());
584 } else {
585 fHistoNEvents[iCut]->GetXaxis()->SetBinLabel(4,"Trigger");
586 }
587 fHistoNEvents[iCut]->GetXaxis()->SetBinLabel(5,"Vertex Z");
588 fHistoNEvents[iCut]->GetXaxis()->SetBinLabel(6,"Cont. Vertex");
589 fHistoNEvents[iCut]->GetXaxis()->SetBinLabel(7,"Pile-Up");
590 fHistoNEvents[iCut]->GetXaxis()->SetBinLabel(8,"no SDD");
591 fHistoNEvents[iCut]->GetXaxis()->SetBinLabel(9,"no V0AND");
592 fESDList[iCut]->Add(fHistoNEvents[iCut]);
593
594 if(fIsHeavyIon == 1) fHistoNGoodESDTracks[iCut] = new TH1I("GoodESDTracks","GoodESDTracks",4000,0,4000);
595 else if(fIsHeavyIon == 2) fHistoNGoodESDTracks[iCut] = new TH1I("GoodESDTracks","GoodESDTracks",400,0,400);
596 else fHistoNGoodESDTracks[iCut] = new TH1I("GoodESDTracks","GoodESDTracks",200,0,200);
597 fESDList[iCut]->Add(fHistoNGoodESDTracks[iCut]);
598 if(fIsHeavyIon == 1) fHistoNGammaCandidates[iCut] = new TH1I("GammaCandidates","GammaCandidates",100,0,100);
599 else if(fIsHeavyIon == 2) fHistoNGammaCandidates[iCut] = new TH1I("GammaCandidates","GammaCandidates",50,0,50);
600 else fHistoNGammaCandidates[iCut] = new TH1I("GammaCandidates","GammaCandidates",50,0,50);
601 fESDList[iCut]->Add(fHistoNGammaCandidates[iCut]);
602 if(fIsHeavyIon == 1) fHistoNGoodESDTracksVsNGammaCanditates[iCut] = new TH2F("GoodESDTracksVsGammaCandidates","GoodESDTracksVsGammaCandidates",4000,0,4000,100,0,100);
603 else if(fIsHeavyIon == 2) fHistoNGoodESDTracksVsNGammaCanditates[iCut] = new TH2F("GoodESDTracksVsGammaCandidates","GoodESDTracksVsGammaCandidates",400,0,400,50,0,50);
604 else fHistoNGoodESDTracksVsNGammaCanditates[iCut] = new TH2F("GoodESDTracksVsGammaCandidates","GoodESDTracksVsGammaCandidates",200,0,200,50,0,50);
605 fESDList[iCut]->Add(fHistoNGoodESDTracksVsNGammaCanditates[iCut]);
606
607
608 if(fIsHeavyIon == 1) fHistoNV0Tracks[iCut] = new TH1I("V0 Multiplicity","V0 Multiplicity",30000,0,30000);
609 else if(fIsHeavyIon == 2) fHistoNV0Tracks[iCut] = new TH1I("V0 Multiplicity","V0 Multiplicity",2500,0,2500);
610 else fHistoNV0Tracks[iCut] = new TH1I("V0 Multiplicity","V0 Multiplicity",1500,0,1500);
611 fESDList[iCut]->Add(fHistoNV0Tracks[iCut]);
612 fProfileEtaShift[iCut] = new TProfile("Eta Shift","Eta Shift",1, -0.5,0.5);
613 fESDList[iCut]->Add(fProfileEtaShift[iCut]);
614 fHistoConvGammaPt[iCut] = new TH1F("ESD_ConvGamma_Pt","ESD_ConvGamma_Pt",250,0,25);
615 fESDList[iCut]->Add(fHistoConvGammaPt[iCut]);
616
617 if (fDoPhotonQA == 2){
618 fPhotonDCAList[iCut] = new TList();
619 fPhotonDCAList[iCut]->SetName(Form("%s_%s_%s Photon DCA tree",cutstring.Data(),cutstringCalo.Data(),cutstringMeson.Data()));
620 fPhotonDCAList[iCut]->SetOwner(kTRUE);
621 fCutFolder[iCut]->Add(fPhotonDCAList[iCut]);
622
623 fTreeConvGammaPtDcazCat[iCut] = new TTree("ESD_ConvGamma_Pt_Dcaz_R_Eta","ESD_ConvGamma_Pt_Dcaz_R_Eta_Cat");
624 fTreeConvGammaPtDcazCat[iCut]->Branch("Pt",&fPtGamma,"fPtGamma/F");
625 fTreeConvGammaPtDcazCat[iCut]->Branch("DcaZPhoton",&fDCAzPhoton,"fDCAzPhoton/F");
626 // fTreeConvGammaPtDcazCat[iCut]->Branch("R",&fRConvPhoton,"fRConvPhoton/F");
627 // fTreeConvGammaPtDcazCat[iCut]->Branch("Eta",&fEtaPhoton,"fEtaPhoton/F");
628
629 fTreeConvGammaPtDcazCat[iCut]->Branch("cat",&fCharCatPhoton,"fCharCatPhoton/b");
630 if(fIsMC){
631 fTreeConvGammaPtDcazCat[iCut]->Branch("photonMCInfo",&fCharPhotonMCInfo,"fCharPhotonMCInfo/b");
632 }
633 fPhotonDCAList[iCut]->Add(fTreeConvGammaPtDcazCat[iCut]);
634 }
635
636 if (fDoPhotonQA > 0){
637 fHistoConvGammaR[iCut] = new TH1F("ESD_ConvGamma_R","ESD_ConvGamma_R",800,0,200);
638 fESDList[iCut]->Add(fHistoConvGammaR[iCut]);
639 fHistoConvGammaEta[iCut] = new TH1F("ESD_ConvGamma_Eta","ESD_ConvGamma_Eta",2000,-2,2);
640 fESDList[iCut]->Add(fHistoConvGammaEta[iCut]);
641 }
642
643 fTagOutputList[iCut] = new TList();
644 fTagOutputList[iCut]->SetName(Form("%s_%s_%s Tagging Output",cutstring.Data(),cutstringCalo.Data(),cutstringMeson.Data()));
645 fTagOutputList[iCut]->SetOwner(1);
646 fCutFolder[iCut]->Add(fTagOutputList[iCut]);
647
648 const Int_t nptbins = 200;
649 const Double_t ptmin = 0.;
650 const Double_t ptmax = 20.;
651
652 const Int_t nmbins = 180;
653 const Double_t mmin = 0.;
654 const Double_t mmax = 0.9;
655
656 // photon candidates
657 // this is maybe not necessary ...
658
659 fHistoConvGammaUntagged[iCut] = new TH1F("ConvGammaUntagged","",nptbins,ptmin,ptmax);
660 fHistoConvGammaUntagged[iCut]->SetXTitle("p_{T} (GeV/c)");
661 fTagOutputList[iCut]->Add(fHistoConvGammaUntagged[iCut]);
662
663 fHistoConvGammaTagged[iCut] = new TH1F("ConvGammaTagged","",nptbins,ptmin,ptmax);
664 fHistoConvGammaTagged[iCut]->SetXTitle("p_{T} (GeV/c)");
665 fTagOutputList[iCut]->Add(fHistoConvGammaTagged[iCut]);
666
667 fHistoConvGammaPi0Tagged[iCut] = new TH1F("ConvGammaPi0Tagged","",nptbins,ptmin,ptmax);
668 fHistoConvGammaPi0Tagged[iCut]->SetXTitle("p_{T} (GeV/c)");
669 fTagOutputList[iCut]->Add(fHistoConvGammaPi0Tagged[iCut]);
670
671 fHistoConvGammaEtaTagged[iCut] = new TH1F("ConvGammaEtaTagged","",nptbins,ptmin,ptmax);
672 fHistoConvGammaEtaTagged[iCut]->SetXTitle("p_{T} (GeV/c)");
673 fTagOutputList[iCut]->Add(fHistoConvGammaEtaTagged[iCut]);
674
675 // pairs
676 fHistoPhotonPairAll[iCut] = new TH2F("PhotonPairAll","",nmbins,mmin,mmax,nptbins,ptmin,ptmax);
677 fHistoPhotonPairAll[iCut]->SetXTitle("M_{inv} (GeV/cc)");
678 fHistoPhotonPairAll[iCut]->SetYTitle("p_{T} (GeV/c)");
679 fTagOutputList[iCut]->Add(fHistoPhotonPairAll[iCut]);
680
681 fHistoPhotonPairAllGam[iCut] = new TH2F("PhotonPairAllGammaConvPt","",nmbins,mmin,mmax,nptbins,ptmin,ptmax);
682 fHistoPhotonPairAllGam[iCut]->SetXTitle("M_{inv} (GeV/cc)");
683 fHistoPhotonPairAllGam[iCut]->SetYTitle("#gamma^{conv} p_{T} (GeV/c)");
684 fTagOutputList[iCut]->Add(fHistoPhotonPairAllGam[iCut]);
685
686 fHistoClusGammaPt[iCut] = new TH1F("ClusGamma_Pt","ClusGamma_Pt",250,0,25);
687 fTagOutputList[iCut]->Add(fHistoClusGammaPt[iCut]);
688
689
690 if(fDoMesonAnalysis){
691 fHistoMotherInvMassPt[iCut] = new TH2F("ESD_Mother_InvMass_Pt","ESD_Mother_InvMass_Pt",800,0,0.8,250,0,25);
692 fESDList[iCut]->Add(fHistoMotherInvMassPt[iCut]);
693 fHistoMotherBackInvMassPt[iCut] = new TH2F("ESD_Background_InvMass_Pt","ESD_Background_InvMass_Pt",800,0,0.8,250,0,25);
694 fESDList[iCut]->Add(fHistoMotherBackInvMassPt[iCut]);
695 fHistoMotherInvMassEalpha[iCut] = new TH2F("ESD_Mother_InvMass_vs_E_alpha","ESD_Mother_InvMass_vs_E_alpha",800,0,0.8,250,0,25);
696 fESDList[iCut]->Add(fHistoMotherInvMassEalpha[iCut]);
697 if (fDoMesonQA == 2){
698 fMesonDCAList[iCut] = new TList();
699 fMesonDCAList[iCut]->SetName(Form("%s_%s_%s Meson DCA tree",cutstring.Data(),cutstringCalo.Data(),cutstringMeson.Data()));
700 fMesonDCAList[iCut]->SetOwner(kTRUE);
701 fCutFolder[iCut]->Add(fMesonDCAList[iCut]);
702
703 fTreeMesonsInvMassPtDcazMinDcazMaxFlag[iCut] = new TTree("ESD_Mesons_InvMass_Pt_DcazMin_DcazMax_Flag","ESD_Mesons_InvMass_Pt_DcazMin_DcazMax_Flag");
704 fTreeMesonsInvMassPtDcazMinDcazMaxFlag[iCut]->Branch("InvMass",&fInvMass,"fInvMass/F");
705 fTreeMesonsInvMassPtDcazMinDcazMaxFlag[iCut]->Branch("Pt",&fPt,"fPt/F");
706 fTreeMesonsInvMassPtDcazMinDcazMaxFlag[iCut]->Branch("DcaZMin",&fDCAzGammaMin,"fDCAzGammaMin/F");
707 fTreeMesonsInvMassPtDcazMinDcazMaxFlag[iCut]->Branch("DcaZMax",&fDCAzGammaMax,"fDCAzGammaMax/F");
708 fTreeMesonsInvMassPtDcazMinDcazMaxFlag[iCut]->Branch("kind",&fCharFlag,"fCharFlag/b");
709 if(fIsMC){
710 fTreeMesonsInvMassPtDcazMinDcazMaxFlag[iCut]->Branch("mesonMCInfo",&fCharMesonMCInfo,"fCharMesonMCInfo/b");
711 }
712 fMesonDCAList[iCut]->Add(fTreeMesonsInvMassPtDcazMinDcazMaxFlag[iCut]);
713
714 }
715 if (fDoMesonQA > 0 ){
716 fHistoMotherPi0PtY[iCut] = new TH2F("ESD_MotherPi0_Pt_Y","ESD_MotherPi0_Pt_Y",150,0.03,15.,150,-1.5,1.5);
717 SetLogBinningXTH2(fHistoMotherPi0PtY[iCut]);
718 fESDList[iCut]->Add(fHistoMotherPi0PtY[iCut]);
719 fHistoMotherEtaPtY[iCut] = new TH2F("ESD_MotherEta_Pt_Y","ESD_MotherEta_Pt_Y",150,0.03,15.,150,-1.5,1.5);
720 SetLogBinningXTH2(fHistoMotherEtaPtY[iCut]);
721 fESDList[iCut]->Add(fHistoMotherEtaPtY[iCut]);
722 fHistoMotherPi0PtAlpha[iCut] = new TH2F("ESD_MotherPi0_Pt_Alpha","ESD_MotherPi0_Pt_Alpha",150,0.03,15.,100,0,1);
723 SetLogBinningXTH2(fHistoMotherPi0PtAlpha[iCut]);
724 fESDList[iCut]->Add(fHistoMotherPi0PtAlpha[iCut]);
725 fHistoMotherEtaPtAlpha[iCut] = new TH2F("ESD_MotherEta_Pt_Alpha","ESD_MotherEta_Pt_Alpha",150,0.03,15.,100,0,1);
726 SetLogBinningXTH2(fHistoMotherEtaPtAlpha[iCut]);
727 fESDList[iCut]->Add(fHistoMotherEtaPtAlpha[iCut]);
728 fHistoMotherPi0PtOpenAngle[iCut] = new TH2F("ESD_MotherPi0_Pt_OpenAngle","ESD_MotherPi0_Pt_OpenAngle",150,0.03,15.,200,0,2*TMath::Pi());
729 SetLogBinningXTH2(fHistoMotherPi0PtOpenAngle[iCut]);
730 fESDList[iCut]->Add(fHistoMotherPi0PtOpenAngle[iCut]);
731 fHistoMotherEtaPtOpenAngle[iCut] = new TH2F("ESD_MotherEta_Pt_OpenAngle","ESD_MotherEta_Pt_OpenAngle",150,0.03,15.,200,0,2*TMath::Pi());
732 SetLogBinningXTH2(fHistoMotherEtaPtOpenAngle[iCut]);
733 fESDList[iCut]->Add(fHistoMotherEtaPtOpenAngle[iCut]);
734 }
735 }
736 }
737 if(fDoMesonAnalysis){
738 InitBack(); // Init Background Handler
739 }
740
741 if(fIsMC){
742 // MC Histogramms
743 fMCList = new TList*[fnCuts];
744 // True Histogramms
745 fTrueList = new TList*[fnCuts];
746 // Selected Header List
747 fHeaderNameList = new TList*[fnCuts];
748 fHistoMCHeaders = new TH1I*[fnCuts];
749 fHistoMCAllGammaPt = new TH1F*[fnCuts];
750 fHistoMCDecayGammaPi0Pt = new TH1F*[fnCuts];
751 fHistoMCDecayGammaRhoPt = new TH1F*[fnCuts];
752 fHistoMCDecayGammaEtaPt = new TH1F*[fnCuts];
753 fHistoMCDecayGammaOmegaPt = new TH1F*[fnCuts];
754 fHistoMCDecayGammaEtapPt = new TH1F*[fnCuts];
755 fHistoMCDecayGammaPhiPt = new TH1F*[fnCuts];
756 fHistoMCDecayGammaSigmaPt = new TH1F*[fnCuts];
757 fHistoMCConvGammaPt = new TH1F*[fnCuts];
758 fHistoTrueConvGammaPt = new TH1F*[fnCuts];
759 fHistoTrueConvPi0GammaPt = new TH1F*[fnCuts];
760
761 fHistoCombinatorialPt = new TH2F*[fnCuts];
762 fHistoTruePrimaryConvGammaPt = new TH1F*[fnCuts];
763 fHistoTruePrimaryConvGammaESDPtMCPt = new TH2F*[fnCuts];
764 fHistoTrueSecondaryConvGammaPt = new TH1F*[fnCuts];
765 fHistoTrueSecondaryConvGammaFromXFromK0sPt = new TH1F*[fnCuts];
766 fHistoTrueSecondaryConvGammaFromXFromLambdaPt = new TH1F*[fnCuts];
767
768 fHistoTrueDalitzPsiPairDeltaPhi= new TH2F*[fnCuts];
769 fHistoTrueGammaPsiPairDeltaPhi= new TH2F*[fnCuts];
770
771 fHistoTrueClusGammaPt = new TH1F*[fnCuts];
772 fHistoTruePrimaryClusGammaPt = new TH1F*[fnCuts];
773 fHistoTruePrimaryClusGammaESDPtMCPt = new TH2F*[fnCuts];
774
775 if (fDoPhotonQA > 0){
776 fHistoMCConvGammaR = new TH1F*[fnCuts];
777 fHistoMCConvGammaEta = new TH1F*[fnCuts];
778 fHistoTrueConvGammaEta = new TH1F*[fnCuts];
779 }
780
781 if(fDoMesonAnalysis){
782 fHistoMCPi0Pt = new TH1F*[fnCuts];
783 fHistoMCPi0WOWeightPt = new TH1F*[fnCuts];
784 fHistoMCEtaPt = new TH1F*[fnCuts];
785 fHistoMCEtaWOWeightPt = new TH1F*[fnCuts];
786 fHistoMCPi0InAccPt = new TH1F*[fnCuts];
787 fHistoMCEtaInAccPt = new TH1F*[fnCuts];
788
789 fHistoTrueMotherInvMassPt = new TH2F*[fnCuts];
790 fHistoTruePrimaryMotherInvMassPt = new TH2F*[fnCuts];
791 fHistoTruePrimaryMotherW0WeightingInvMassPt = new TH2F*[fnCuts];
792 fProfileTruePrimaryMotherWeightsInvMassPt = new TProfile2D*[fnCuts];
793 fHistoTrueSecondaryMotherInvMassPt = new TH2F*[fnCuts];
794 fHistoTrueSecondaryMotherFromK0sInvMassPt = new TH2F*[fnCuts];
795 fHistoTrueSecondaryMotherFromEtaInvMassPt = new TH2F*[fnCuts];
796 fHistoTrueSecondaryMotherFromLambdaInvMassPt = new TH2F*[fnCuts];
797 fHistoTrueMotherDalitzInvMassPt = new TH2F*[fnCuts];
798 if (fDoMesonQA > 0){
799 fHistoMCPi0PtY = new TH2F*[fnCuts];
800 fHistoMCEtaPtY = new TH2F*[fnCuts];
801 fHistoMCK0sPt = new TH1F*[fnCuts];
802 fHistoMCK0sWOWeightPt = new TH1F*[fnCuts];
803 fHistoMCK0sPtY = new TH2F*[fnCuts];
804 fHistoMCSecPi0PtvsSource= new TH2F*[fnCuts];
805 fHistoMCSecPi0Source = new TH1F*[fnCuts];
806 fHistoMCSecEtaPt = new TH1F*[fnCuts];
807 fHistoMCSecEtaSource = new TH1F*[fnCuts];
808 fHistoTruePrimaryPi0MCPtResolPt = new TH2F*[fnCuts];
809 fHistoTruePrimaryEtaMCPtResolPt = new TH2F*[fnCuts];
810 fHistoTrueK0sWithPi0DaughterMCPt = new TH1F*[fnCuts];
811 fHistoTrueEtaWithPi0DaughterMCPt = new TH1F*[fnCuts];
812 fHistoTrueLambdaWithPi0DaughterMCPt = new TH1F*[fnCuts];
813 fHistoTrueBckGGInvMassPt = new TH2F*[fnCuts];
814 fHistoTrueBckContInvMassPt = new TH2F*[fnCuts];
815 fHistoTruePi0PtY = new TH2F*[fnCuts];
816 fHistoTrueEtaPtY = new TH2F*[fnCuts];
817 fHistoTruePi0PtAlpha = new TH2F*[fnCuts];
818 fHistoTrueEtaPtAlpha = new TH2F*[fnCuts];
819 fHistoTruePi0PtOpenAngle = new TH2F*[fnCuts];
820 fHistoTrueEtaPtOpenAngle = new TH2F*[fnCuts];
821 }
822 }
823
824
825
826 for(Int_t iCut = 0; iCut<fnCuts;iCut++){
827 TString cutstring = ((AliConversionCuts*)fCutArray->At(iCut))->GetCutNumber();
828 TString cutstringCalo = ((AliCaloPhotonCuts*)fClusterCutArray->At(iCut))->GetCutNumber();
829 TString cutstringMeson = "NoMesonCut";
830 if(fDoMesonAnalysis)cutstringMeson = ((AliConversionMesonCuts*)fMesonCutArray->At(iCut))->GetCutNumber();
831
832 fMCList[iCut] = new TList();
833 fMCList[iCut]->SetName(Form("%s_%s_%s MC histograms",cutstring.Data(),cutstringCalo.Data(),cutstringMeson.Data()));
834 fMCList[iCut]->SetOwner(kTRUE);
835 fCutFolder[iCut]->Add(fMCList[iCut]);
836 fHistoMCHeaders[iCut] = new TH1I("MC_Headers","MC_Headers",20,0,20);
837 fMCList[iCut]->Add(fHistoMCHeaders[iCut]);
838 fHistoMCAllGammaPt[iCut] = new TH1F("MC_AllGamma_Pt","MC_AllGamma_Pt",250,0,25);
839 fMCList[iCut]->Add(fHistoMCAllGammaPt[iCut]);
840 fHistoMCDecayGammaPi0Pt[iCut] = new TH1F("MC_DecayGammaPi0_Pt","MC_DecayGammaPi0_Pt",250,0,25);
841 fMCList[iCut]->Add(fHistoMCDecayGammaPi0Pt[iCut]);
842 fHistoMCDecayGammaRhoPt[iCut] = new TH1F("MC_DecayGammaRho_Pt","MC_DecayGammaRho_Pt",250,0,25);
843 fMCList[iCut]->Add(fHistoMCDecayGammaRhoPt[iCut]);
844 fHistoMCDecayGammaEtaPt[iCut] = new TH1F("MC_DecayGammaEta_Pt","MC_DecayGammaEta_Pt",250,0,25);
845 fMCList[iCut]->Add(fHistoMCDecayGammaEtaPt[iCut]);
846 fHistoMCDecayGammaOmegaPt[iCut] = new TH1F("MC_DecayGammaOmega_Pt","MC_DecayGammaOmmega_Pt",250,0,25);
847 fMCList[iCut]->Add(fHistoMCDecayGammaOmegaPt[iCut]);
848 fHistoMCDecayGammaEtapPt[iCut] = new TH1F("MC_DecayGammaEtap_Pt","MC_DecayGammaEtap_Pt",250,0,25);
849 fMCList[iCut]->Add(fHistoMCDecayGammaEtapPt[iCut]);
850 fHistoMCDecayGammaPhiPt[iCut] = new TH1F("MC_DecayGammaPhi_Pt","MC_DecayGammaPhi_Pt",250,0,25);
851 fMCList[iCut]->Add(fHistoMCDecayGammaPhiPt[iCut]);
852 fHistoMCDecayGammaSigmaPt[iCut] = new TH1F("MC_DecayGammaSigma_Pt","MC_DecayGammaSigma_Pt",250,0,25);
853 fMCList[iCut]->Add(fHistoMCDecayGammaSigmaPt[iCut]);
854 fHistoMCConvGammaPt[iCut] = new TH1F("MC_ConvGamma_Pt","MC_ConvGamma_Pt",250,0,25);
855 fMCList[iCut]->Add(fHistoMCConvGammaPt[iCut]);
856
857 if (fDoPhotonQA > 0){
858 fHistoMCConvGammaR[iCut] = new TH1F("MC_ConvGamma_R","MC_ConvGamma_R",800,0,200);
859 fMCList[iCut]->Add(fHistoMCConvGammaR[iCut]);
860 fHistoMCConvGammaEta[iCut] = new TH1F("MC_ConvGamma_Eta","MC_ConvGamma_Eta",2000,-2,2);
861 fMCList[iCut]->Add(fHistoMCConvGammaEta[iCut]);
862 }
863
864 if(fDoMesonAnalysis){
865 fHistoMCPi0Pt[iCut] = new TH1F("MC_Pi0_Pt","MC_Pi0_Pt",250,0,25);
866 fHistoMCPi0Pt[iCut]->Sumw2();
867 fMCList[iCut]->Add(fHistoMCPi0Pt[iCut]);
868 fHistoMCPi0WOWeightPt[iCut] = new TH1F("MC_Pi0_WOWeights_Pt","MC_Pi0_WOWeights_Pt",250,0,25);
869 fHistoMCPi0WOWeightPt[iCut]->Sumw2();
870 fMCList[iCut]->Add(fHistoMCPi0WOWeightPt[iCut]);
871
872 fHistoMCEtaPt[iCut] = new TH1F("MC_Eta_Pt","MC_Eta_Pt",250,0,25);
873 fHistoMCEtaPt[iCut]->Sumw2();
874 fMCList[iCut]->Add(fHistoMCEtaPt[iCut]);
875 fHistoMCEtaWOWeightPt[iCut] = new TH1F("MC_Eta_WOWeights_Pt","MC_Eta_WOWeights_Pt",250,0,25);
876 fHistoMCEtaWOWeightPt[iCut]->Sumw2();
877 fMCList[iCut]->Add(fHistoMCEtaWOWeightPt[iCut]);
878 fHistoMCPi0InAccPt[iCut] = new TH1F("MC_Pi0InAcc_Pt","MC_Pi0InAcc_Pt",250,0,25);
879 fHistoMCPi0InAccPt[iCut]->Sumw2();
880 fMCList[iCut]->Add(fHistoMCPi0InAccPt[iCut]);
881 fHistoMCEtaInAccPt[iCut] = new TH1F("MC_EtaInAcc_Pt","MC_EtaInAcc_Pt",250,0,25);
882 fHistoMCEtaInAccPt[iCut]->Sumw2();
883 fMCList[iCut]->Add(fHistoMCEtaInAccPt[iCut]);
884 if (fDoMesonQA > 0){
885 fHistoMCPi0PtY[iCut] = new TH2F("MC_Pi0_Pt_Y","MC_Pi0_Pt_Y",150,0.03,15.,150,-1.5,1.5);
886 fHistoMCPi0PtY[iCut]->Sumw2();
887 SetLogBinningXTH2(fHistoMCPi0PtY[iCut]);
888 fMCList[iCut]->Add(fHistoMCPi0PtY[iCut]);
889 fHistoMCEtaPtY[iCut] = new TH2F("MC_Eta_Pt_Y","MC_Eta_Pt_Y",150,0.03,15.,150,-1.5,1.5);
890 fHistoMCEtaPtY[iCut]->Sumw2();
891 SetLogBinningXTH2(fHistoMCEtaPtY[iCut]);
892 fMCList[iCut]->Add(fHistoMCEtaPtY[iCut]);
893 fHistoMCK0sPt[iCut] = new TH1F("MC_K0s_Pt","MC_K0s_Pt",150,0,15);
894 fHistoMCK0sPt[iCut]->Sumw2();
895 fMCList[iCut]->Add(fHistoMCK0sPt[iCut]);
896 fHistoMCK0sWOWeightPt[iCut] = new TH1F("MC_K0s_WOWeights_Pt","MC_K0s_WOWeights_Pt",150,0,15);
897 fHistoMCK0sWOWeightPt[iCut]->Sumw2();
898 fMCList[iCut]->Add(fHistoMCK0sWOWeightPt[iCut]);
899 fHistoMCK0sPtY[iCut] = new TH2F("MC_K0s_Pt_Y","MC_K0s_Pt_Y",150,0.03,15.,150,-1.5,1.5);
900 fHistoMCK0sPtY[iCut]->Sumw2();
901 SetLogBinningXTH2(fHistoMCK0sPtY[iCut]);
902 fMCList[iCut]->Add(fHistoMCK0sPtY[iCut]);
903
904 fHistoMCSecPi0Source[iCut] = new TH1F("MC_SecPi0_Source","MC_SecPi0_Source",5000,0.,5000);
905 fMCList[iCut]->Add(fHistoMCSecPi0Source[iCut]);
906 fHistoMCSecEtaSource[iCut] = new TH1F("MC_SecEta_Source","MC_SecEta_Source",5000,0,5000);
907 fMCList[iCut]->Add(fHistoMCSecEtaSource[iCut]);
908 fHistoMCSecPi0PtvsSource[iCut] = new TH2F("MC_SecPi0_Pt_Source","MC_SecPi0_Pt_Source",250,0.0,25.,16,-0.5,15.5);
909 fHistoMCSecPi0PtvsSource[iCut]->Sumw2();
910 fMCList[iCut]->Add(fHistoMCSecPi0PtvsSource[iCut]);
911 fHistoMCSecEtaPt[iCut] = new TH1F("MC_SecEta_Pt","MC_SecEta_Pt",250,0,25);
912 fHistoMCSecEtaPt[iCut]->Sumw2();
913 fMCList[iCut]->Add(fHistoMCSecEtaPt[iCut]);
914 }
915
916 }
917 fTrueList[iCut] = new TList();
918 fTrueList[iCut]->SetName(Form("%s_%s_%s True histograms",cutstring.Data(),cutstringCalo.Data(),cutstringMeson.Data()));
919 fTrueList[iCut]->SetOwner(kTRUE);
920 fCutFolder[iCut]->Add(fTrueList[iCut]);
921
922 fHistoTrueConvGammaPt[iCut] = new TH1F("ESD_TrueConvGamma_Pt","ESD_TrueConvGamma_Pt",250,0,25);
923 fTrueList[iCut]->Add(fHistoTrueConvGammaPt[iCut]);
924
925 fHistoTrueConvPi0GammaPt[iCut] = new TH1F("ESD_TrueConvGamma_Pt","ESD_TrueConvGamma_Pt",250,0,25);
926 fTrueList[iCut]->Add(fHistoTrueConvPi0GammaPt[iCut]);
927
928 fHistoCombinatorialPt[iCut] = new TH2F("ESD_TrueCombinatorial_Pt","ESD_TrueCombinatorial_Pt",250,0,25,16,-0.5,15.5);
929 fHistoCombinatorialPt[iCut]->GetYaxis()->SetBinLabel( 1,"Elec+Elec");
930 fHistoCombinatorialPt[iCut]->GetYaxis()->SetBinLabel( 2,"Elec+Pion");
931 fHistoCombinatorialPt[iCut]->GetYaxis()->SetBinLabel( 3,"Elec+Kaon");
932 fHistoCombinatorialPt[iCut]->GetYaxis()->SetBinLabel( 4,"Elec+Proton");
933 fHistoCombinatorialPt[iCut]->GetYaxis()->SetBinLabel( 5,"Elec+Muon");
934 fHistoCombinatorialPt[iCut]->GetYaxis()->SetBinLabel( 6,"Pion+Pion");
935 fHistoCombinatorialPt[iCut]->GetYaxis()->SetBinLabel( 7,"Pion+Kaon");
936 fHistoCombinatorialPt[iCut]->GetYaxis()->SetBinLabel( 8,"Pion+Proton");
937 fHistoCombinatorialPt[iCut]->GetYaxis()->SetBinLabel( 9,"Pion+Muon");
938 fHistoCombinatorialPt[iCut]->GetYaxis()->SetBinLabel(10,"Kaon+Kaon");
939 fHistoCombinatorialPt[iCut]->GetYaxis()->SetBinLabel(11,"Kaon+Proton");
940 fHistoCombinatorialPt[iCut]->GetYaxis()->SetBinLabel(12,"Kaon+Muon");
941 fHistoCombinatorialPt[iCut]->GetYaxis()->SetBinLabel(13,"Proton+Proton");
942 fHistoCombinatorialPt[iCut]->GetYaxis()->SetBinLabel(14,"Proton+Muon");
943 fHistoCombinatorialPt[iCut]->GetYaxis()->SetBinLabel(15,"Muon+Muon");
944 fHistoCombinatorialPt[iCut]->GetYaxis()->SetBinLabel(16,"Rest");
945 fTrueList[iCut]->Add(fHistoCombinatorialPt[iCut]);
946 fHistoTruePrimaryConvGammaPt[iCut] = new TH1F("ESD_TruePrimaryConvGamma_Pt","ESD_TruePrimaryConvGamma_Pt",250,0,25);
947 fTrueList[iCut]->Add(fHistoTruePrimaryConvGammaPt[iCut]);
948 fHistoTrueSecondaryConvGammaPt[iCut] = new TH1F("ESD_TrueSecondaryConvGamma_Pt","ESD_TrueSecondaryConvGamma_Pt",250,0,25);
949 fTrueList[iCut]->Add(fHistoTrueSecondaryConvGammaPt[iCut]);
950
951 fHistoTrueSecondaryConvGammaFromXFromK0sPt[iCut] = new TH1F("ESD_TrueSecondaryConvGammaFromXFromK0s_Pt", "ESD_TrueSecondaryConvGammaFromXFromK0s_Pt",250,0,25);
952 fTrueList[iCut]->Add(fHistoTrueSecondaryConvGammaFromXFromK0sPt[iCut]);
953 fHistoTrueSecondaryConvGammaFromXFromLambdaPt[iCut] = new TH1F("ESD_TrueSecondaryConvGammaFromXFromLambda_Pt", "ESD_TrueSecondaryConvGammaFromXFromLambda_Pt",250,0,25);
954 fTrueList[iCut]->Add(fHistoTrueSecondaryConvGammaFromXFromLambdaPt[iCut]);
955
956 fHistoTrueDalitzPsiPairDeltaPhi[iCut] = new TH2F("ESD_TrueDalitzPsiPairDeltaPhi_Pt", "ESD_TrueDalitzPsiPairDeltaPhi_Pt",100,-0.5,2,100,-0.5,0.5);
957 fTrueList[iCut]->Add(fHistoTrueDalitzPsiPairDeltaPhi[iCut]);
958
959 fHistoTrueGammaPsiPairDeltaPhi[iCut] = new TH2F("ESD_TrueGammaPsiPairDeltaPhi_Pt", "ESD_TrueGammaPsiPairDeltaPhi_Pt",100,-0.5,2,100,-0.5,0.5);
960 fTrueList[iCut]->Add(fHistoTrueGammaPsiPairDeltaPhi[iCut]);
961
962 fHistoTruePrimaryConvGammaESDPtMCPt[iCut] = new TH2F("ESD_TruePrimaryConvGammaESD_PtMCPt", "ESD_TruePrimaryConvGammaESD_PtMCPt",250,0,25,250,0,25);
963 fTrueList[iCut]->Add(fHistoTruePrimaryConvGammaESDPtMCPt[iCut]);
964
965 fHistoTrueClusGammaPt[iCut] = new TH1F("TrueClusGamma_Pt","ESD_TruePrimaryClusGamma_Pt",250,0,25);
966 fTagOutputList[iCut]->Add(fHistoTrueClusGammaPt[iCut]);
967 fHistoTruePrimaryClusGammaPt[iCut] = new TH1F("TruePrimaryClusGamma_Pt","ESD_TruePrimaryClusGamma_Pt",250,0,25);
968 fTagOutputList[iCut]->Add(fHistoTruePrimaryClusGammaPt[iCut]);
969 fHistoTruePrimaryClusGammaESDPtMCPt[iCut] = new TH2F("TruePrimaryClusGamma_Pt_MCPt","ESD_TruePrimaryClusGamma_MCPt",250,0,25,250,0,25);
970 fTagOutputList[iCut]->Add(fHistoTruePrimaryClusGammaESDPtMCPt[iCut]);
971
972 if(fDoMesonAnalysis){
973 fHistoTrueMotherInvMassPt[iCut] = new TH2F("ESD_TrueMother_InvMass_Pt","ESD_TrueMother_InvMass_Pt",800,0,0.8,250,0,25);
974 fTrueList[iCut]->Add(fHistoTrueMotherInvMassPt[iCut]);
975 fHistoTruePrimaryMotherInvMassPt[iCut] = new TH2F("ESD_TruePrimaryMother_InvMass_Pt", "ESD_TruePrimaryMother_InvMass_Pt", 800,0,0.8,250,0,25);
976 fHistoTruePrimaryMotherInvMassPt[iCut]->Sumw2();
977 fTrueList[iCut]->Add(fHistoTruePrimaryMotherInvMassPt[iCut]);
978 fHistoTruePrimaryMotherW0WeightingInvMassPt[iCut] = new TH2F("ESD_TruePrimaryMotherW0Weights_InvMass_Pt", "ESD_TruePrimaryMotherW0Weights_InvMass_Pt", 800,0,0.8,250,0,25);
979 fHistoTruePrimaryMotherW0WeightingInvMassPt[iCut]->Sumw2();
980 fTrueList[iCut]->Add(fHistoTruePrimaryMotherW0WeightingInvMassPt[iCut]);
981 fProfileTruePrimaryMotherWeightsInvMassPt[iCut] = new TProfile2D("ESD_TruePrimaryMotherWeights_InvMass_Pt", "ESD_TruePrimaryMotherWeights_InvMass_Pt", 800,0,0.8,250,0,25);
982 fProfileTruePrimaryMotherWeightsInvMassPt[iCut]->Sumw2();
983 fTrueList[iCut]->Add(fProfileTruePrimaryMotherWeightsInvMassPt[iCut]);
984 fHistoTrueSecondaryMotherInvMassPt[iCut] = new TH2F("ESD_TrueSecondaryMother_InvMass_Pt", "ESD_TrueSecondaryMother_InvMass_Pt", 800,0,0.8,250,0,25);
985 fHistoTrueSecondaryMotherInvMassPt[iCut]->Sumw2();
986 fTrueList[iCut]->Add(fHistoTrueSecondaryMotherInvMassPt[iCut]);
987 fHistoTrueSecondaryMotherFromK0sInvMassPt[iCut] = new TH2F("ESD_TrueSecondaryMotherFromK0s_InvMass_Pt","ESD_TrueSecondaryMotherFromK0s_InvMass_Pt",800,0,0.8,250,0,25);
988 fHistoTrueSecondaryMotherFromK0sInvMassPt[iCut]->Sumw2();
989 fTrueList[iCut]->Add(fHistoTrueSecondaryMotherFromK0sInvMassPt[iCut]);
990 fHistoTrueSecondaryMotherFromEtaInvMassPt[iCut] = new TH2F("ESD_TrueSecondaryMotherFromEta_InvMass_Pt","ESD_TrueSecondaryMotherFromEta_InvMass_Pt",800,0,0.8,250,0,25);
991 fTrueList[iCut]->Add(fHistoTrueSecondaryMotherFromEtaInvMassPt[iCut]);
992 fHistoTrueSecondaryMotherFromLambdaInvMassPt[iCut] = new TH2F("ESD_TrueSecondaryMotherFromLambda_InvMass_Pt","ESD_TrueSecondaryMotherFromLambda_InvMass_Pt",800,0,0.8,250,0,25);
993 fTrueList[iCut]->Add(fHistoTrueSecondaryMotherFromLambdaInvMassPt[iCut]);
994 fHistoTrueMotherDalitzInvMassPt[iCut] = new TH2F("ESD_TrueDalitz_InvMass_Pt","ESD_TrueDalitz_InvMass_Pt",800,0,0.8,250,0,25);
995 fTrueList[iCut]->Add(fHistoTrueMotherDalitzInvMassPt[iCut]);
996 if (fDoMesonQA > 0){
997 fHistoTruePrimaryPi0MCPtResolPt[iCut] = new TH2F("ESD_TruePrimaryPi0_MCPt_ResolPt","ESD_TruePrimaryPi0_ResolPt_MCPt",500,0.03,25,1000,-1.,1.);
998 fHistoTruePrimaryPi0MCPtResolPt[iCut]->Sumw2();
999 SetLogBinningXTH2(fHistoTruePrimaryPi0MCPtResolPt[iCut]);
1000 fTrueList[iCut]->Add(fHistoTruePrimaryPi0MCPtResolPt[iCut]);
1001 fHistoTruePrimaryEtaMCPtResolPt[iCut] = new TH2F("ESD_TruePrimaryEta_MCPt_ResolPt","ESD_TruePrimaryEta_ResolPt_MCPt",500,0.03,25,1000,-1.,1.);
1002 fHistoTruePrimaryEtaMCPtResolPt[iCut]->Sumw2();
1003 SetLogBinningXTH2(fHistoTruePrimaryEtaMCPtResolPt[iCut]);
1004 fTrueList[iCut]->Add(fHistoTruePrimaryEtaMCPtResolPt[iCut]);
1005 fHistoTrueBckGGInvMassPt[iCut] = new TH2F("ESD_TrueBckGG_InvMass_Pt","ESD_TrueBckGG_InvMass_Pt",800,0,0.8,250,0,25);
1006 fTrueList[iCut]->Add(fHistoTrueBckGGInvMassPt[iCut]);
1007 fHistoTrueBckContInvMassPt[iCut] = new TH2F("ESD_TrueBckCont_InvMass_Pt","ESD_TrueBckCont_InvMass_Pt",800,0,0.8,250,0,25);
1008 fTrueList[iCut]->Add(fHistoTrueBckContInvMassPt[iCut]);
1009 fHistoTrueK0sWithPi0DaughterMCPt[iCut] = new TH1F("ESD_TrueK0sWithPi0Daughter_MCPt","ESD_TrueK0sWithPi0Daughter_MCPt",250,0,25);
1010 fTrueList[iCut]->Add(fHistoTrueK0sWithPi0DaughterMCPt[iCut]);
1011 fHistoTrueEtaWithPi0DaughterMCPt[iCut] = new TH1F("ESD_TrueEtaWithPi0Daughter_MCPt","ESD_TrueEtaWithPi0Daughter_MCPt",250,0,25);
1012 fTrueList[iCut]->Add(fHistoTrueEtaWithPi0DaughterMCPt[iCut]);
1013 fHistoTrueLambdaWithPi0DaughterMCPt[iCut] = new TH1F("ESD_TrueLambdaWithPi0Daughter_MCPt","ESD_TrueLambdaWithPi0Daughter_MCPt",250,0,25);
1014 fTrueList[iCut]->Add(fHistoTrueLambdaWithPi0DaughterMCPt[iCut]);
1015
1016 fHistoTruePi0PtY[iCut] = new TH2F("ESD_TruePi0_Pt_Y","ESD_TruePi0_Pt_Y",150,0.03,15.,150,-1.5,1.5);
1017 SetLogBinningXTH2(fHistoTruePi0PtY[iCut]);
1018 fTrueList[iCut]->Add(fHistoTruePi0PtY[iCut]);
1019 fHistoTrueEtaPtY[iCut] = new TH2F("ESD_TrueEta_Pt_Y","ESD_TrueEta_Pt_Y",150,0.03,15.,150,-1.5,1.5);
1020 SetLogBinningXTH2(fHistoTrueEtaPtY[iCut]);
1021 fTrueList[iCut]->Add(fHistoTrueEtaPtY[iCut]);
1022 fHistoTruePi0PtAlpha[iCut] = new TH2F("ESD_TruePi0_Pt_Alpha","ESD_TruePi0_Pt_Alpha",150,0.03,15.,100,0,1);
1023 SetLogBinningXTH2(fHistoTruePi0PtAlpha[iCut]);
1024 fTrueList[iCut]->Add(fHistoTruePi0PtAlpha[iCut]);
1025 fHistoTrueEtaPtAlpha[iCut] = new TH2F("ESD_TrueEta_Pt_Alpha","ESD_TrueEta_Pt_Alpha",150,0.03,15.,100,0,1);
1026 SetLogBinningXTH2(fHistoTrueEtaPtAlpha[iCut]);
1027 fTrueList[iCut]->Add(fHistoTrueEtaPtAlpha[iCut]);
1028
1029 fHistoTruePi0PtOpenAngle[iCut] = new TH2F("ESD_TruePi0_Pt_OpenAngle","ESD_TruePi0_Pt_OpenAngle",150,0.03,15.,200,0,2*TMath::Pi());
1030 SetLogBinningXTH2(fHistoTruePi0PtOpenAngle[iCut]);
1031 fTrueList[iCut]->Add(fHistoTruePi0PtOpenAngle[iCut]);
1032 fHistoTrueEtaPtOpenAngle[iCut] = new TH2F("ESD_TrueEta_Pt_OpenAngle","ESD_TrueEta_Pt_OpenAngle",150,0.03,15.,200,0,2*TMath::Pi());
1033 SetLogBinningXTH2(fHistoTrueEtaPtOpenAngle[iCut]);
1034 fTrueList[iCut]->Add(fHistoTrueEtaPtOpenAngle[iCut]);
1035
1036 fHistoTrueConvGammaEta[iCut] = new TH1F("ESD_TrueConvGamma_Eta","ESD_TrueConvGamma_Eta",2000,-2,2);
1037 fTrueList[iCut]->Add(fHistoTrueConvGammaEta[iCut]);
1038
1039 }
1040 }
1041 }
1042 }
1043
1044 fV0Reader=(AliV0ReaderV1*)AliAnalysisManager::GetAnalysisManager()->GetTask("V0ReaderV1");
1045 if(!fV0Reader){printf("Error: No V0 Reader");return;} // GetV0Reader
1046
1047 if(fV0Reader)
1048 if((AliConversionCuts*)fV0Reader->GetConversionCuts())
1049 if(((AliConversionCuts*)fV0Reader->GetConversionCuts())->GetCutHistograms())
1050 fOutputContainer->Add(((AliConversionCuts*)fV0Reader->GetConversionCuts())->GetCutHistograms());
1051
1052 for(Int_t iCut = 0; iCut<fnCuts;iCut++){
1053 if(!((AliConversionCuts*)fCutArray->At(iCut))) continue;
1054 if(((AliConversionCuts*)fCutArray->At(iCut))->GetCutHistograms()){
1055 fCutFolder[iCut]->Add(((AliConversionCuts*)fCutArray->At(iCut))->GetCutHistograms());
1056 }
1057 if(!((AliCaloPhotonCuts*)fClusterCutArray->At(iCut))) continue;
1058 if(((AliCaloPhotonCuts*)fClusterCutArray->At(iCut))->GetCutHistograms()){
1059 fCutFolder[iCut]->Add(((AliCaloPhotonCuts*)fClusterCutArray->At(iCut))->GetCutHistograms());
1060 }
1061 if(fDoMesonAnalysis){
1062 if(!((AliConversionMesonCuts*)fMesonCutArray->At(iCut))) continue;
1063 if(((AliConversionMesonCuts*)fMesonCutArray->At(iCut))->GetCutHistograms()){
1064 fCutFolder[iCut]->Add(((AliConversionMesonCuts*)fMesonCutArray->At(iCut))->GetCutHistograms());
1065 }
1066 }
1067 }
1068 PostData(1, fOutputContainer);
1069}
1070//_____________________________________________________________________________
1071Bool_t AliAnalysisTaskGammaConvCalo::Notify()
1072{
1073 for(Int_t iCut = 0; iCut<fnCuts;iCut++){
1074 if(!((AliConversionCuts*)fCutArray->At(iCut))->GetDoEtaShift()){
1075 fProfileEtaShift[iCut]->Fill(0.,(((AliConversionCuts*)fCutArray->At(iCut))->GetEtaShift()));
1076 continue; // No Eta Shift requested, continue
1077 }
1078 if(((AliConversionCuts*)fCutArray->At(iCut))->GetEtaShift() == 0.0){ // Eta Shift requested but not set, get shift automatically
1079 ((AliConversionCuts*)fCutArray->At(iCut))->GetCorrectEtaShiftFromPeriod(fV0Reader->GetPeriodName());
1080 fProfileEtaShift[iCut]->Fill(0.,(((AliConversionCuts*)fCutArray->At(iCut))->GetEtaShift()));
1081 ((AliConversionCuts*)fCutArray->At(iCut))->DoEtaShift(kFALSE); // Eta Shift Set, make sure that it is called only once
1082 continue;
1083 }
1084 else{
1085 printf(" Gamma Conversion Task %s :: Eta Shift Manually Set to %f \n\n",
1086 (((AliConversionCuts*)fCutArray->At(iCut))->GetCutNumber()).Data(),((AliConversionCuts*)fCutArray->At(iCut))->GetEtaShift());
1087 fProfileEtaShift[iCut]->Fill(0.,(((AliConversionCuts*)fCutArray->At(iCut))->GetEtaShift()));
1088 ((AliConversionCuts*)fCutArray->At(iCut))->DoEtaShift(kFALSE); // Eta Shift Set, make sure that it is called only once
1089 }
1090 }
1091
1092 return kTRUE;
1093}
1094//_____________________________________________________________________________
1095void AliAnalysisTaskGammaConvCalo::UserExec(Option_t *)
1096{
1097 //
1098 // Called for each event
1099 //
1100 Int_t eventQuality = ((AliConversionCuts*)fV0Reader->GetConversionCuts())->GetEventQuality();
1101 if(eventQuality == 2 || eventQuality == 3){// Event Not Accepted due to MC event missing or wrong trigger for V0ReaderV1
1102 for(Int_t iCut = 0; iCut<fnCuts; iCut++){
1103 fHistoNEvents[iCut]->Fill(eventQuality);
1104 }
1105 return;
1106 }
1107
1108 if(fIsMC) fMCEvent = MCEvent();
1109 if(fMCEvent == NULL) fIsMC = kFALSE;
1110
1111 fInputEvent = InputEvent();
1112
1113 if(fIsMC && fInputEvent->IsA()==AliESDEvent::Class()){
1114 fMCStack = fMCEvent->Stack();
1115 if(fMCStack == NULL) fIsMC = kFALSE;
1116 }
1117
1118 fReaderGammas = fV0Reader->GetReconstructedGammas(); // Gammas from default Cut
1119
1120 // ------------------- BeginEvent ----------------------------
1121
1122 AliEventplane *EventPlane = fInputEvent->GetEventplane();
1123 if(fIsHeavyIon ==1)fEventPlaneAngle = EventPlane->GetEventplane("V0",fInputEvent,2);
1124 else fEventPlaneAngle=0.0;
1125
1126 if(fIsMC && fInputEvent->IsA()==AliAODEvent::Class() && !(fV0Reader->AreAODsRelabeled())){
1127 RelabelAODPhotonCandidates(kTRUE); // In case of AODMC relabeling MC
1128 fV0Reader->RelabelAODs(kTRUE);
1129 }
1130
1131
1132 for(Int_t iCut = 0; iCut<fnCuts; iCut++){
1133
1134 fiCut = iCut;
1135 Int_t eventNotAccepted = ((AliConversionCuts*)fCutArray->At(iCut))->IsEventAcceptedByConversionCut(fV0Reader->GetConversionCuts(),fInputEvent,fMCEvent,fIsHeavyIon);
1136
1137 if(eventNotAccepted){
1138 // cout << "event rejected due to wrong trigger: " <<eventNotAccepted << endl;
1139 fHistoNEvents[iCut]->Fill(eventNotAccepted); // Check Centrality, PileUp, SDD and V0AND --> Not Accepted => eventQuality = 1
1140 continue;
1141 }
1142
1143 if(eventQuality != 0){// Event Not Accepted
1144 //cout << "event rejected due to: " <<eventQuality << endl;
1145 fHistoNEvents[iCut]->Fill(eventQuality);
1146 continue;
1147 }
1148
1149 fHistoNEvents[iCut]->Fill(eventQuality); // Should be 0 here
1150 fHistoNGoodESDTracks[iCut]->Fill(fV0Reader->GetNumberOfPrimaryTracks());
1151 if(((AliConversionCuts*)fCutArray->At(iCut))->IsHeavyIon() == 2) fHistoNV0Tracks[iCut]->Fill(fInputEvent->GetVZEROData()->GetMTotV0A());
1152 else fHistoNV0Tracks[iCut]->Fill(fInputEvent->GetVZEROData()->GetMTotV0A()+fInputEvent->GetVZEROData()->GetMTotV0C());
1153
1154 if(fIsMC){
1155 // Process MC Particle
1156 if(((AliConversionCuts*)fCutArray->At(iCut))->GetSignalRejection() != 0){
1157 if(fInputEvent->IsA()==AliESDEvent::Class()){
1158 ((AliConversionCuts*)fCutArray->At(iCut))->GetNotRejectedParticles(((AliConversionCuts*)fCutArray->At(iCut))->GetSignalRejection(),
1159 ((AliConversionCuts*)fCutArray->At(iCut))->GetAcceptedHeader(),
1160 fMCEvent);
1161 }
1162 else if(fInputEvent->IsA()==AliAODEvent::Class()){
1163 ((AliConversionCuts*)fCutArray->At(iCut))->GetNotRejectedParticles(((AliConversionCuts*)fCutArray->At(iCut))->GetSignalRejection(),
1164 ((AliConversionCuts*)fCutArray->At(iCut))->GetAcceptedHeader(),
1165 fInputEvent);
1166 }
1167
1168 if(((AliConversionCuts*)fCutArray->At(iCut))->GetAcceptedHeader()){
1169 for(Int_t i = 0;i<(((AliConversionCuts*)fCutArray->At(iCut))->GetAcceptedHeader())->GetEntries();i++){
1170 TString nameBin= fHistoMCHeaders[iCut]->GetXaxis()->GetBinLabel(i+1);
1171 if (nameBin.CompareTo("")== 0){
1172 TString nameHeader = ((TObjString*)((TList*)((AliConversionCuts*)fCutArray->At(iCut))
1173 ->GetAcceptedHeader())->At(i))->GetString();
1174 fHistoMCHeaders[iCut]->GetXaxis()->SetBinLabel(i+1,nameHeader.Data());
1175 }
1176 }
1177 }
1178 }
1179 }
1180 if(fIsMC){
1181 if(fInputEvent->IsA()==AliESDEvent::Class())
1182 ProcessMCParticles();
1183 if(fInputEvent->IsA()==AliAODEvent::Class())
1184 ProcessAODMCParticles();
1185 }
1186
1187 // it is in the loop to have the same conversion cut string (used also for MC stuff that should be same for V0 and Cluster)
1188 ProcessClusters(); // process calo clusters
1189 ProcessPhotonCandidates(); // Process this cuts gammas
1190
1191 fHistoNGammaCandidates[iCut]->Fill(fGammaCandidates->GetEntries());
1192 fHistoNGoodESDTracksVsNGammaCanditates[iCut]->Fill(fV0Reader->GetNumberOfPrimaryTracks(),fGammaCandidates->GetEntries());
1193 if(fDoMesonAnalysis){ // Meson Analysis
1194 if(((AliConversionMesonCuts*)fMesonCutArray->At(iCut))->UseMCPSmearing() && fIsMC){
1195 fUnsmearedPx = new Double_t[fGammaCandidates->GetEntries()]; // Store unsmeared Momenta
1196 fUnsmearedPy = new Double_t[fGammaCandidates->GetEntries()];
1197 fUnsmearedPz = new Double_t[fGammaCandidates->GetEntries()];
1198 fUnsmearedE = new Double_t[fGammaCandidates->GetEntries()];
1199
1200 for(Int_t gamma=0;gamma<fGammaCandidates->GetEntries();gamma++){ // Smear the AODPhotons in MC
1201 fUnsmearedPx[gamma] = ((AliAODConversionPhoton*)fGammaCandidates->At(gamma))->Px();
1202 fUnsmearedPy[gamma] = ((AliAODConversionPhoton*)fGammaCandidates->At(gamma))->Py();
1203 fUnsmearedPz[gamma] = ((AliAODConversionPhoton*)fGammaCandidates->At(gamma))->Pz();
1204 fUnsmearedE[gamma] = ((AliAODConversionPhoton*)fGammaCandidates->At(gamma))->E();
1205 ((AliConversionMesonCuts*)fMesonCutArray->At(iCut))->SmearParticle(dynamic_cast<AliAODConversionPhoton*>(fGammaCandidates->At(gamma)));
1206 }
1207 }
1208
1209 PhotonTagging(); // tag PCM photons with calorimeter
1210
1211 CalculatePi0Candidates(); // Combine Gammas from conversion and from calo
1212
1213 if(((AliConversionMesonCuts*)fMesonCutArray->At(iCut))->DoBGCalculation()){
1214 if(((AliConversionMesonCuts*)fMesonCutArray->At(iCut))->BackgroundHandlerType() == 0){
1215 CalculateBackground(); // Combinatorial Background
1216 UpdateEventByEventData(); // Store Event for mixed Events
1217 }
1218 else{
1219 CalculateBackgroundRP(); // Combinatorial Background
1220 fBGHandlerRP[iCut]->AddEvent(fGammaCandidates,fInputEvent); // Store Event for mixed Events
1221 fBGClusHandlerRP[iCut]->AddEvent(fClusterCandidates,fInputEvent); // Store Event for mixed Events
1222 }
1223 }
1224 if(((AliConversionMesonCuts*)fMesonCutArray->At(iCut))->UseMCPSmearing() && fIsMC){
1225 for(Int_t gamma=0;gamma<fGammaCandidates->GetEntries();gamma++){ // Smear the AODPhotons in MC
1226 ((AliAODConversionPhoton*)fGammaCandidates->At(gamma))->SetPx(fUnsmearedPx[gamma]); // Reset Unsmeared Momenta
1227 ((AliAODConversionPhoton*)fGammaCandidates->At(gamma))->SetPy(fUnsmearedPy[gamma]);
1228 ((AliAODConversionPhoton*)fGammaCandidates->At(gamma))->SetPz(fUnsmearedPz[gamma]);
1229 ((AliAODConversionPhoton*)fGammaCandidates->At(gamma))->SetE(fUnsmearedE[gamma]);
1230 }
1231 delete[] fUnsmearedPx; fUnsmearedPx = 0x0;
1232 delete[] fUnsmearedPy; fUnsmearedPy = 0x0;
1233 delete[] fUnsmearedPz; fUnsmearedPz = 0x0;
1234 delete[] fUnsmearedE; fUnsmearedE = 0x0;
1235 }
1236 }
1237
1238 fGammaCandidates->Clear(); // delete this cuts good gammas
1239 fClusterCandidates->Clear(); // delete cluster candidates
1240 }
1241
1242 if(fIsMC && fInputEvent->IsA()==AliAODEvent::Class() && !(fV0Reader->AreAODsRelabeled())){
1243 RelabelAODPhotonCandidates(kFALSE); // Back to ESDMC Label
1244 fV0Reader->RelabelAODs(kFALSE);
1245 }
1246
1247 PostData(1, fOutputContainer);
1248}
1249
1250//________________________________________________________________________
1251void AliAnalysisTaskGammaConvCalo::ProcessClusters()
1252{
1253
1254 Int_t nclus = 0;
1255 nclus = fInputEvent->GetNumberOfCaloClusters();
1256
1257// cout << nclus << endl;
1258
1259 if(nclus == 0) return;
1260
1261 // vertex
1262 Double_t vertex[3] = {0};
1263 InputEvent()->GetPrimaryVertex()->GetXYZ(vertex);
1264
1265 // Loop over EMCal clusters
1266 for(Int_t i = 0; i < nclus; i++){
1267
1268 AliVCluster* clus = fInputEvent->GetCaloCluster(i);
1269 if (!clus) continue;
1270 if(!((AliCaloPhotonCuts*)fClusterCutArray->At(fiCut))->ClusterIsSelected(clus,fInputEvent,fIsMC)) continue;
1271 // TLorentzvector with cluster
1272 TLorentzVector clusterVector;
1273 clus->GetMomentum(clusterVector,vertex);
1274
1275 TLorentzVector* tmpvec = new TLorentzVector();
1276 tmpvec->SetPxPyPzE(clusterVector.Px(),clusterVector.Py(),clusterVector.Pz(),clusterVector.E());
1277
1278 // convert to AODConversionPhoton
1279 AliAODConversionPhoton *PhotonCandidate=new AliAODConversionPhoton(tmpvec);
1280 if(!PhotonCandidate) continue;
1281
1282 // get MC label
1283 Int_t mclab[2] = {0,0};
1284 if(fIsMC){
1285 mclab[0] = clus->GetLabel();
1286 mclab[1] = clus->GetLabel();
1287// cout << "mclabels: " << mclab[0] << " " << mclab[1] << endl;
1288 }
1289
1290 PhotonCandidate->SetMCLabel(mclab);
1291
1292 // fIsFromMBHeader = kTRUE;
1293 // if(fIsMC){
1294 // Int_t isPosFromMBHeader
1295 // = ((AliConversionCuts*)fCutArray->At(fiCut))->IsParticleFromBGEvent(mclab, fMCStack, fInputEvent);
1296 // if(isPosFromMBHeader == 0 && ((AliConversionCuts*)fCutArray->At(fiCut))->GetSignalRejection() != 3) continue;
1297 //
1298 // Int_t isNegFromMBHeader
1299 // = ((AliConversionCuts*)fCutArray->At(fiCut))->IsParticleFromBGEvent(mclab, fMCStack, fInputEvent);
1300 // if(isNegFromMBHeader == 0 && ((AliConversionCuts*)fCutArray->At(fiCut))->GetSignalRejection() != 3) continue;
1301 //
1302 // if( (isNegFromMBHeader+isPosFromMBHeader) != 4) fIsFromMBHeader = kFALSE;
1303 // }
1304
1305 fHistoClusGammaPt[fiCut]->Fill(PhotonCandidate->Pt());
1306 fClusterCandidates->Add(PhotonCandidate); // if no second loop is required add to events good gammas
1307
1308 // if(fIsFromMBHeader){
1309 // fHistoConvGammaPt[fiCut]->Fill(PhotonCandidate->Pt());
1310 // if (fDoPhotonQA > 0){
1311 // fHistoConvGammaR[fiCut]->Fill(PhotonCandidate->GetConversionRadius());
1312 // fHistoConvGammaEta[fiCut]->Fill(PhotonCandidate->Eta());
1313 // }
1314 // }
1315
1316 if(fIsMC){
1317 ProcessTrueClusterCandidates(PhotonCandidate);
1318 }
1319
1320 // if (fIsFromMBHeader && fDoPhotonQA == 2){
1321 // if (fIsHeavyIon == 1 && PhotonCandidate->Pt() > 0.399 && PhotonCandidate->Pt() < 12.){
1322 // fPtGamma = PhotonCandidate->Pt();
1323 // fDCAzPhoton = PhotonCandidate->GetDCAzToPrimVtx();
1324 // fRConvPhoton = PhotonCandidate->GetConversionRadius();
1325 // fEtaPhoton = PhotonCandidate->GetPhotonEta();
1326 // fCharCatPhoton = PhotonCandidate->GetPhotonQuality();
1327 // fTreeConvGammaPtDcazCat[fiCut]->Fill();
1328 // } else if ( PhotonCandidate->Pt() > 0.299 && PhotonCandidate->Pt() < 16.){
1329 // fPtGamma = PhotonCandidate->Pt();
1330 // fDCAzPhoton = PhotonCandidate->GetDCAzToPrimVtx();
1331 // fRConvPhoton = PhotonCandidate->GetConversionRadius();
1332 // fEtaPhoton = PhotonCandidate->GetPhotonEta();
1333 // fCharCatPhoton = PhotonCandidate->GetPhotonQuality();
1334 // fTreeConvGammaPtDcazCat[fiCut]->Fill();
1335 // }
1336 // }
1337
1338 }
1339
1340}
1341
1342//________________________________________________________________________
1343void AliAnalysisTaskGammaConvCalo::ProcessTrueClusterCandidates(AliAODConversionPhoton *TruePhotonCandidate)
1344{
1345
1346 TParticle *Photon = TruePhotonCandidate->GetMCParticle(fMCStack);
1347
1348 if(Photon == NULL){
1349 // cout << "no photon" << endl;
1350 return;
1351 }
1352
1353 if(Photon->GetPdgCode() != 22){
1354 // cout << "pdg code: " << Photon->GetPdgCode() << endl;
1355 return; // Particle is no Photon
1356 }
1357
1358 // True Photon
1359 if(fIsFromMBHeader){
1360 fHistoTrueClusGammaPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1361 if (fDoPhotonQA > 0) fHistoTrueConvGammaEta[fiCut]->Fill(TruePhotonCandidate->Eta());
1362 }
1363
1364 if(Photon->GetMother(0) <= fMCStack->GetNprimary()){
1365 // Count just primary MC Gammas as true --> For Ratio esdtruegamma / mcconvgamma
1366 if(fIsFromMBHeader){
1367 fHistoTruePrimaryClusGammaPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1368 fHistoTruePrimaryClusGammaESDPtMCPt[fiCut]->Fill(TruePhotonCandidate->Pt(),Photon->Pt()); // Allways Filled
1369 }
1370 }
1371
1372 // maybe later
1373 // else{
1374 // if(fIsFromMBHeader){
1375 // fCharPhotonMCInfo = 2;
1376 // fHistoTrueSecondaryConvGammaPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1377 // if(fMCStack->Particle(Photon->GetMother(0))->GetMother(0) > -1 &&
1378 // fMCStack->Particle(fMCStack->Particle(Photon->GetMother(0))->GetMother(0))->GetPdgCode() == 3122){
1379 // fHistoTrueSecondaryConvGammaFromXFromLambdaPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1380 // fCharPhotonMCInfo = 5;
1381 // }
1382 // if(fMCStack->Particle(Photon->GetMother(0))->GetMother(0) > -1 &&
1383 // fMCStack->Particle(fMCStack->Particle(Photon->GetMother(0))->GetMother(0))->GetPdgCode() == 310){
1384 // fHistoTrueSecondaryConvGammaFromXFromK0sPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1385 // fCharPhotonMCInfo = 4;
1386 // }
1387 // if(fMCStack->Particle(Photon->GetMother(0))->GetMother(0) > -1 &&
1388 // fMCStack->Particle(fMCStack->Particle(Photon->GetMother(0))->GetMother(0))->GetPdgCode() == 221){
1389 // fCharPhotonMCInfo = 3;
1390 // }
1391 // }
1392 // }
1393}
1394
1395
1396//________________________________________________________________________
1397void AliAnalysisTaskGammaConvCalo::ProcessPhotonCandidates()
1398{
1399 Int_t nV0 = 0;
1400 TList *GammaCandidatesStepOne = new TList();
1401 TList *GammaCandidatesStepTwo = new TList();
1402 // Loop over Photon Candidates allocated by ReaderV1
1403 for(Int_t i = 0; i < fReaderGammas->GetEntriesFast(); i++){
1404 AliAODConversionPhoton* PhotonCandidate = (AliAODConversionPhoton*) fReaderGammas->At(i);
1405 if(!PhotonCandidate) continue;
1406 fIsFromMBHeader = kTRUE;
1407 if(fIsMC && ((AliConversionCuts*)fCutArray->At(fiCut))->GetSignalRejection() != 0){
1408 Int_t isPosFromMBHeader
1409 = ((AliConversionCuts*)fCutArray->At(fiCut))->IsParticleFromBGEvent(PhotonCandidate->GetMCLabelPositive(), fMCStack, fInputEvent);
1410 if(isPosFromMBHeader == 0 && ((AliConversionCuts*)fCutArray->At(fiCut))->GetSignalRejection() != 3) continue;
1411 Int_t isNegFromMBHeader
1412 = ((AliConversionCuts*)fCutArray->At(fiCut))->IsParticleFromBGEvent(PhotonCandidate->GetMCLabelNegative(), fMCStack, fInputEvent);
1413 if(isNegFromMBHeader == 0 && ((AliConversionCuts*)fCutArray->At(fiCut))->GetSignalRejection() != 3) continue;
1414
1415 if( (isNegFromMBHeader+isPosFromMBHeader) != 4) fIsFromMBHeader = kFALSE;
1416 }
1417
1418 if(!((AliConversionCuts*)fCutArray->At(fiCut))->PhotonIsSelected(PhotonCandidate,fInputEvent)) continue;
1419 if(!((AliConversionCuts*)fCutArray->At(fiCut))->InPlaneOutOfPlaneCut(PhotonCandidate->GetPhotonPhi(),fEventPlaneAngle)) continue;
1420 if(!((AliConversionCuts*)fCutArray->At(fiCut))->UseElecSharingCut() &&
1421 !((AliConversionCuts*)fCutArray->At(fiCut))->UseToCloseV0sCut()){
1422 fGammaCandidates->Add(PhotonCandidate); // if no second loop is required add to events good gammas
1423
1424 if(fIsFromMBHeader){
1425 fHistoConvGammaPt[fiCut]->Fill(PhotonCandidate->Pt());
1426 if (fDoPhotonQA > 0){
1427 fHistoConvGammaR[fiCut]->Fill(PhotonCandidate->GetConversionRadius());
1428 fHistoConvGammaEta[fiCut]->Fill(PhotonCandidate->Eta());
1429 }
1430 }
1431 if(fIsMC){
1432 if(fInputEvent->IsA()==AliESDEvent::Class())
1433 ProcessTruePhotonCandidates(PhotonCandidate);
1434 if(fInputEvent->IsA()==AliAODEvent::Class())
1435 ProcessTruePhotonCandidatesAOD(PhotonCandidate);
1436 }
1437 if (fIsFromMBHeader && fDoPhotonQA == 2){
1438 if (fIsHeavyIon == 1 && PhotonCandidate->Pt() > 0.399 && PhotonCandidate->Pt() < 12.){
1439 fPtGamma = PhotonCandidate->Pt();
1440 fDCAzPhoton = PhotonCandidate->GetDCAzToPrimVtx();
1441 fRConvPhoton = PhotonCandidate->GetConversionRadius();
1442 fEtaPhoton = PhotonCandidate->GetPhotonEta();
1443 fCharCatPhoton = PhotonCandidate->GetPhotonQuality();
1444 fTreeConvGammaPtDcazCat[fiCut]->Fill();
1445 } else if ( PhotonCandidate->Pt() > 0.299 && PhotonCandidate->Pt() < 16.){
1446 fPtGamma = PhotonCandidate->Pt();
1447 fDCAzPhoton = PhotonCandidate->GetDCAzToPrimVtx();
1448 fRConvPhoton = PhotonCandidate->GetConversionRadius();
1449 fEtaPhoton = PhotonCandidate->GetPhotonEta();
1450 fCharCatPhoton = PhotonCandidate->GetPhotonQuality();
1451 fTreeConvGammaPtDcazCat[fiCut]->Fill();
1452 }
1453 }
1454 } else if(((AliConversionCuts*)fCutArray->At(fiCut))->UseElecSharingCut()){ // if Shared Electron cut is enabled, Fill array, add to step one
1455 ((AliConversionCuts*)fCutArray->At(fiCut))->FillElectonLabelArray(PhotonCandidate,nV0);
1456 nV0++;
1457 GammaCandidatesStepOne->Add(PhotonCandidate);
1458 } else if(!((AliConversionCuts*)fCutArray->At(fiCut))->UseElecSharingCut() &&
1459 ((AliConversionCuts*)fCutArray->At(fiCut))->UseToCloseV0sCut()){ // shared electron is disabled, step one not needed -> step two
1460 GammaCandidatesStepTwo->Add(PhotonCandidate);
1461 }
1462 }
1463 if(((AliConversionCuts*)fCutArray->At(fiCut))->UseElecSharingCut()){
1464 for(Int_t i = 0;i<GammaCandidatesStepOne->GetEntries();i++){
1465 AliAODConversionPhoton *PhotonCandidate= (AliAODConversionPhoton*) GammaCandidatesStepOne->At(i);
1466 if(!PhotonCandidate) continue;
1467 fIsFromMBHeader = kTRUE;
1468 if(fMCStack && ((AliConversionCuts*)fCutArray->At(fiCut))->GetSignalRejection() != 0){
1469 Int_t isPosFromMBHeader
1470 = ((AliConversionCuts*)fCutArray->At(fiCut))->IsParticleFromBGEvent(PhotonCandidate->GetMCLabelPositive(), fMCStack, fInputEvent);
1471 Int_t isNegFromMBHeader
1472 = ((AliConversionCuts*)fCutArray->At(fiCut))->IsParticleFromBGEvent(PhotonCandidate->GetMCLabelNegative(), fMCStack, fInputEvent);
1473 if( (isNegFromMBHeader+isPosFromMBHeader) != 4) fIsFromMBHeader = kFALSE;
1474 }
1475 if(!((AliConversionCuts*)fCutArray->At(fiCut))->RejectSharedElectronV0s(PhotonCandidate,i,GammaCandidatesStepOne->GetEntries())) continue;
1476 if(!((AliConversionCuts*)fCutArray->At(fiCut))->UseToCloseV0sCut()){ // To Colse v0s cut diabled, step two not needed
1477 fGammaCandidates->Add(PhotonCandidate);
1478 if(fIsFromMBHeader){
1479 fHistoConvGammaPt[fiCut]->Fill(PhotonCandidate->Pt());
1480 if (fDoPhotonQA > 0){
1481 fHistoConvGammaR[fiCut]->Fill(PhotonCandidate->GetConversionRadius());
1482 fHistoConvGammaEta[fiCut]->Fill(PhotonCandidate->Eta());
1483 }
1484 }
1485 }
1486 if(fIsMC){
1487 if(fInputEvent->IsA()==AliESDEvent::Class())
1488 ProcessTruePhotonCandidates(PhotonCandidate);
1489 if(fInputEvent->IsA()==AliAODEvent::Class())
1490 ProcessTruePhotonCandidatesAOD(PhotonCandidate);
1491 } else GammaCandidatesStepTwo->Add(PhotonCandidate); // Close v0s cut enabled -> add to list two
1492
1493 if (fIsFromMBHeader && fDoPhotonQA == 2){
1494 if (fIsHeavyIon ==1 && PhotonCandidate->Pt() > 0.399 && PhotonCandidate->Pt() < 12.){
1495 fPtGamma = PhotonCandidate->Pt();
1496 fDCAzPhoton = PhotonCandidate->GetDCAzToPrimVtx();
1497 fRConvPhoton = PhotonCandidate->GetConversionRadius();
1498 fEtaPhoton = PhotonCandidate->GetPhotonEta();
1499 fCharCatPhoton = PhotonCandidate->GetPhotonQuality();
1500 fTreeConvGammaPtDcazCat[fiCut]->Fill();
1501 } else if ( PhotonCandidate->Pt() > 0.299 && PhotonCandidate->Pt() < 16.){
1502 fPtGamma = PhotonCandidate->Pt();
1503 fDCAzPhoton = PhotonCandidate->GetDCAzToPrimVtx();
1504 fRConvPhoton = PhotonCandidate->GetConversionRadius();
1505 fEtaPhoton = PhotonCandidate->GetPhotonEta();
1506 fCharCatPhoton = PhotonCandidate->GetPhotonQuality();
1507 fTreeConvGammaPtDcazCat[fiCut]->Fill();
1508 }
1509 }
1510 }
1511 }
1512 if(((AliConversionCuts*)fCutArray->At(fiCut))->UseToCloseV0sCut()){
1513 for(Int_t i = 0;i<GammaCandidatesStepTwo->GetEntries();i++){
1514 AliAODConversionPhoton* PhotonCandidate = (AliAODConversionPhoton*) GammaCandidatesStepTwo->At(i);
1515 if(!PhotonCandidate) continue;
1516 fIsFromMBHeader = kTRUE;
1517 if(fMCStack && ((AliConversionCuts*)fCutArray->At(fiCut))->GetSignalRejection() != 0){
1518 Int_t isPosFromMBHeader
1519 = ((AliConversionCuts*)fCutArray->At(fiCut))->IsParticleFromBGEvent(PhotonCandidate->GetMCLabelPositive(), fMCStack, fInputEvent);
1520 Int_t isNegFromMBHeader
1521 = ((AliConversionCuts*)fCutArray->At(fiCut))->IsParticleFromBGEvent(PhotonCandidate->GetMCLabelNegative(), fMCStack, fInputEvent);
1522 if( (isNegFromMBHeader+isPosFromMBHeader) != 4) fIsFromMBHeader = kFALSE;
1523 }
1524 if(!((AliConversionCuts*)fCutArray->At(fiCut))->RejectToCloseV0s(PhotonCandidate,GammaCandidatesStepTwo,i)) continue;
1525 fGammaCandidates->Add(PhotonCandidate); // Add gamma to current cut TList
1526 if(fIsFromMBHeader){
1527 fHistoConvGammaPt[fiCut]->Fill(PhotonCandidate->Pt());
1528 if (fDoPhotonQA > 0){
1529 fHistoConvGammaR[fiCut]->Fill(PhotonCandidate->GetConversionRadius());
1530 fHistoConvGammaEta[fiCut]->Fill(PhotonCandidate->Eta());
1531 }
1532 }
1533 if(fIsMC){
1534 if(fInputEvent->IsA()==AliESDEvent::Class())
1535 ProcessTruePhotonCandidates(PhotonCandidate);
1536 if(fInputEvent->IsA()==AliAODEvent::Class())
1537 ProcessTruePhotonCandidatesAOD(PhotonCandidate);
1538 }
1539 if (fIsFromMBHeader){
1540 if (fIsHeavyIon == 1 && PhotonCandidate->Pt() > 0.399 && PhotonCandidate->Pt() < 12.){
1541 fPtGamma = PhotonCandidate->Pt();
1542 fDCAzPhoton = PhotonCandidate->GetDCAzToPrimVtx();
1543 fRConvPhoton = PhotonCandidate->GetConversionRadius();
1544 fEtaPhoton = PhotonCandidate->GetPhotonEta();
1545 fCharCatPhoton = PhotonCandidate->GetPhotonQuality();
1546 fTreeConvGammaPtDcazCat[fiCut]->Fill();
1547 } else if ( PhotonCandidate->Pt() > 0.299 && PhotonCandidate->Pt() < 16.){
1548 fPtGamma = PhotonCandidate->Pt();
1549 fDCAzPhoton = PhotonCandidate->GetDCAzToPrimVtx();
1550 fRConvPhoton = PhotonCandidate->GetConversionRadius();
1551 fEtaPhoton = PhotonCandidate->GetPhotonEta();
1552 fCharCatPhoton = PhotonCandidate->GetPhotonQuality();
1553 fTreeConvGammaPtDcazCat[fiCut]->Fill();
1554 }
1555 }
1556 }
1557 }
1558
1559 delete GammaCandidatesStepOne;
1560 GammaCandidatesStepOne = 0x0;
1561 delete GammaCandidatesStepTwo;
1562 GammaCandidatesStepTwo = 0x0;
1563
1564}
1565//________________________________________________________________________
1566void AliAnalysisTaskGammaConvCalo::ProcessTruePhotonCandidatesAOD(AliAODConversionPhoton *TruePhotonCandidate)
1567{
1568
1569 Double_t magField = fInputEvent->GetMagneticField();
1570 if( magField < 0.0 ){
1571 magField = 1.0;
1572 }
1573 else {
1574 magField = -1.0;
1575 }
1576
1577 TClonesArray *AODMCTrackArray = dynamic_cast<TClonesArray*>(fInputEvent->FindListObject(AliAODMCParticle::StdBranchName()));
1578 AliAODMCParticle *posDaughter = (AliAODMCParticle*) AODMCTrackArray->At(TruePhotonCandidate->GetMCLabelPositive());
1579 AliAODMCParticle *negDaughter = (AliAODMCParticle*) AODMCTrackArray->At(TruePhotonCandidate->GetMCLabelNegative());
1580 fCharPhotonMCInfo = 0;
1581
1582 if(posDaughter == NULL || negDaughter == NULL) return; // One particle does not exist
1583 Int_t pdgCode[2] = {abs(posDaughter->GetPdgCode()),abs(negDaughter->GetPdgCode())};
1584
1585 if(posDaughter->GetMother() != negDaughter->GetMother()){
1586 FillPhotonCombinatorialBackgroundHist(TruePhotonCandidate, pdgCode);
1587 fCharPhotonMCInfo = 1;
1588 return;
1589 }
1590 else if(posDaughter->GetMother() == -1){
1591 FillPhotonCombinatorialBackgroundHist(TruePhotonCandidate, pdgCode);
1592 fCharPhotonMCInfo = 1;
1593 return;
1594 }
1595
1596 if(pdgCode[0]!=11 || pdgCode[1]!=11){
1597 fCharPhotonMCInfo = 1;
1598 return; //One Particle is not a electron
1599 }
1600 if(posDaughter->GetPdgCode()==negDaughter->GetPdgCode()){
1601 fCharPhotonMCInfo = 1;
1602 return; // Same Charge
1603 }
1604
1605 AliAODMCParticle *Photon = (AliAODMCParticle*) AODMCTrackArray->At(posDaughter->GetMother());
1606 AliVTrack * electronCandidate = ((AliConversionCuts*)fCutArray->At(fiCut))->GetTrack(fInputEvent,TruePhotonCandidate->GetTrackLabelNegative() );
1607 AliVTrack * positronCandidate = ((AliConversionCuts*)fCutArray->At(fiCut))->GetTrack(fInputEvent,TruePhotonCandidate->GetTrackLabelPositive() );
1608 Double_t deltaPhi = magField * TVector2::Phi_mpi_pi( electronCandidate->Phi()-positronCandidate->Phi());
1609
1610 if(Photon->GetPdgCode() != 22){
1611 fHistoTrueDalitzPsiPairDeltaPhi[fiCut]->Fill(deltaPhi,TruePhotonCandidate->GetPsiPair());
1612 fCharPhotonMCInfo = 1;
1613 return; // Mother is no Photon
1614 }
1615
1616 if(((posDaughter->GetMCProcessCode())) != 5 || ((negDaughter->GetMCProcessCode())) != 5){
1617 fCharPhotonMCInfo = 1;
1618 return;// check if the daughters come from a conversion
1619 }
1620 // STILL A BUG IN ALIROOT >>8 HAS TPO BE REMOVED AFTER FIX
1621
1622
1623
1624 // True Photon
1625 if(fIsFromMBHeader){
1626 fHistoTrueConvGammaPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1627 if (fDoPhotonQA > 0) fHistoTrueConvGammaEta[fiCut]->Fill(TruePhotonCandidate->Eta());
1628 }
1629 fHistoTrueGammaPsiPairDeltaPhi[fiCut]->Fill(deltaPhi,TruePhotonCandidate->GetPsiPair());
1630 if(Photon->IsPrimary()){
1631 // Count just primary MC Gammas as true --> For Ratio esdtruegamma / mcconvgamma
1632 if(fIsFromMBHeader){
1633 fCharPhotonMCInfo = 6;
1634 fHistoTruePrimaryConvGammaPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1635 fHistoTruePrimaryConvGammaESDPtMCPt[fiCut]->Fill(TruePhotonCandidate->Pt(),Photon->Pt()); // Allways Filled
1636 }
1637 // (Not Filled for i6, Extra Signal Gamma (parambox) are secondary)
1638 } else {
1639 if(fIsFromMBHeader){
1640 fHistoTrueSecondaryConvGammaPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1641 fCharPhotonMCInfo = 2;
1642 if(((AliAODMCParticle*)AODMCTrackArray->At(Photon->GetMother()))->GetMother() > -1 &&
1643 ((AliAODMCParticle*)AODMCTrackArray->At(((AliAODMCParticle*)AODMCTrackArray->At(Photon->GetMother()))->GetMother()))->GetPdgCode() == 3122){
1644 fCharPhotonMCInfo = 5;
1645 fHistoTrueSecondaryConvGammaFromXFromLambdaPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1646 }
1647 if(((AliAODMCParticle*)AODMCTrackArray->At(Photon->GetMother()))->GetMother() > -1 &&
1648 ((AliAODMCParticle*)AODMCTrackArray->At(((AliAODMCParticle*)AODMCTrackArray->At(Photon->GetMother()))->GetMother()))->GetPdgCode() == 310){
1649 fCharPhotonMCInfo = 4;
1650 fHistoTrueSecondaryConvGammaFromXFromK0sPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1651 }
1652 if(((AliAODMCParticle*)AODMCTrackArray->At(Photon->GetMother()))->GetMother() > -1 &&
1653 ((AliAODMCParticle*)AODMCTrackArray->At(((AliAODMCParticle*)AODMCTrackArray->At(Photon->GetMother()))->GetMother()))->GetPdgCode() == 221){
1654 fCharPhotonMCInfo = 3;
1655 }
1656 }
1657 }
1658
1659}
1660//________________________________________________________________________
1661void AliAnalysisTaskGammaConvCalo::ProcessTruePhotonCandidates(AliAODConversionPhoton *TruePhotonCandidate)
1662{
1663
1664 Double_t magField = fInputEvent->GetMagneticField();
1665 if( magField < 0.0 ){
1666 magField = 1.0;
1667 }
1668 else {
1669 magField = -1.0;
1670 }
1671
1672 // Process True Photons
1673 TParticle *posDaughter = TruePhotonCandidate->GetPositiveMCDaughter(fMCStack);
1674 TParticle *negDaughter = TruePhotonCandidate->GetNegativeMCDaughter(fMCStack);
1675 fCharPhotonMCInfo = 0;
1676
1677 if(posDaughter == NULL || negDaughter == NULL) return; // One particle does not exist
1678 Int_t pdgCode[2] = {abs(posDaughter->GetPdgCode()),abs(negDaughter->GetPdgCode())};
1679 fCharPhotonMCInfo = 1;
1680 if(posDaughter->GetMother(0) != negDaughter->GetMother(0)){
1681 FillPhotonCombinatorialBackgroundHist(TruePhotonCandidate, pdgCode);
1682 return;
1683 }
1684 else if(posDaughter->GetMother(0) == -1){
1685 FillPhotonCombinatorialBackgroundHist(TruePhotonCandidate, pdgCode);
1686 return;
1687 }
1688
1689 if(pdgCode[0]!=11 || pdgCode[1]!=11) return; //One Particle is not a electron
1690
1691 if(posDaughter->GetPdgCode()==negDaughter->GetPdgCode()) return; // Same Charge
1692
1693 TParticle *Photon = TruePhotonCandidate->GetMCParticle(fMCStack);
1694 AliVTrack * electronCandidate = ((AliConversionCuts*)fCutArray->At(fiCut))->GetTrack(fInputEvent,TruePhotonCandidate->GetTrackLabelNegative() );
1695 AliVTrack * positronCandidate = ((AliConversionCuts*)fCutArray->At(fiCut))->GetTrack(fInputEvent,TruePhotonCandidate->GetTrackLabelPositive() );
1696 Double_t deltaPhi = magField * TVector2::Phi_mpi_pi( electronCandidate->Phi()-positronCandidate->Phi());
1697
1698 if(Photon->GetPdgCode() != 22){
1699 fHistoTrueDalitzPsiPairDeltaPhi[fiCut]->Fill(deltaPhi,TruePhotonCandidate->GetPsiPair());
1700 return; // Mother is no Photon
1701 }
1702
1703 if(posDaughter->GetUniqueID() != 5 || negDaughter->GetUniqueID() !=5) return;// check if the daughters come from a conversion
1704
1705
1706
1707 // True Photon
1708 if(fIsFromMBHeader){
1709 fHistoTrueConvGammaPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1710 if (fDoPhotonQA > 0) fHistoTrueConvGammaEta[fiCut]->Fill(TruePhotonCandidate->Eta());
1711 }
1712 fHistoTrueGammaPsiPairDeltaPhi[fiCut]->Fill(deltaPhi,TruePhotonCandidate->GetPsiPair());
1713 if(posDaughter->GetMother(0) <= fMCStack->GetNprimary()){
1714 // Count just primary MC Gammas as true --> For Ratio esdtruegamma / mcconvgamma
1715 if(fIsFromMBHeader){
1716 fCharPhotonMCInfo = 6;
1717 fHistoTruePrimaryConvGammaPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1718 fHistoTruePrimaryConvGammaESDPtMCPt[fiCut]->Fill(TruePhotonCandidate->Pt(),Photon->Pt()); // Allways Filled
1719
1720 }
1721 // (Not Filled for i6, Extra Signal Gamma (parambox) are secondary)
1722 }
1723 else{
1724 if(fIsFromMBHeader){
1725 fCharPhotonMCInfo = 2;
1726 fHistoTrueSecondaryConvGammaPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1727 if(fMCStack->Particle(Photon->GetMother(0))->GetMother(0) > -1 &&
1728 fMCStack->Particle(fMCStack->Particle(Photon->GetMother(0))->GetMother(0))->GetPdgCode() == 3122){
1729 fHistoTrueSecondaryConvGammaFromXFromLambdaPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1730 fCharPhotonMCInfo = 5;
1731 }
1732 if(fMCStack->Particle(Photon->GetMother(0))->GetMother(0) > -1 &&
1733 fMCStack->Particle(fMCStack->Particle(Photon->GetMother(0))->GetMother(0))->GetPdgCode() == 310){
1734 fHistoTrueSecondaryConvGammaFromXFromK0sPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1735 fCharPhotonMCInfo = 4;
1736 }
1737 if(fMCStack->Particle(Photon->GetMother(0))->GetMother(0) > -1 &&
1738 fMCStack->Particle(fMCStack->Particle(Photon->GetMother(0))->GetMother(0))->GetPdgCode() == 221){
1739 fCharPhotonMCInfo = 3;
1740 }
1741 }
1742 }
1743
1744 // pi0 photon
1745 //Bool_t bpi0 = 0;
1746 Int_t imother = Photon->GetMother(0);
1747 AliMCParticle *McMother = static_cast<AliMCParticle*>(fMCEvent->GetTrack(imother));
1748 if(McMother->PdgCode() == 111){
1749 fHistoTrueConvPi0GammaPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1750 }
1751
1752}
1753//________________________________________________________________________
1754void AliAnalysisTaskGammaConvCalo::ProcessAODMCParticles()
1755{
1756
1757 TClonesArray *AODMCTrackArray = dynamic_cast<TClonesArray*>(fInputEvent->FindListObject(AliAODMCParticle::StdBranchName()));
1758
1759 // Loop over all primary MC particle
1760 for(Int_t i = 0; i < AODMCTrackArray->GetEntriesFast(); i++) {
1761
1762 AliAODMCParticle* particle = static_cast<AliAODMCParticle*>(AODMCTrackArray->At(i));
1763 if (!particle) continue;
1764 if (!particle->IsPrimary()) continue;
1765
1766 Int_t isMCFromMBHeader = -1;
1767 if(((AliConversionCuts*)fCutArray->At(fiCut))->GetSignalRejection() != 0){
1768 isMCFromMBHeader
1769 = ((AliConversionCuts*)fCutArray->At(fiCut))->IsParticleFromBGEvent(i, fMCStack, fInputEvent);
1770 if(isMCFromMBHeader == 0 && ((AliConversionCuts*)fCutArray->At(fiCut))->GetSignalRejection() != 3) continue;
1771 }
1772
1773 if(!((AliConversionCuts*)fCutArray->At(fiCut))->InPlaneOutOfPlaneCut(particle->Phi(),fEventPlaneAngle,kFALSE)) continue;
1774 if(((AliConversionCuts*)fCutArray->At(fiCut))->PhotonIsSelectedAODMC(particle,AODMCTrackArray,kFALSE)){
1775 fHistoMCAllGammaPt[fiCut]->Fill(particle->Pt()); // All MC Gamma
1776 if(particle->GetMother() >-1){ // Meson Decay Gamma
1777 switch((static_cast<AliAODMCParticle*>(AODMCTrackArray->At(particle->GetMother())))->GetPdgCode()){
1778 case 111: // Pi0
1779 fHistoMCDecayGammaPi0Pt[fiCut]->Fill(particle->Pt());
1780 break;
1781 case 113: // Rho0
1782 fHistoMCDecayGammaRhoPt[fiCut]->Fill(particle->Pt());
1783 break;
1784 case 221: // Eta
1785 fHistoMCDecayGammaEtaPt[fiCut]->Fill(particle->Pt());
1786 break;
1787 case 223: // Omega
1788 fHistoMCDecayGammaOmegaPt[fiCut]->Fill(particle->Pt());
1789 break;
1790 case 331: // Eta'
1791 fHistoMCDecayGammaEtapPt[fiCut]->Fill(particle->Pt());
1792 break;
1793 case 333: // Phi
1794 fHistoMCDecayGammaPhiPt[fiCut]->Fill(particle->Pt());
1795 break;
1796 case 3212: // Sigma
1797 fHistoMCDecayGammaSigmaPt[fiCut]->Fill(particle->Pt());
1798 break;
1799 }
1800 }
1801 }
1802 if(((AliConversionCuts*)fCutArray->At(fiCut))->PhotonIsSelectedAODMC(particle,AODMCTrackArray,kTRUE)){
1803 Double_t rConv = 0;
1804 for(Int_t daughterIndex=particle->GetDaughter(0);daughterIndex<=particle->GetDaughter(1);daughterIndex++){
1805 AliAODMCParticle *tmpDaughter = static_cast<AliAODMCParticle*>(AODMCTrackArray->At(daughterIndex));
1806 if(!tmpDaughter) continue;
1807 if(abs(tmpDaughter->GetPdgCode()) == 11){
1808 rConv = sqrt( (tmpDaughter->Xv()*tmpDaughter->Xv()) + (tmpDaughter->Yv()*tmpDaughter->Yv()) );
1809 }
1810 }
1811 fHistoMCConvGammaPt[fiCut]->Fill(particle->Pt());
1812 if (fDoPhotonQA > 0){
1813 fHistoMCConvGammaR[fiCut]->Fill(rConv);
1814 fHistoMCConvGammaEta[fiCut]->Fill(particle->Eta());
1815 }
1816 }
1817 // Converted MC Gamma
1818 if(fDoMesonAnalysis){
1819 if(particle->GetPdgCode() == 310 && fDoMesonQA > 0){
1820 Double_t mesonY = 10.;
1821 if(particle->E() - particle->Pz() == 0 || particle->E() + particle->Pz() == 0){
1822 mesonY=10.-((AliConversionCuts*)fCutArray->At(fiCut))->GetEtaShift();
1823 } else {
1824 mesonY = 0.5*(TMath::Log((particle->E()+particle->Pz()) / (particle->E()-particle->Pz())))-((AliConversionCuts*)fCutArray->At(fiCut))->GetEtaShift();
1825 }
1826 Float_t weightedK0s= 1;
1827 if(((AliConversionCuts*)fCutArray->At(fiCut))->IsParticleFromBGEvent(i, fMCStack, fInputEvent)){
1828 if (particle->Pt()>0.005){
1829 weightedK0s= ((AliConversionCuts*)fCutArray->At(fiCut))->GetWeightForMeson(fV0Reader->GetPeriodName(),i, 0x0, fInputEvent);
1830 //cout << "MC input \t"<<i << "\t" << particle->Pt()<<"\t"<<weighted << endl;
1831 }
1832 }
1833 fHistoMCK0sPt[fiCut]->Fill(particle->Pt(),weightedK0s);
1834 fHistoMCK0sWOWeightPt[fiCut]->Fill(particle->Pt());
1835 fHistoMCK0sPtY[fiCut]->Fill(particle->Pt(),mesonY,weightedK0s);
1836 }
1837 if(((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))
1838 ->MesonIsSelectedAODMC(particle,AODMCTrackArray,((AliConversionCuts*)fCutArray->At(fiCut))->GetEtaShift())){
1839 AliAODMCParticle* daughter0 = static_cast<AliAODMCParticle*>(AODMCTrackArray->At(particle->GetDaughter(0)));
1840 AliAODMCParticle* daughter1 = static_cast<AliAODMCParticle*>(AODMCTrackArray->At(particle->GetDaughter(1)));
1841 Float_t weighted= 1;
1842 if(((AliConversionCuts*)fCutArray->At(fiCut))->IsParticleFromBGEvent(i, fMCStack, fInputEvent)){
1843 if (particle->Pt()>0.005){
1844 weighted= ((AliConversionCuts*)fCutArray->At(fiCut))->GetWeightForMeson(fV0Reader->GetPeriodName(),i, 0x0, fInputEvent);
1845 // if(particle->GetPdgCode() == 221){
1846 // cout << "MC input \t"<<i << "\t" << particle->Pt()<<"\t"<<weighted << endl;
1847 // }
1848 }
1849 }
1850 Double_t mesonY = 10.;
1851 if(particle->E() - particle->Pz() == 0 || particle->E() + particle->Pz() == 0){
1852 mesonY=10.-((AliConversionCuts*)fCutArray->At(fiCut))->GetEtaShift();
1853 } else{
1854 mesonY = 0.5*(TMath::Log((particle->E()+particle->Pz()) / (particle->E()-particle->Pz())))-((AliConversionCuts*)fCutArray->At(fiCut))->GetEtaShift();
1855 }
1856
1857 if(particle->GetPdgCode() == 111){
1858 fHistoMCPi0Pt[fiCut]->Fill(particle->Pt(),weighted); // All MC Pi0
1859 fHistoMCPi0WOWeightPt[fiCut]->Fill(particle->Pt());
1860 if (fDoMesonQA > 0) fHistoMCPi0PtY[fiCut]->Fill(particle->Pt(),mesonY,weighted); // All MC Pi0
1861 } else if(particle->GetPdgCode() == 221){
1862 fHistoMCEtaPt[fiCut]->Fill(particle->Pt(),weighted); // All MC Eta
1863 fHistoMCEtaWOWeightPt[fiCut]->Fill(particle->Pt());
1864 if (fDoMesonQA > 0) fHistoMCEtaPtY[fiCut]->Fill(particle->Pt(),mesonY,weighted); // All MC Pi0
1865 }
1866
1867 // Check the acceptance for both gammas
1868 if(((AliConversionCuts*)fCutArray->At(fiCut))->PhotonIsSelectedAODMC(daughter0,AODMCTrackArray,kFALSE) &&
1869 ((AliConversionCuts*)fCutArray->At(fiCut))->PhotonIsSelectedAODMC(daughter1,AODMCTrackArray,kFALSE) &&
1870 ((AliConversionCuts*)fCutArray->At(fiCut))->InPlaneOutOfPlaneCut(daughter0->Phi(),fEventPlaneAngle,kFALSE) &&
1871 ((AliConversionCuts*)fCutArray->At(fiCut))->InPlaneOutOfPlaneCut(daughter1->Phi(),fEventPlaneAngle,kFALSE)){
1872
1873 if(particle->GetPdgCode() == 111){
1874 fHistoMCPi0InAccPt[fiCut]->Fill(particle->Pt(),weighted); // MC Pi0 with gamma in acc
1875 } else if(particle->GetPdgCode() == 221){
1876 fHistoMCEtaInAccPt[fiCut]->Fill(particle->Pt(),weighted); // MC Eta with gamma in acc
1877 }
1878 }
1879 }
1880 }
1881 }
1882
1883}
1884//________________________________________________________________________
1885void AliAnalysisTaskGammaConvCalo::ProcessMCParticles()
1886{
1887 // Loop over all primary MC particle
1888 for(Int_t i = 0; i < fMCStack->GetNprimary(); i++) {
1889 TParticle* particle = (TParticle *)fMCStack->Particle(i);
1890 if (!particle) continue;
1891
1892 Int_t isMCFromMBHeader = -1;
1893 if(((AliConversionCuts*)fCutArray->At(fiCut))->GetSignalRejection() != 0){
1894 isMCFromMBHeader
1895 = ((AliConversionCuts*)fCutArray->At(fiCut))->IsParticleFromBGEvent(i, fMCStack, fInputEvent);
1896 if(isMCFromMBHeader == 0 && ((AliConversionCuts*)fCutArray->At(fiCut))->GetSignalRejection() != 3) continue;
1897 }
1898
1899 if(!((AliConversionCuts*)fCutArray->At(fiCut))->InPlaneOutOfPlaneCut(particle->Phi(),fEventPlaneAngle,kFALSE)) continue;
1900 if(((AliConversionCuts*)fCutArray->At(fiCut))->PhotonIsSelectedMC(particle,fMCStack,kFALSE)){
1901 fHistoMCAllGammaPt[fiCut]->Fill(particle->Pt()); // All MC Gamma
1902 if(particle->GetMother(0) >-1){ // Meson Decay Gamma
1903 switch(fMCStack->Particle(particle->GetMother(0))->GetPdgCode()){
1904 case 111: // Pi0
1905 fHistoMCDecayGammaPi0Pt[fiCut]->Fill(particle->Pt());
1906 break;
1907 case 113: // Rho0
1908 fHistoMCDecayGammaRhoPt[fiCut]->Fill(particle->Pt());
1909 break;
1910 case 221: // Eta
1911 fHistoMCDecayGammaEtaPt[fiCut]->Fill(particle->Pt());
1912 break;
1913 case 223: // Omega
1914 fHistoMCDecayGammaOmegaPt[fiCut]->Fill(particle->Pt());
1915 break;
1916 case 331: // Eta'
1917 fHistoMCDecayGammaEtapPt[fiCut]->Fill(particle->Pt());
1918 break;
1919 case 333: // Phi
1920 fHistoMCDecayGammaPhiPt[fiCut]->Fill(particle->Pt());
1921 break;
1922 case 3212: // Sigma
1923 fHistoMCDecayGammaSigmaPt[fiCut]->Fill(particle->Pt());
1924 break;
1925 }
1926 }
1927 }
1928 if(((AliConversionCuts*)fCutArray->At(fiCut))->PhotonIsSelectedMC(particle,fMCStack,kTRUE)){
1929 fHistoMCConvGammaPt[fiCut]->Fill(particle->Pt());
1930 if (fDoPhotonQA > 0){
1931 fHistoMCConvGammaR[fiCut]->Fill(((TParticle*)fMCStack->Particle(particle->GetFirstDaughter()))->R());
1932 fHistoMCConvGammaEta[fiCut]->Fill(particle->Eta());
1933 }
1934 } // Converted MC Gamma
1935 if(fDoMesonAnalysis){
1936 if(particle->GetPdgCode() == 310 && fDoMesonQA > 0){
1937 Double_t mesonY = 10.;
1938 if(particle->Energy() - particle->Pz() == 0 || particle->Energy() + particle->Pz() == 0){
1939 mesonY=10.-((AliConversionCuts*)fCutArray->At(fiCut))->GetEtaShift();
1940 } else{
1941 mesonY = 0.5*(TMath::Log((particle->Energy()+particle->Pz()) / (particle->Energy()-particle->Pz())))-((AliConversionCuts*)fCutArray->At(fiCut))->GetEtaShift();
1942 }
1943 Float_t weightedK0s= 1;
1944 if(((AliConversionCuts*)fCutArray->At(fiCut))->IsParticleFromBGEvent(i, fMCStack, fInputEvent)){
1945 if (particle->Pt()>0.005){
1946 weightedK0s= ((AliConversionCuts*)fCutArray->At(fiCut))->GetWeightForMeson(fV0Reader->GetPeriodName(),i, fMCStack, fInputEvent);
1947 //cout << "MC input \t"<<i << "\t" << particle->Pt()<<"\t"<<weighted << endl;
1948 }
1949 }
1950 if (fMCStack->IsPhysicalPrimary(i)){
1951 fHistoMCK0sPt[fiCut]->Fill(particle->Pt(),weightedK0s);
1952 fHistoMCK0sWOWeightPt[fiCut]->Fill(particle->Pt());
1953 fHistoMCK0sPtY[fiCut]->Fill(particle->Pt(),mesonY,weightedK0s);
1954 }
1955 }
1956 if(((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))
1957 ->MesonIsSelectedMC(particle,fMCStack,((AliConversionCuts*)fCutArray->At(fiCut))->GetEtaShift())){
1958 TParticle* daughter0 = (TParticle*)fMCStack->Particle(particle->GetFirstDaughter());
1959 TParticle* daughter1 = (TParticle*)fMCStack->Particle(particle->GetLastDaughter());
1960
1961 Float_t weighted= 1;
1962 if(((AliConversionCuts*)fCutArray->At(fiCut))->IsParticleFromBGEvent(i, fMCStack, fInputEvent)){
1963 if (particle->Pt()>0.005){
1964 weighted= ((AliConversionCuts*)fCutArray->At(fiCut))->GetWeightForMeson(fV0Reader->GetPeriodName(),i, fMCStack, fInputEvent);
1965 // if(particle->GetPdgCode() == 221){
1966 // cout << "MC input \t"<<i << "\t" << particle->Pt()<<"\t"<<weighted << endl;
1967 // }
1968 }
1969 }
1970 Double_t mesonY = 10.;
1971 if(particle->Energy() - particle->Pz() == 0 || particle->Energy() + particle->Pz() == 0){
1972 mesonY=10.-((AliConversionCuts*)fCutArray->At(fiCut))->GetEtaShift();
1973 } else{
1974 mesonY = 0.5*(TMath::Log((particle->Energy()+particle->Pz()) / (particle->Energy()-particle->Pz())))-((AliConversionCuts*)fCutArray->At(fiCut))->GetEtaShift();
1975 }
1976
1977 if(particle->GetPdgCode() == 111){
1978 fHistoMCPi0Pt[fiCut]->Fill(particle->Pt(),weighted); // All MC Pi0
1979 fHistoMCPi0WOWeightPt[fiCut]->Fill(particle->Pt());
1980 if (fDoMesonQA > 0) fHistoMCPi0PtY[fiCut]->Fill(particle->Pt(),mesonY,weighted); // All MC Pi0
1981 } else if(particle->GetPdgCode() == 221){
1982 fHistoMCEtaPt[fiCut]->Fill(particle->Pt(),weighted); // All MC Eta
1983 fHistoMCEtaWOWeightPt[fiCut]->Fill(particle->Pt());
1984 if (fDoMesonQA > 0) fHistoMCEtaPtY[fiCut]->Fill(particle->Pt(),mesonY,weighted); // All MC Pi0
1985 }
1986
1987 // Check the acceptance for both gammas
1988 if(((AliConversionCuts*)fCutArray->At(fiCut))->PhotonIsSelectedMC(daughter0,fMCStack,kFALSE) &&
1989 ((AliConversionCuts*)fCutArray->At(fiCut))->PhotonIsSelectedMC(daughter1,fMCStack,kFALSE) &&
1990 ((AliConversionCuts*)fCutArray->At(fiCut))->InPlaneOutOfPlaneCut(daughter0->Phi(),fEventPlaneAngle,kFALSE) &&
1991 ((AliConversionCuts*)fCutArray->At(fiCut))->InPlaneOutOfPlaneCut(daughter1->Phi(),fEventPlaneAngle,kFALSE)){
1992
1993 if(particle->GetPdgCode() == 111){
1994 fHistoMCPi0InAccPt[fiCut]->Fill(particle->Pt(),weighted); // MC Pi0 with gamma in acc
1995 } else if(particle->GetPdgCode() == 221){
1996 fHistoMCEtaInAccPt[fiCut]->Fill(particle->Pt(),weighted); // MC Eta with gamma in acc
1997 }
1998 }
1999 }
2000 }
2001 }
2002
2003 if (fDoMesonQA){
2004 for(Int_t i = fMCStack->GetNprimary(); i < fMCStack->GetNtrack(); i++) {
2005 TParticle* particle = (TParticle *)fMCStack->Particle(i);
2006 if (!particle) continue;
2007
2008 Int_t isMCFromMBHeader = -1;
2009 if(((AliConversionCuts*)fCutArray->At(fiCut))->GetSignalRejection() != 0){
2010 isMCFromMBHeader
2011 = ((AliConversionCuts*)fCutArray->At(fiCut))->IsParticleFromBGEvent(i, fMCStack, fInputEvent);
2012 if(isMCFromMBHeader == 0 && ((AliConversionCuts*)fCutArray->At(fiCut))->GetSignalRejection() != 3) continue;
2013 }
2014
2015 if(fDoMesonAnalysis){
2016 if(((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))
2017 ->MesonIsSelectedMC(particle,fMCStack,((AliConversionCuts*)fCutArray->At(fiCut))->GetEtaShift())){
2018 Float_t weighted= 1;
2019 if(((AliConversionCuts*)fCutArray->At(fiCut))->IsParticleFromBGEvent(i, fMCStack, fInputEvent)){
2020 if (particle->Pt()>0.005){
2021 weighted= ((AliConversionCuts*)fCutArray->At(fiCut))->GetWeightForMeson(fV0Reader->GetPeriodName(),i, fMCStack, fInputEvent);
2022 // if(particle->GetPdgCode() == 221){
2023 // cout << "MC input \t"<<i << "\t" << particle->Pt()<<"\t"<<weighted << endl;
2024 // }
2025 }
2026 }
2027
2028 if(particle->GetPdgCode() == 111){
2029 Int_t pdgCode = ((TParticle*)fMCStack->Particle( particle->GetFirstMother() ))->GetPdgCode();
2030 Int_t source = GetSourceClassification(111,pdgCode);
2031 fHistoMCSecPi0PtvsSource[fiCut]->Fill(particle->Pt(),source,weighted); // All MC Pi0
2032 fHistoMCSecPi0Source[fiCut]->Fill(pdgCode);
2033 } else if(particle->GetPdgCode() == 221){
2034 Int_t pdgCode = ((TParticle*)fMCStack->Particle( particle->GetFirstMother() ))->GetPdgCode();
2035 fHistoMCSecEtaPt[fiCut]->Fill(particle->Pt(),weighted); // All MC Pi0
2036 fHistoMCSecEtaSource[fiCut]->Fill(pdgCode);
2037 }
2038 }
2039 }
2040 }
2041 }
2042}
2043
2044//________________________________________________________________________
2045void AliAnalysisTaskGammaConvCalo::PhotonTagging(){
2046
2047 // Conversion Gammas
2048 if(fGammaCandidates->GetEntries()>0){
2049
2050 for(Int_t firstGammaIndex=0;firstGammaIndex<fGammaCandidates->GetEntries()-1;firstGammaIndex++){
2051
2052 // get conversion photon
2053 AliAODConversionPhoton *gamma0=dynamic_cast<AliAODConversionPhoton*>(fGammaCandidates->At(firstGammaIndex));
2054 if (gamma0==NULL) continue;
2055
2056 TLorentzVector photonVector;
2057 photonVector.SetPxPyPzE(gamma0->GetPx(),gamma0->GetPy(),gamma0->GetPz(),gamma0->GetPhotonP());
2058
2059 Bool_t btagpi0 = 0;
2060 Bool_t btageta = 0;
2061
2062 // loop over clusters
2063 for(Int_t secondGammaIndex = 0; secondGammaIndex<fClusterCandidates->GetEntries(); ++secondGammaIndex) {
2064
2065 AliAODConversionPhoton *gamma1=dynamic_cast<AliAODConversionPhoton*>(fClusterCandidates->At(secondGammaIndex));
2066 if (gamma1==NULL) continue;
2067
2068 TLorentzVector clusterVector;
2069 clusterVector.SetPxPyPzE(gamma1->GetPx(),gamma1->GetPy(),gamma1->GetPz(),gamma1->GetPhotonP());
2070
2071 // do the tagging
2072 TLorentzVector pairVector = photonVector+clusterVector;
2073
2074 // see if pi0?
2075 if((pairVector.M() > 0.11 && pairVector.M() < 0.15)){
2076 btagpi0 = 1;
2077 }
2078 // or eta
2079 if((pairVector.M() > 0.50 && pairVector.M() < 0.6)){
2080 btageta = 1;
2081 }
2082 }// end loop over clusters
2083
2084 if(btagpi0 && btageta)
2085 fHistoConvGammaTagged[fiCut]->Fill(photonVector.Pt());
2086 else if(btagpi0 && !btageta)
2087 fHistoConvGammaPi0Tagged[fiCut]->Fill(photonVector.Pt());
2088 else if(btageta && !btagpi0)
2089 fHistoConvGammaEtaTagged[fiCut]->Fill(photonVector.Pt());
2090 else
2091 fHistoConvGammaUntagged[fiCut]->Fill(photonVector.Pt());
2092
2093 }// end loop over gammas
2094 }// end if
2095 return;
2096}
2097
2098//________________________________________________________________________
2099void AliAnalysisTaskGammaConvCalo::CalculatePi0Candidates(){
2100
2101 // Conversion Gammas
2102 if(fGammaCandidates->GetEntries()>0){
2103
2104 // vertex
2105 Double_t vertex[3] = {0};
2106 InputEvent()->GetPrimaryVertex()->GetXYZ(vertex);
2107
2108 for(Int_t firstGammaIndex=0;firstGammaIndex<fGammaCandidates->GetEntries();firstGammaIndex++){
2109 AliAODConversionPhoton *gamma0=dynamic_cast<AliAODConversionPhoton*>(fGammaCandidates->At(firstGammaIndex));
2110 if (gamma0==NULL) continue;
2111
2112 for(Int_t secondGammaIndex=0;secondGammaIndex<fClusterCandidates->GetEntries();secondGammaIndex++){
2113
2114 AliAODConversionPhoton *gamma1=dynamic_cast<AliAODConversionPhoton*>(fClusterCandidates->At(secondGammaIndex));
2115 if (gamma1==NULL) continue;
2116
2117 AliAODConversionMother *pi0cand = new AliAODConversionMother(gamma0,gamma1);
2118 pi0cand->SetLabels(firstGammaIndex,secondGammaIndex);
2119
2120 if((((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->MesonIsSelected(pi0cand,kTRUE,((AliConversionCuts*)fCutArray->At(fiCut))->GetEtaShift()))){
2121 fHistoMotherInvMassPt[fiCut]->Fill(pi0cand->M(),pi0cand->Pt());
2122
2123 // fill new histograms
2124 fHistoPhotonPairAll[fiCut]->Fill(pi0cand->M(),pi0cand->Pt());
2125 fHistoPhotonPairAllGam[fiCut]->Fill(pi0cand->M(),gamma0->Pt());
2126
2127 if(pi0cand->GetAlpha()<0.1)
2128 fHistoMotherInvMassEalpha[fiCut]->Fill(pi0cand->M(),pi0cand->E());
2129
2130 if (fDoMesonQA > 0){
2131 if ( pi0cand->M() > 0.05 && pi0cand->M() < 0.17){
2132 fHistoMotherPi0PtY[fiCut]->Fill(pi0cand->Pt(),pi0cand->Rapidity()-((AliConversionCuts*)fCutArray->At(fiCut))->GetEtaShift());
2133 fHistoMotherPi0PtAlpha[fiCut]->Fill(pi0cand->Pt(),pi0cand->GetAlpha());
2134 fHistoMotherPi0PtOpenAngle[fiCut]->Fill(pi0cand->Pt(),pi0cand->GetOpeningAngle());
2135 }
2136 if ( pi0cand->M() > 0.45 && pi0cand->M() < 0.65){
2137 fHistoMotherEtaPtY[fiCut]->Fill(pi0cand->Pt(),pi0cand->Rapidity()-((AliConversionCuts*)fCutArray->At(fiCut))->GetEtaShift());
2138 fHistoMotherEtaPtAlpha[fiCut]->Fill(pi0cand->Pt(),pi0cand->GetAlpha());
2139 fHistoMotherEtaPtOpenAngle[fiCut]->Fill(pi0cand->Pt(),pi0cand->GetOpeningAngle());
2140 }
2141 }
2142 if(((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->DoBGCalculation()){
2143 Int_t zbin = 0;
2144 Int_t mbin = 0;
2145
2146 if(((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->BackgroundHandlerType() == 0){
2147 zbin = fBGHandler[fiCut]->GetZBinIndex(fInputEvent->GetPrimaryVertex()->GetZ());
2148 if(((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->UseTrackMultiplicity()){
2149 mbin = fBGHandler[fiCut]->GetMultiplicityBinIndex(fV0Reader->GetNumberOfPrimaryTracks());
2150 } else {
2151 mbin = fBGHandler[fiCut]->GetMultiplicityBinIndex(fGammaCandidates->GetEntries());
2152 }
2153 } else{
2154 zbin = fBGHandlerRP[fiCut]->GetZBinIndex(fInputEvent->GetPrimaryVertex()->GetZ());
2155 if(((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->UseTrackMultiplicity()){
2156 mbin = fBGHandlerRP[fiCut]->GetMultiplicityBinIndex(fV0Reader->GetNumberOfPrimaryTracks());
2157 } else {
2158 mbin = fBGHandlerRP[fiCut]->GetMultiplicityBinIndex(fGammaCandidates->GetEntries());
2159 }
2160 }
2161 Double_t sparesFill[4] = {pi0cand->M(),pi0cand->Pt(),(Double_t)zbin,(Double_t)mbin};
2162 fSparseMotherInvMassPtZM[fiCut]->Fill(sparesFill,1);
2163 }
2164
2165
2166 if(fIsMC){
2167 if(fInputEvent->IsA()==AliESDEvent::Class())
2168 ProcessTrueMesonCandidates(pi0cand,gamma0,gamma1);
2169 if(fInputEvent->IsA()==AliAODEvent::Class())
2170 ProcessTrueMesonCandidatesAOD(pi0cand,gamma0,gamma1);
2171 }
2172 if (fDoMesonQA == 2){
2173 fInvMass = pi0cand->M();
2174 fPt = pi0cand->Pt();
2175 if (abs(gamma0->GetDCAzToPrimVtx()) < abs(gamma1->GetDCAzToPrimVtx())){
2176 fDCAzGammaMin = gamma0->GetDCAzToPrimVtx();
2177 fDCAzGammaMax = gamma1->GetDCAzToPrimVtx();
2178 } else {
2179 fDCAzGammaMin = gamma1->GetDCAzToPrimVtx();
2180 fDCAzGammaMax = gamma0->GetDCAzToPrimVtx();
2181 }
2182 fCharFlag = pi0cand->GetMesonQuality();
2183 // cout << "gamma 0: " << gamma0->GetV0Index()<< "\t" << gamma0->GetPx() << "\t" << gamma0->GetPy() << "\t" << gamma0->GetPz() << "\t" << endl;
2184 // cout << "gamma 1: " << gamma1->GetV0Index()<< "\t"<< gamma1->GetPx() << "\t" << gamma1->GetPy() << "\t" << gamma1->GetPz() << "\t" << endl;
2185 // cout << "pi0: "<<fInvMass << "\t" << fPt <<"\t" << fDCAzGammaMin << "\t" << fDCAzGammaMax << "\t" << (Int_t)fCharFlag << "\t" << (Int_t)fCharMesonMCInfo <<endl;
2186 if (fIsHeavyIon == 1 && fPt > 0.399 && fPt < 20. ) {
2187 if (fInvMass > 0.08 && fInvMass < 0.2) fTreeMesonsInvMassPtDcazMinDcazMaxFlag[fiCut]->Fill();
2188 if ((fInvMass > 0.45 && fInvMass < 0.6) && (fPt > 0.999 && fPt < 20.) )fTreeMesonsInvMassPtDcazMinDcazMaxFlag[fiCut]->Fill();
2189 } else if (fPt > 0.299 && fPt < 20. ) {
2190 if ( (fInvMass > 0.08 && fInvMass < 0.2) || (fInvMass > 0.45 && fInvMass < 0.6)) fTreeMesonsInvMassPtDcazMinDcazMaxFlag[fiCut]->Fill();
2191 }
2192 }
2193 }
2194 delete pi0cand;
2195 pi0cand=0x0;
2196 }
2197 }
2198 }
2199}
2200//______________________________________________________________________
2201void AliAnalysisTaskGammaConvCalo::ProcessTrueMesonCandidates(AliAODConversionMother *Pi0Candidate, AliAODConversionPhoton *TrueGammaCandidate0, AliAODConversionPhoton *TrueGammaCandidate1)
2202{
2203 // Process True Mesons
2204 AliStack *MCStack = fMCEvent->Stack();
2205 fCharMesonMCInfo = 0;
2206 if(TrueGammaCandidate0->GetV0Index()<fInputEvent->GetNumberOfV0s()){
2207 Bool_t isTruePi0 = kFALSE;
2208 Bool_t isTrueEta = kFALSE;
2209 Bool_t isTruePi0Dalitz = kFALSE;
2210 Bool_t isTrueEtaDalitz = kFALSE;
2211 Bool_t gamma0DalitzCand = kFALSE;
2212 Bool_t gamma1DalitzCand = kFALSE;
2213 Int_t gamma0MCLabel = TrueGammaCandidate0->GetMCParticleLabel(MCStack);
2214 Int_t gamma0MotherLabel = -1;
2215 if(gamma0MCLabel != -1){ // Gamma is Combinatorial; MC Particles don't belong to the same Mother
2216 // Daughters Gamma 0
2217 TParticle * negativeMC = (TParticle*)TrueGammaCandidate0->GetNegativeMCDaughter(MCStack);
2218 TParticle * positiveMC = (TParticle*)TrueGammaCandidate0->GetPositiveMCDaughter(MCStack);
2219 TParticle * gammaMC0 = (TParticle*)MCStack->Particle(gamma0MCLabel);
2220 if(abs(negativeMC->GetPdgCode())==11 && abs(positiveMC->GetPdgCode())==11){ // Electrons ...
2221 if(negativeMC->GetUniqueID() == 5 && positiveMC->GetUniqueID() ==5){ // ... From Conversion ...
2222 if(gammaMC0->GetPdgCode() == 22){ // ... with Gamma Mother
2223 gamma0MotherLabel=gammaMC0->GetFirstMother();
2224 }
2225 }
2226 if(gammaMC0->GetPdgCode() ==111){ // Dalitz candidate
2227 gamma0DalitzCand = kTRUE;
2228 gamma0MotherLabel=-111;
2229 }
2230 if(gammaMC0->GetPdgCode() ==221){ // Dalitz candidate
2231 gamma0DalitzCand = kTRUE;
2232 gamma0MotherLabel=-221;
2233 }
2234 }
2235 }
2236 if(TrueGammaCandidate1->GetV0Index()<fInputEvent->GetNumberOfV0s()){
2237 Int_t gamma1MCLabel = TrueGammaCandidate1->GetMCParticleLabel(MCStack);
2238 Int_t gamma1MotherLabel = -1;
2239 if(gamma1MCLabel != -1){ // Gamma is Combinatorial; MC Particles don't belong to the same Mother
2240 // Daughters Gamma 1
2241 TParticle * negativeMC = (TParticle*)TrueGammaCandidate1->GetNegativeMCDaughter(MCStack);
2242 TParticle * positiveMC = (TParticle*)TrueGammaCandidate1->GetPositiveMCDaughter(MCStack);
2243 TParticle * gammaMC1 = (TParticle*)MCStack->Particle(gamma1MCLabel);
2244 if(abs(negativeMC->GetPdgCode())==11 && abs(positiveMC->GetPdgCode())==11){ // Electrons ...
2245 if(negativeMC->GetUniqueID() == 5 && positiveMC->GetUniqueID() ==5){ // ... From Conversion ...
2246 if(gammaMC1->GetPdgCode() == 22){ // ... with Gamma Mother
2247 gamma1MotherLabel=gammaMC1->GetFirstMother();
2248 }
2249 }
2250 if(gammaMC1->GetPdgCode() ==111 ){ // Dalitz candidate
2251 gamma1DalitzCand = kTRUE;
2252 gamma1MotherLabel=-111;
2253 }
2254 if(gammaMC1->GetPdgCode() ==221){ // Dalitz candidate
2255 gamma1DalitzCand = kTRUE;
2256 gamma1MotherLabel=-221;
2257 }
2258 }
2259 }
2260 if(gamma0MotherLabel>=0 && gamma0MotherLabel==gamma1MotherLabel){
2261 if(((TParticle*)MCStack->Particle(gamma1MotherLabel))->GetPdgCode() == 111){
2262 isTruePi0=kTRUE;
2263 }
2264 if(((TParticle*)MCStack->Particle(gamma1MotherLabel))->GetPdgCode() == 221){
2265 isTrueEta=kTRUE;
2266 }
2267 }
2268
2269 //Identify Dalitz candidate
2270 if (gamma1DalitzCand || gamma0DalitzCand){
2271 if (gamma0DalitzCand && gamma0MCLabel >=0 && gamma0MCLabel==gamma1MotherLabel){
2272 if (gamma0MotherLabel == -111) isTruePi0Dalitz = kTRUE;
2273 if (gamma0MotherLabel == -221) isTrueEtaDalitz = kTRUE;
2274 }
2275 if (gamma1DalitzCand && gamma1MCLabel >=0 && gamma1MCLabel==gamma0MotherLabel){
2276 if (gamma1MotherLabel == -111) isTruePi0Dalitz = kTRUE;
2277 if (gamma1MotherLabel == -221) isTrueEtaDalitz = kTRUE;
2278 }
2279 }
2280
2281
2282 if(isTruePi0 || isTrueEta){// True Pion or Eta
2283 fHistoTrueMotherInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt());
2284 if (fDoMesonQA > 0){
2285 if (isTruePi0){
2286 if ( Pi0Candidate->M() > 0.05 && Pi0Candidate->M() < 0.17){
2287 fHistoTruePi0PtY[fiCut]->Fill(Pi0Candidate->Pt(),Pi0Candidate->Rapidity()-((AliConversionCuts*)fCutArray->At(fiCut))->GetEtaShift());
2288 fHistoTruePi0PtAlpha[fiCut]->Fill(Pi0Candidate->Pt(),Pi0Candidate->GetAlpha());
2289 fHistoTruePi0PtOpenAngle[fiCut]->Fill(Pi0Candidate->Pt(),Pi0Candidate->GetOpeningAngle());
2290 }
2291 } else if (isTrueEta){
2292 if ( Pi0Candidate->M() > 0.45 && Pi0Candidate->M() < 0.65){
2293 fHistoTrueEtaPtY[fiCut]->Fill(Pi0Candidate->Pt(),Pi0Candidate->Rapidity()-((AliConversionCuts*)fCutArray->At(fiCut))->GetEtaShift());
2294 fHistoTrueEtaPtAlpha[fiCut]->Fill(Pi0Candidate->Pt(),Pi0Candidate->GetAlpha());
2295 fHistoTrueEtaPtOpenAngle[fiCut]->Fill(Pi0Candidate->Pt(),Pi0Candidate->GetOpeningAngle());
2296 }
2297 }
2298 }
2299 if(gamma0MotherLabel >= MCStack->GetNprimary()){ // Secondary Meson
2300 Int_t secMotherLabel = ((TParticle*)MCStack->Particle(gamma1MotherLabel))->GetMother(0);
2301 Float_t weightedSec= 1;
2302 if(((AliConversionCuts*)fCutArray->At(fiCut))->IsParticleFromBGEvent(secMotherLabel, fMCStack, fInputEvent) && MCStack->Particle(secMotherLabel)->GetPdgCode()==310){
2303 weightedSec= ((AliConversionCuts*)fCutArray->At(fiCut))->GetWeightForMeson(fV0Reader->GetPeriodName(),secMotherLabel, fMCStack, fInputEvent)/2.; //invariant mass is additive thus the weight for the daughters has to be devide by two for the K0s at a certain pt
2304 //cout << "MC input \t"<<i << "\t" << particle->Pt()<<"\t"<<weighted << endl;
2305 }
2306 fHistoTrueSecondaryMotherInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt(),weightedSec);
2307 fCharMesonMCInfo = 2;
2308 if (secMotherLabel >-1){
2309 if(MCStack->Particle(secMotherLabel)->GetPdgCode()==310){
2310 fCharMesonMCInfo = 4;
2311 fHistoTrueSecondaryMotherFromK0sInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt(),weightedSec);
2312 if (fDoMesonQA > 0)fHistoTrueK0sWithPi0DaughterMCPt[fiCut]->Fill(MCStack->Particle(secMotherLabel)->Pt());
2313 }
2314 if(MCStack->Particle(secMotherLabel)->GetPdgCode()==221){
2315 fCharMesonMCInfo = 3;
2316 fHistoTrueSecondaryMotherFromEtaInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt(),weightedSec);
2317 if (fDoMesonQA > 0)fHistoTrueEtaWithPi0DaughterMCPt[fiCut]->Fill(MCStack->Particle(secMotherLabel)->Pt());
2318 }
2319 if(MCStack->Particle(secMotherLabel)->GetPdgCode()==3122){
2320 fCharMesonMCInfo = 7;
2321 fHistoTrueSecondaryMotherFromLambdaInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt(),weightedSec);
2322 if (fDoMesonQA > 0)fHistoTrueLambdaWithPi0DaughterMCPt[fiCut]->Fill(MCStack->Particle(secMotherLabel)->Pt());
2323 }
2324 }
2325 } else { // Only primary pi0 for efficiency calculation
2326 fCharMesonMCInfo = 6;
2327 Float_t weighted= 1;
2328 if(((AliConversionCuts*)fCutArray->At(fiCut))->IsParticleFromBGEvent(gamma1MotherLabel, fMCStack, fInputEvent)){
2329 if (((TParticle*)MCStack->Particle(gamma1MotherLabel))->Pt()>0.005){
2330 weighted= ((AliConversionCuts*)fCutArray->At(fiCut))->GetWeightForMeson(fV0Reader->GetPeriodName(),gamma1MotherLabel, fMCStack, fInputEvent);
2331 // cout << "rec \t " <<gamma1MotherLabel << "\t" << weighted << endl;
2332 }
2333 }
2334 fHistoTruePrimaryMotherInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt(),weighted);
2335 fHistoTruePrimaryMotherW0WeightingInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt());
2336 fProfileTruePrimaryMotherWeightsInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt(),weighted);
2337
2338
2339 if (fDoMesonQA > 0){
2340 if(isTruePi0){ // Only primary pi0 for resolution
2341 fHistoTruePrimaryPi0MCPtResolPt[fiCut]->Fill(((TParticle*)MCStack->Particle(gamma1MotherLabel))->Pt(),(Pi0Candidate->Pt()-((TParticle*)MCStack->Particle(gamma1MotherLabel))->Pt())/((TParticle*)MCStack->Particle(gamma1MotherLabel))->Pt(),weighted);
2342 }
2343 if (isTrueEta){ // Only primary eta for resolution
2344 fHistoTruePrimaryEtaMCPtResolPt[fiCut]->Fill(((TParticle*)MCStack->Particle(gamma1MotherLabel))->Pt(),(Pi0Candidate->Pt()-((TParticle*)MCStack->Particle(gamma1MotherLabel))->Pt())/((TParticle*)MCStack->Particle(gamma1MotherLabel))->Pt(),weighted);
2345 }
2346 }
2347 }
2348 } else if(!isTruePi0 && !isTrueEta){ // Background
2349 if (fDoMesonQA > 0){
2350 if(gamma0MotherLabel>-1 && gamma1MotherLabel>-1){ // Both Tracks are Photons and have a mother but not Pi0 or Eta
2351 fHistoTrueBckGGInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt());
2352 fCharMesonMCInfo = 1;
2353 } else { // No photon or without mother
2354 fHistoTrueBckContInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt());
2355 }
2356 }
2357 if( isTruePi0Dalitz || isTrueEtaDalitz ){
2358 // Dalitz
2359 fCharMesonMCInfo = 5;
2360 fHistoTrueMotherDalitzInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt());
2361 } else if (gamma0DalitzCand || gamma1DalitzCand){
2362 if (fDoMesonQA > 0)fHistoTrueBckContInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt());
2363 }
2364 }
2365 }
2366 }
2367}
2368//______________________________________________________________________
2369void AliAnalysisTaskGammaConvCalo::ProcessTrueMesonCandidatesAOD(AliAODConversionMother *Pi0Candidate, AliAODConversionPhoton *TrueGammaCandidate0, AliAODConversionPhoton *TrueGammaCandidate1)
2370{
2371
2372 // Process True Mesons
2373 TClonesArray *AODMCTrackArray = dynamic_cast<TClonesArray*>(fInputEvent->FindListObject(AliAODMCParticle::StdBranchName()));
2374 Bool_t isTruePi0 = kFALSE;
2375 Bool_t isTrueEta = kFALSE;
2376 Bool_t isTruePi0Dalitz = kFALSE;
2377 Bool_t isTrueEtaDalitz = kFALSE;
2378 Bool_t gamma0DalitzCand = kFALSE;
2379 Bool_t gamma1DalitzCand = kFALSE;
2380
2381 AliAODMCParticle *positiveMC = static_cast<AliAODMCParticle*>(AODMCTrackArray->At(TrueGammaCandidate0->GetMCLabelPositive()));
2382 AliAODMCParticle *negativeMC = static_cast<AliAODMCParticle*>(AODMCTrackArray->At(TrueGammaCandidate0->GetMCLabelNegative()));
2383
2384 fCharMesonMCInfo = 0;
2385 Int_t gamma0MCLabel = -1;
2386 Int_t gamma0MotherLabel = -1;
2387 if(!positiveMC||!negativeMC)
2388 return;
2389
2390 if(positiveMC->GetMother()>-1&&(negativeMC->GetMother() == positiveMC->GetMother())){
2391 gamma0MCLabel = positiveMC->GetMother();
2392 }
2393
2394 if(gamma0MCLabel != -1){ // Gamma is Combinatorial; MC Particles don't belong to the same Mother
2395 // Daughters Gamma 0
2396 AliAODMCParticle * gammaMC0 = static_cast<AliAODMCParticle*>(AODMCTrackArray->At(gamma0MCLabel));
2397 if(abs(negativeMC->GetPdgCode())==11 && abs(positiveMC->GetPdgCode())==11){ // Electrons ...
2398 if(((positiveMC->GetMCProcessCode())) == 5 && ((negativeMC->GetMCProcessCode())) == 5){ // ... From Conversion ...
2399 if(gammaMC0->GetPdgCode() == 22){ // ... with Gamma Mother
2400 gamma0MotherLabel=gammaMC0->GetMother();
2401 }
2402 }
2403 if(gammaMC0->GetPdgCode() ==111){ // Dalitz candidate
2404 gamma0DalitzCand = kTRUE;
2405 gamma0MotherLabel=-111;
2406 }
2407 if(gammaMC0->GetPdgCode() ==221){ // Dalitz candidate
2408 gamma0DalitzCand = kTRUE;
2409 gamma0MotherLabel=-221;
2410 }
2411 }
2412 }
2413 positiveMC = static_cast<AliAODMCParticle*>(AODMCTrackArray->At(TrueGammaCandidate1->GetMCLabelPositive()));
2414 negativeMC = static_cast<AliAODMCParticle*>(AODMCTrackArray->At(TrueGammaCandidate1->GetMCLabelNegative()));
2415
2416 Int_t gamma1MCLabel = -1;
2417 Int_t gamma1MotherLabel = -1;
2418 if(!positiveMC||!negativeMC)
2419 return;
2420
2421 if(positiveMC->GetMother()>-1&&(negativeMC->GetMother() == positiveMC->GetMother())){
2422 gamma1MCLabel = positiveMC->GetMother();
2423 }
2424 if(gamma1MCLabel != -1){ // Gamma is Combinatorial; MC Particles don't belong to the same Mother
2425 // Daughters Gamma 1
2426 AliAODMCParticle * gammaMC1 = static_cast<AliAODMCParticle*>(AODMCTrackArray->At(gamma1MCLabel));
2427 if(abs(negativeMC->GetPdgCode())==11 && abs(positiveMC->GetPdgCode())==11){ // Electrons ...
2428 if(((positiveMC->GetMCProcessCode())) == 5 && ((negativeMC->GetMCProcessCode())) == 5){ // ... From Conversion ...
2429 if(gammaMC1->GetPdgCode() == 22){ // ... with Gamma Mother
2430 gamma1MotherLabel=gammaMC1->GetMother();
2431 }
2432 }
2433 if(gammaMC1->GetPdgCode() ==111 ){ // Dalitz candidate
2434 gamma1DalitzCand = kTRUE;
2435 gamma1MotherLabel=-111;
2436 }
2437 if(gammaMC1->GetPdgCode() ==221){ // Dalitz candidate
2438 gamma1DalitzCand = kTRUE;
2439 gamma1MotherLabel=-221;
2440 }
2441 }
2442 }
2443 if(gamma0MotherLabel>=0 && gamma0MotherLabel==gamma1MotherLabel){
2444 if(static_cast<AliAODMCParticle*>(AODMCTrackArray->At(gamma1MotherLabel))->GetPdgCode() == 111){
2445 isTruePi0=kTRUE;
2446 }
2447 if(static_cast<AliAODMCParticle*>(AODMCTrackArray->At(gamma1MotherLabel))->GetPdgCode() == 221){
2448 isTrueEta=kTRUE;
2449 }
2450 }
2451
2452 //Identify Dalitz candidate
2453 if (gamma1DalitzCand || gamma0DalitzCand){
2454 if (gamma0DalitzCand && gamma0MCLabel >=0 && gamma0MCLabel==gamma1MotherLabel){
2455 if (gamma0MotherLabel == -111) isTruePi0Dalitz = kTRUE;
2456 if (gamma0MotherLabel == -221) isTrueEtaDalitz = kTRUE;
2457 }
2458 if (gamma1DalitzCand && gamma1MCLabel >=0 && gamma1MCLabel==gamma0MotherLabel){
2459 if (gamma1MotherLabel == -111) isTruePi0Dalitz = kTRUE;
2460 if (gamma1MotherLabel == -221) isTrueEtaDalitz = kTRUE;
2461 }
2462 }
2463
2464 if(isTruePi0 || isTrueEta){// True Pion or Eta
2465 fHistoTrueMotherInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt());
2466 if (fDoMesonQA > 0){
2467 if (isTruePi0){
2468 if ( Pi0Candidate->M() > 0.05 && Pi0Candidate->M() < 0.17){
2469 fHistoTruePi0PtY[fiCut]->Fill(Pi0Candidate->Pt(),Pi0Candidate->Rapidity()-((AliConversionCuts*)fCutArray->At(fiCut))->GetEtaShift());
2470 fHistoTruePi0PtAlpha[fiCut]->Fill(Pi0Candidate->Pt(),Pi0Candidate->GetAlpha());
2471 fHistoTruePi0PtOpenAngle[fiCut]->Fill(Pi0Candidate->Pt(),Pi0Candidate->GetOpeningAngle());
2472 }
2473 } else if (isTrueEta){
2474 if ( Pi0Candidate->M() > 0.45 && Pi0Candidate->M() < 0.65){
2475 fHistoTrueEtaPtY[fiCut]->Fill(Pi0Candidate->Pt(),Pi0Candidate->Rapidity()-((AliConversionCuts*)fCutArray->At(fiCut))->GetEtaShift());
2476 fHistoTrueEtaPtAlpha[fiCut]->Fill(Pi0Candidate->Pt(),Pi0Candidate->GetAlpha());
2477 fHistoTrueEtaPtOpenAngle[fiCut]->Fill(Pi0Candidate->Pt(),Pi0Candidate->GetOpeningAngle());
2478 }
2479 }
2480 }
2481 if(!(static_cast<AliAODMCParticle*>(AODMCTrackArray->At(gamma0MotherLabel))->IsPrimary())){ // Secondary Meson
2482 Int_t secMotherLabel = static_cast<AliAODMCParticle*>(AODMCTrackArray->At(gamma1MotherLabel))->GetMother();
2483 Float_t weightedSec= 1;
2484 if(((AliConversionCuts*)fCutArray->At(fiCut))->IsParticleFromBGEvent(secMotherLabel, 0x0, fInputEvent) && static_cast<AliAODMCParticle*>(AODMCTrackArray->At(secMotherLabel))->GetPdgCode()==310){
2485 weightedSec= ((AliConversionCuts*)fCutArray->At(fiCut))->GetWeightForMeson(fV0Reader->GetPeriodName(),secMotherLabel, 0x0, fInputEvent)/2.; //invariant mass is additive thus the weight for the daughters has to be devide by two for the K0s at a certain pt
2486 //cout << "MC input \t"<<i << "\t" << particle->Pt()<<"\t"<<weighted << endl;
2487 }
2488 fHistoTrueSecondaryMotherInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt(),weightedSec);
2489 fCharMesonMCInfo = 2;
2490 if (secMotherLabel >-1){
2491 if(static_cast<AliAODMCParticle*>(AODMCTrackArray->At(secMotherLabel))->GetPdgCode()==310){
2492 fCharMesonMCInfo = 4;
2493 fHistoTrueSecondaryMotherFromK0sInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt(),weightedSec);
2494 if (fDoMesonQA > 0)fHistoTrueK0sWithPi0DaughterMCPt[fiCut]->Fill(static_cast<AliAODMCParticle*>(AODMCTrackArray->At(secMotherLabel))->Pt());
2495 }
2496 if(static_cast<AliAODMCParticle*>(AODMCTrackArray->At(secMotherLabel))->GetPdgCode()==221){
2497 fCharMesonMCInfo = 3;
2498 fHistoTrueSecondaryMotherFromEtaInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt(),weightedSec);
2499 if (fDoMesonQA > 0)fHistoTrueEtaWithPi0DaughterMCPt[fiCut]->Fill(static_cast<AliAODMCParticle*>(AODMCTrackArray->At(secMotherLabel))->Pt());
2500 }
2501 if(static_cast<AliAODMCParticle*>(AODMCTrackArray->At(secMotherLabel))->GetPdgCode()==3122){
2502 fCharMesonMCInfo = 7;
2503 fHistoTrueSecondaryMotherFromLambdaInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt(),weightedSec);
2504 if (fDoMesonQA > 0)fHistoTrueLambdaWithPi0DaughterMCPt[fiCut]->Fill(static_cast<AliAODMCParticle*>(AODMCTrackArray->At(secMotherLabel))->Pt());
2505 }
2506 }
2507 }else{ // Only primary pi0 for efficiency calculation
2508 Float_t weighted= 1;
2509 fCharMesonMCInfo = 6;
2510 if(((AliConversionCuts*)fCutArray->At(fiCut))->IsParticleFromBGEvent(gamma1MotherLabel, 0x0, fInputEvent)){
2511 if (static_cast<AliAODMCParticle*>(AODMCTrackArray->At(gamma1MotherLabel))->Pt()>0.005){
2512 weighted= ((AliConversionCuts*)fCutArray->At(fiCut))->GetWeightForMeson(fV0Reader->GetPeriodName(),gamma1MotherLabel, 0x0, fInputEvent);
2513 // cout << "rec \t " <<gamma1MotherLabel << "\t" << weighted << endl;
2514 }
2515 }
2516 fHistoTruePrimaryMotherInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt(),weighted);
2517 fHistoTruePrimaryMotherW0WeightingInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt());
2518 fProfileTruePrimaryMotherWeightsInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt(),weighted);
2519
2520 if (fDoMesonQA > 0){
2521 if(isTruePi0){ // Only primary pi0 for resolution
2522 fHistoTruePrimaryPi0MCPtResolPt[fiCut]->Fill(static_cast<AliAODMCParticle*>(AODMCTrackArray->At(gamma1MotherLabel))->Pt(),
2523 (Pi0Candidate->Pt()-static_cast<AliAODMCParticle*>(AODMCTrackArray->At(gamma1MotherLabel))->Pt())/static_cast<AliAODMCParticle*>(AODMCTrackArray->At(gamma1MotherLabel))->Pt(),weighted);
2524
2525 }
2526 if (isTrueEta){ // Only primary eta for resolution
2527 fHistoTruePrimaryEtaMCPtResolPt[fiCut]->Fill(static_cast<AliAODMCParticle*>(AODMCTrackArray->At(gamma1MotherLabel))->Pt(),
2528 (Pi0Candidate->Pt()-static_cast<AliAODMCParticle*>(AODMCTrackArray->At(gamma1MotherLabel))->Pt())/static_cast<AliAODMCParticle*>(AODMCTrackArray->At(gamma1MotherLabel))->Pt(),weighted);
2529 }
2530 }
2531 }
2532 } else if(!isTruePi0 && !isTrueEta) { // Background
2533 if (fDoMesonQA > 0){
2534 if(gamma0MotherLabel>-1 && gamma1MotherLabel>-1){ // Both Tracks are Photons and have a mother but not Pi0 or Eta
2535 fHistoTrueBckGGInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt());
2536 fCharMesonMCInfo = 1;
2537 } else { // No photon or without mother
2538 fHistoTrueBckContInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt());
2539 }
2540 }
2541 if( isTruePi0Dalitz || isTrueEtaDalitz ){
2542 // Dalitz
2543 fCharMesonMCInfo = 5;
2544 fHistoTrueMotherDalitzInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt());
2545 } else if (gamma0DalitzCand || gamma1DalitzCand){
2546 if (fDoMesonQA > 0)fHistoTrueBckContInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt());
2547 }
2548 }
2549}
2550//________________________________________________________________________
2551void AliAnalysisTaskGammaConvCalo::CalculateBackground(){
2552
2553 Int_t zbin= fBGClusHandler[fiCut]->GetZBinIndex(fInputEvent->GetPrimaryVertex()->GetZ());
2554 Int_t mbin = 0;
2555
2556 if(((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->UseTrackMultiplicity()){
2557 mbin = fBGClusHandler[fiCut]->GetMultiplicityBinIndex(fV0Reader->GetNumberOfPrimaryTracks());
2558 } else {
2559 mbin = fBGClusHandler[fiCut]->GetMultiplicityBinIndex(fGammaCandidates->GetEntries());
2560 }
2561
2562 AliGammaConversionAODBGHandler::GammaConversionVertex *bgEventVertex = NULL;
2563 if(((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->UseTrackMultiplicity()){
2564 for(Int_t nEventsInBG=0;nEventsInBG<fBGClusHandler[fiCut]->GetNBGEvents();nEventsInBG++){
2565 AliGammaConversionAODVector *previousEventV0s = fBGClusHandler[fiCut]->GetBGGoodV0s(zbin,mbin,nEventsInBG);
2566 if(fMoveParticleAccordingToVertex == kTRUE){
2567 bgEventVertex = fBGClusHandler[fiCut]->GetBGEventVertex(zbin,mbin,nEventsInBG);
2568 }
2569
2570 for(Int_t iCurrent=0;iCurrent<fGammaCandidates->GetEntries();iCurrent++){
2571 AliAODConversionPhoton currentEventGoodV0 = *(AliAODConversionPhoton*)(fGammaCandidates->At(iCurrent));
2572 for(UInt_t iPrevious=0;iPrevious<previousEventV0s->size();iPrevious++){
2573 AliAODConversionPhoton previousGoodV0 = (AliAODConversionPhoton)(*(previousEventV0s->at(iPrevious)));
2574 if(fMoveParticleAccordingToVertex == kTRUE){
2575 MoveParticleAccordingToVertex(&previousGoodV0,bgEventVertex);
2576 }
2577 if(((AliConversionCuts*)fCutArray->At(fiCut))->GetInPlaneOutOfPlaneCut() != 0){
2578 RotateParticleAccordingToEP(&previousGoodV0,bgEventVertex->fEP,fEventPlaneAngle);
2579 }
2580
2581 AliAODConversionMother *backgroundCandidate = new AliAODConversionMother(&currentEventGoodV0,&previousGoodV0);
2582 backgroundCandidate->CalculateDistanceOfClossetApproachToPrimVtx(fInputEvent->GetPrimaryVertex());
2583 if((((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))
2584 ->MesonIsSelected(backgroundCandidate,kFALSE,((AliConversionCuts*)fCutArray->At(fiCut))->GetEtaShift()))){
2585 fHistoMotherBackInvMassPt[fiCut]->Fill(backgroundCandidate->M(),backgroundCandidate->Pt());
2586 Double_t sparesFill[4] = {backgroundCandidate->M(),backgroundCandidate->Pt(),(Double_t)zbin,(Double_t)mbin};
2587 fSparseMotherBackInvMassPtZM[fiCut]->Fill(sparesFill,1);
2588 }
2589 delete backgroundCandidate;
2590 backgroundCandidate = 0x0;
2591 }
2592 }
2593 }
2594 } else {
2595 for(Int_t nEventsInBG=0;nEventsInBG <fBGClusHandler[fiCut]->GetNBGEvents();nEventsInBG++){
2596 AliGammaConversionAODVector *previousEventV0s = fBGClusHandler[fiCut]->GetBGGoodV0s(zbin,mbin,nEventsInBG);
2597 if(previousEventV0s){
2598 if(fMoveParticleAccordingToVertex == kTRUE){
2599 bgEventVertex = fBGClusHandler[fiCut]->GetBGEventVertex(zbin,mbin,nEventsInBG);
2600 }
2601 for(Int_t iCurrent=0;iCurrent<fGammaCandidates->GetEntries();iCurrent++){
2602 AliAODConversionPhoton currentEventGoodV0 = *(AliAODConversionPhoton*)(fGammaCandidates->At(iCurrent));
2603 for(UInt_t iPrevious=0;iPrevious<previousEventV0s->size();iPrevious++){
2604
2605 AliAODConversionPhoton previousGoodV0 = (AliAODConversionPhoton)(*(previousEventV0s->at(iPrevious)));
2606
2607 if(fMoveParticleAccordingToVertex == kTRUE){
2608 MoveParticleAccordingToVertex(&previousGoodV0,bgEventVertex);
2609 }
2610 if(((AliConversionCuts*)fCutArray->At(fiCut))->GetInPlaneOutOfPlaneCut() != 0){
2611 RotateParticleAccordingToEP(&previousGoodV0,bgEventVertex->fEP,fEventPlaneAngle);
2612 }
2613
2614
2615 AliAODConversionMother *backgroundCandidate = new AliAODConversionMother(&currentEventGoodV0,&previousGoodV0);
2616 backgroundCandidate->CalculateDistanceOfClossetApproachToPrimVtx(fInputEvent->GetPrimaryVertex());
2617 if((((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->MesonIsSelected(backgroundCandidate,kFALSE,((AliConversionCuts*)fCutArray->At(fiCut))->GetEtaShift()))){
2618 fHistoMotherBackInvMassPt[fiCut]->Fill(backgroundCandidate->M(),backgroundCandidate->Pt());
2619 Double_t sparesFill[4] = {backgroundCandidate->M(),backgroundCandidate->Pt(),(Double_t)zbin,(Double_t)mbin};
2620 fSparseMotherBackInvMassPtZM[fiCut]->Fill(sparesFill,1);
2621 }
2622 delete backgroundCandidate;
2623 backgroundCandidate = 0x0;
2624 }
2625 }
2626 }
2627 }
2628 }
2629}
2630
2631//________________________________________________________________________
2632void AliAnalysisTaskGammaConvCalo::CalculateBackgroundRP(){
2633
2634 Int_t zbin= fBGHandlerRP[fiCut]->GetZBinIndex(fInputEvent->GetPrimaryVertex()->GetZ());
2635 Int_t mbin = 0;
2636 if(((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->UseTrackMultiplicity()){
2637 mbin = fBGHandlerRP[fiCut]->GetMultiplicityBinIndex(fV0Reader->GetNumberOfPrimaryTracks());
2638 } else {
2639 mbin = fBGHandlerRP[fiCut]->GetMultiplicityBinIndex(fGammaCandidates->GetEntries());
2640 }
2641
2642
2643 //Rotation Method
2644 if(((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->UseRotationMethod()){
2645 // Correct for the number of rotations
2646 // BG is for rotation the same, except for factor NRotations
2647 Double_t weight=1./Double_t(((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->GetNumberOfBGEvents());
2648
2649 for(Int_t firstGammaIndex=0;firstGammaIndex<fGammaCandidates->GetEntries();firstGammaIndex++){
2650
2651 AliAODConversionPhoton *gamma0=dynamic_cast<AliAODConversionPhoton*>(fGammaCandidates->At(firstGammaIndex));
2652 if (gamma0==NULL) continue;
2653 for(Int_t secondGammaIndex=firstGammaIndex+1;secondGammaIndex<fGammaCandidates->GetEntries();secondGammaIndex++){
2654 AliAODConversionPhoton *gamma1=dynamic_cast<AliAODConversionPhoton*>(fGammaCandidates->At(secondGammaIndex));
2655 if (gamma1 == NULL) continue;
2656 if(!((AliConversionCuts*)fCutArray->At(fiCut))->PhotonIsSelected(gamma1,fInputEvent))continue;
2657 for(Int_t nRandom=0;nRandom<((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->GetNumberOfBGEvents();nRandom++){
2658
2659 RotateParticle(gamma1);
2660 AliAODConversionMother backgroundCandidate(gamma0,gamma1);
2661 backgroundCandidate.CalculateDistanceOfClossetApproachToPrimVtx(fInputEvent->GetPrimaryVertex());
2662 if(((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))
2663 ->MesonIsSelected(&backgroundCandidate,kFALSE,((AliConversionCuts*)fCutArray->At(fiCut))->GetEtaShift())){
2664 fHistoMotherBackInvMassPt[fiCut]->Fill(backgroundCandidate.M(),backgroundCandidate.Pt());
2665 Double_t sparesFill[4] = {backgroundCandidate.M(),backgroundCandidate.Pt(),(Double_t)zbin,(Double_t)mbin};
2666 fSparseMotherBackInvMassPtZM[fiCut]->Fill(sparesFill,weight);
2667 }
2668 }
2669 }
2670 }
2671 }
2672 else{
2673 // Do Event Mixing
2674 for(Int_t nEventsInBG=0;nEventsInBG <fBGHandlerRP[fiCut]->GetNBGEvents(fGammaCandidates,fInputEvent);nEventsInBG++){
2675
2676 AliGammaConversionPhotonVector *previousEventGammas = fBGHandlerRP[fiCut]->GetBGGoodGammas(fGammaCandidates,fInputEvent,nEventsInBG);
2677
2678 if(previousEventGammas){
2679 // test weighted background
2680 Double_t weight=1.0;
2681 // Correct for the number of eventmixing:
2682 // N gammas -> (N-1) + (N-2) +(N-3) ...+ (N-(N-1)) using sum formula sum(i)=N*(N-1)/2 -> N*(N-1)/2
2683 // real combinations (since you cannot combine a photon with its own)
2684 // but BG leads to N_{a}*N_{b} combinations
2685 weight*=0.5*(Double_t(fGammaCandidates->GetEntries()-1))/Double_t(previousEventGammas->size());
2686
2687 for(Int_t iCurrent=0;iCurrent<fGammaCandidates->GetEntries();iCurrent++){
2688 AliAODConversionPhoton *gamma0 = (AliAODConversionPhoton*)(fGammaCandidates->At(iCurrent));
2689 for(UInt_t iPrevious=0;iPrevious<previousEventGammas->size();iPrevious++){
2690
2691 AliAODConversionPhoton *gamma1 = (AliAODConversionPhoton*)(previousEventGammas->at(iPrevious));
2692
2693 AliAODConversionMother backgroundCandidate(gamma0,gamma1);
2694 backgroundCandidate.CalculateDistanceOfClossetApproachToPrimVtx(fInputEvent->GetPrimaryVertex());
2695 if(((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))
2696 ->MesonIsSelected(&backgroundCandidate,kFALSE,((AliConversionCuts*)fCutArray->At(fiCut))->GetEtaShift())){
2697 fHistoMotherBackInvMassPt[fiCut]->Fill(backgroundCandidate.M(),backgroundCandidate.Pt());
2698 Double_t sparesFill[4] = {backgroundCandidate.M(),backgroundCandidate.Pt(),(Double_t)zbin,(Double_t)mbin};
2699 fSparseMotherBackInvMassPtZM[fiCut]->Fill(sparesFill,weight);
2700 }
2701 }
2702 }
2703 }
2704 }
2705 }
2706}
2707//________________________________________________________________________
2708void AliAnalysisTaskGammaConvCalo::RotateParticle(AliAODConversionPhoton *gamma){
2709 Int_t fNDegreesPMBackground= ((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->NDegreesRotation();
2710 Double_t nRadiansPM = fNDegreesPMBackground*TMath::Pi()/180;
2711 Double_t rotationValue = fRandom.Rndm()*2*nRadiansPM + TMath::Pi()-nRadiansPM;
2712 gamma->RotateZ(rotationValue);
2713}
2714
2715//________________________________________________________________________
2716void AliAnalysisTaskGammaConvCalo::RotateParticleAccordingToEP(AliAODConversionPhoton *gamma, Double_t previousEventEP, Double_t thisEventEP){
2717
2718 previousEventEP=previousEventEP+TMath::Pi();
2719 thisEventEP=thisEventEP+TMath::Pi();
2720 Double_t rotationValue= thisEventEP-previousEventEP;
2721 gamma->RotateZ(rotationValue);
2722}
2723
2724//________________________________________________________________________
2725void AliAnalysisTaskGammaConvCalo::MoveParticleAccordingToVertex(AliAODConversionPhoton* particle,const AliGammaConversionAODBGHandler::GammaConversionVertex *vertex){
2726 //see header file for documentation
2727
2728 Double_t dx = vertex->fX - fInputEvent->GetPrimaryVertex()->GetX();
2729 Double_t dy = vertex->fY - fInputEvent->GetPrimaryVertex()->GetY();
2730 Double_t dz = vertex->fZ - fInputEvent->GetPrimaryVertex()->GetZ();
2731
2732 Double_t movedPlace[3] = {particle->GetConversionX() - dx,particle->GetConversionY() - dy,particle->GetConversionZ() - dz};
2733 particle->SetConversionPoint(movedPlace);
2734}
2735//________________________________________________________________________
2736void AliAnalysisTaskGammaConvCalo::UpdateEventByEventData(){
2737 //see header file for documentation
2738 if(fGammaCandidates->GetEntries() >0 ){
2739 if(((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->UseTrackMultiplicity()){
2740 fBGHandler[fiCut]->AddEvent(fGammaCandidates,fInputEvent->GetPrimaryVertex()->GetX(),fInputEvent->GetPrimaryVertex()->GetY(),fInputEvent->GetPrimaryVertex()->GetZ(),fV0Reader->GetNumberOfPrimaryTracks(),fEventPlaneAngle);
2741 fBGClusHandler[fiCut]->AddEvent(fClusterCandidates,fInputEvent->GetPrimaryVertex()->GetX(),fInputEvent->GetPrimaryVertex()->GetY(),fInputEvent->GetPrimaryVertex()->GetZ(),fV0Reader->GetNumberOfPrimaryTracks(),fEventPlaneAngle);
2742 } else { // means we use #V0s for multiplicity
2743 fBGHandler[fiCut]->AddEvent(fGammaCandidates,fInputEvent->GetPrimaryVertex()->GetX(),fInputEvent->GetPrimaryVertex()->GetY(),fInputEvent->GetPrimaryVertex()->GetZ(),fGammaCandidates->GetEntries(),fEventPlaneAngle);
2744 fBGClusHandler[fiCut]->AddEvent(fClusterCandidates,fInputEvent->GetPrimaryVertex()->GetX(),fInputEvent->GetPrimaryVertex()->GetY(),fInputEvent->GetPrimaryVertex()->GetZ(),fGammaCandidates->GetEntries(),fEventPlaneAngle);
2745 }
2746 }
2747}
2748
2749
2750//________________________________________________________________________
2751void AliAnalysisTaskGammaConvCalo::FillPhotonCombinatorialBackgroundHist(AliAODConversionPhoton *TruePhotonCandidate, Int_t pdgCode[])
2752{
2753 // Combinatorial Bck = 0 ee, 1 ep,i 2 ek, 3 ep, 4 emu, 5 pipi, 6 pik, 7 pip, 8 pimu, 9 kk, 10 kp, 11 kmu, 12 pp, 13 pmu, 14 mumu, 15 Rest
2754 if(pdgCode[0]==11 && pdgCode[1]==11){if(fIsFromMBHeader)fHistoCombinatorialPt[fiCut]->Fill(TruePhotonCandidate->Pt(),0);}
2755 else if( (pdgCode[0]==11 && pdgCode[1]==211) || (pdgCode[0]==211 && pdgCode[1]==11) )
2756 {if(fIsFromMBHeader)fHistoCombinatorialPt[fiCut]->Fill(TruePhotonCandidate->Pt(),1);}
2757 else if( (pdgCode[0]==11 && pdgCode[1]==321) || (pdgCode[0]==321 && pdgCode[1]==11) )
2758 {if(fIsFromMBHeader)fHistoCombinatorialPt[fiCut]->Fill(TruePhotonCandidate->Pt(),2);}
2759 else if( (pdgCode[0]==11 && pdgCode[1]==2212) || (pdgCode[0]==2212 && pdgCode[1]==11) )
2760 {if(fIsFromMBHeader)fHistoCombinatorialPt[fiCut]->Fill(TruePhotonCandidate->Pt(),3);}
2761 else if( (pdgCode[0]==11 && pdgCode[1]==13) || (pdgCode[0]==13 && pdgCode[1]==11) )
2762 {if(fIsFromMBHeader)fHistoCombinatorialPt[fiCut]->Fill(TruePhotonCandidate->Pt(),4);}
2763 else if( pdgCode[0]==211 && pdgCode[1]==211 ){if(fIsFromMBHeader)fHistoCombinatorialPt[fiCut]->Fill(TruePhotonCandidate->Pt(),5);}
2764 else if( (pdgCode[0]==211 && pdgCode[1]==321) || (pdgCode[0]==321 && pdgCode[1]==211) )
2765 {if(fIsFromMBHeader)fHistoCombinatorialPt[fiCut]->Fill(TruePhotonCandidate->Pt(),6);}
2766 else if( (pdgCode[0]==211 && pdgCode[1]==2212) || (pdgCode[0]==2212 && pdgCode[1]==211) )
2767 {if(fIsFromMBHeader)fHistoCombinatorialPt[fiCut]->Fill(TruePhotonCandidate->Pt(),7);}
2768 else if( (pdgCode[0]==211 && pdgCode[1]==13) || (pdgCode[0]==13 && pdgCode[1]==211) )
2769 {if(fIsFromMBHeader)fHistoCombinatorialPt[fiCut]->Fill(TruePhotonCandidate->Pt(),8);}
2770 else if( pdgCode[0]==321 && pdgCode[1]==321 ){if(fIsFromMBHeader)fHistoCombinatorialPt[fiCut]->Fill(TruePhotonCandidate->Pt(),9);}
2771 else if( (pdgCode[0]==321 && pdgCode[1]==2212) || (pdgCode[0]==2212 && pdgCode[1]==321) )
2772 {if(fIsFromMBHeader)fHistoCombinatorialPt[fiCut]->Fill(TruePhotonCandidate->Pt(),10);}
2773 else if( (pdgCode[0]==321 && pdgCode[1]==13) || (pdgCode[0]==13 && pdgCode[1]==321) )
2774 {if(fIsFromMBHeader)fHistoCombinatorialPt[fiCut]->Fill(TruePhotonCandidate->Pt(),11);}
2775 else if( pdgCode[0]==2212 && pdgCode[1]==2212 ){if(fIsFromMBHeader)fHistoCombinatorialPt[fiCut]->Fill(TruePhotonCandidate->Pt(),12);}
2776 else if( (pdgCode[0]==2212 && pdgCode[1]==13) || (pdgCode[0]==13 && pdgCode[1]==2212) )
2777 {if(fIsFromMBHeader)fHistoCombinatorialPt[fiCut]->Fill(TruePhotonCandidate->Pt(),13);}
2778 else if( pdgCode[0]==13 && pdgCode[1]==13 ){if(fIsFromMBHeader)fHistoCombinatorialPt[fiCut]->Fill(TruePhotonCandidate->Pt(),14);}
2779 else {if(fIsFromMBHeader)fHistoCombinatorialPt[fiCut]->Fill(TruePhotonCandidate->Pt(),15);}
2780}
2781//________________________________________________________________________
2782void AliAnalysisTaskGammaConvCalo::RelabelAODPhotonCandidates(Bool_t mode){
2783
2784 // Relabeling For AOD Event
2785 // ESDiD -> AODiD
2786 // MCLabel -> AODMCLabel
2787
2788 if(mode){
2789 fMCStackPos = new Int_t[fReaderGammas->GetEntries()];
2790 fMCStackNeg = new Int_t[fReaderGammas->GetEntries()];
2791 fESDArrayPos = new Int_t[fReaderGammas->GetEntries()];
2792 fESDArrayNeg = new Int_t[fReaderGammas->GetEntries()];
2793 }
2794
2795 for(Int_t iGamma = 0;iGamma<fReaderGammas->GetEntries();iGamma++){
2796 AliAODConversionPhoton* PhotonCandidate = (AliAODConversionPhoton*) fReaderGammas->At(iGamma);
2797 if(!PhotonCandidate) continue;
2798 if(!mode){// Back to ESD Labels
2799 PhotonCandidate->SetMCLabelPositive(fMCStackPos[iGamma]);
2800 PhotonCandidate->SetMCLabelNegative(fMCStackNeg[iGamma]);
2801 PhotonCandidate->SetLabelPositive(fESDArrayPos[iGamma]);
2802 PhotonCandidate->SetLabelNegative(fESDArrayNeg[iGamma]);
2803 continue;
2804 }
2805 fMCStackPos[iGamma] = PhotonCandidate->GetMCLabelPositive();
2806 fMCStackNeg[iGamma] = PhotonCandidate->GetMCLabelNegative();
2807 fESDArrayPos[iGamma] = PhotonCandidate->GetTrackLabelPositive();
2808 fESDArrayNeg[iGamma] = PhotonCandidate->GetTrackLabelNegative();
2809
2810 Bool_t AODLabelPos = kFALSE;
2811 Bool_t AODLabelNeg = kFALSE;
2812
2813 for(Int_t i = 0; i<fInputEvent->GetNumberOfTracks();i++){
2814 AliAODTrack *tempDaughter = static_cast<AliAODTrack*>(fInputEvent->GetTrack(i));
2815 if(!AODLabelPos){
2816 if( tempDaughter->GetID() == PhotonCandidate->GetTrackLabelPositive() ){
2817 PhotonCandidate->SetMCLabelPositive(abs(tempDaughter->GetLabel()));
2818 PhotonCandidate->SetLabelPositive(i);
2819 AODLabelPos = kTRUE;
2820 }
2821 }
2822 if(!AODLabelNeg){
2823 if( tempDaughter->GetID() == PhotonCandidate->GetTrackLabelNegative()){
2824 PhotonCandidate->SetMCLabelNegative(abs(tempDaughter->GetLabel()));
2825 PhotonCandidate->SetLabelNegative(i);
2826 AODLabelNeg = kTRUE;
2827 }
2828 }
2829 if(AODLabelNeg && AODLabelPos){
2830 break;
2831 }
2832 }
2833 if(!AODLabelPos || !AODLabelNeg){
2834 cout<<"WARNING!!! AOD TRACKS NOT FOUND FOR"<<endl;
2835 }
2836 }
2837
2838
2839 if(!mode){
2840 delete[] fMCStackPos;
2841 delete[] fMCStackNeg;
2842 delete[] fESDArrayPos;
2843 delete[] fESDArrayNeg;
2844 }
2845}
2846
2847void AliAnalysisTaskGammaConvCalo::SetLogBinningXTH2(TH2* histoRebin){
2848 TAxis *axisafter = histoRebin->GetXaxis();
2849 Int_t bins = axisafter->GetNbins();
2850 Double_t from = axisafter->GetXmin();
2851 Double_t to = axisafter->GetXmax();
2852 Double_t *newbins = new Double_t[bins+1];
2853 newbins[0] = from;
2854 Double_t factor = TMath::Power(to/from, 1./bins);
2855 for(Int_t i=1; i<=bins; ++i) newbins[i] = factor * newbins[i-1];
2856 axisafter->Set(bins, newbins);
2857 delete [] newbins;
2858}
2859
2860//________________________________________________________________________
2861void AliAnalysisTaskGammaConvCalo::Terminate(const Option_t *)
2862{
2863
2864 //fOutputContainer->Print(); // Will crash on GRID
2865}
2866
2867//________________________________________________________________________
2868Int_t AliAnalysisTaskGammaConvCalo::GetSourceClassification(Int_t daughter, Int_t pdgCode){
2869
2870 if (daughter == 111) {
2871 if (abs(pdgCode) == 310) return 1; // k0s
2872 else if (abs(pdgCode) == 3122) return 2; // Lambda
2873 else if (abs(pdgCode) == 130) return 3; // K0L
2874 else if (abs(pdgCode) == 2212) return 4; // proton
2875 else if (abs(pdgCode) == 2112) return 5; // neutron
2876 else if (abs(pdgCode) == 211) return 6; // pion
2877 else if (abs(pdgCode) == 321) return 7; // kaon
2878 else if (abs(pdgCode) == 113 || abs(pdgCode) == 213 ) return 8; // rho 0,+,-
2879 else if (abs(pdgCode) == 3222 || abs(pdgCode) == 3212 || abs(pdgCode) == 3112 ) return 9; // Sigma
2880 else if (abs(pdgCode) == 2224 || abs(pdgCode) == 2214 || abs(pdgCode) == 2114 || abs(pdgCode) == 1114 ) return 10; // Delta
2881 else if (abs(pdgCode) == 313 || abs(pdgCode) == 323 ) return 11; // K*
2882 else return 15;
2883 }
2884 return 15;
2885
2886}
2887
2888// //________________________________________________________________________
2889// Double_t AliAnalysisTaskGammaConvCalo::GetMaxCellEnergy(const AliVCluster *cluster, Short_t &id) const
2890// {
2891// // Get maximum energy of attached cell.
2892//
2893// id = -1;
2894// Double_t maxe = 0;
2895// Int_t ncells = cluster->GetNCells();
2896// if (fEsdCells) {
2897// for (Int_t i=0; i<ncells; i++) {
2898// Double_t e = fEsdCells->GetCellAmplitude(TMath::Abs(cluster->GetCellAbsId(i)));
2899// if (e>maxe) {
2900// maxe = e;
2901// id = cluster->GetCellAbsId(i);
2902// }
2903// }
2904// } else {
2905// for (Int_t i=0; i<ncells; i++) {
2906// Double_t e = fAodCells->GetCellAmplitude(TMath::Abs(cluster->GetCellAbsId(i)));
2907// if (e>maxe)
2908// maxe = e;
2909// id = cluster->GetCellAbsId(i);
2910// }
2911// }
2912// return maxe;
2913// }