d31d7f3628a7827ca2c1b39a63c8bbd256ffe346
[u/mrichter/AliRoot.git] / PWGGA / GammaConv / AliAnalysisTaskGammaCalo.cxx
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 "AliAnalysisTaskGammaCalo.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 "AliAODMCParticle.h"
50 #include "AliAODMCHeader.h"
51 #include "AliEventplane.h"
52 #include "AliAnalysisTaskEMCALClusterizeFast.h"
53 #include "AliAODEvent.h"
54 #include "AliAODInputHandler.h"
55 #include "AliESDEvent.h"
56 #include "AliESDInputHandler.h"
57 #include "AliInputEventHandler.h"
58
59 ClassImp(AliAnalysisTaskGammaCalo)
60
61 //________________________________________________________________________
62 AliAnalysisTaskGammaCalo::AliAnalysisTaskGammaCalo(): AliAnalysisTaskSE(),
63         fV0Reader(NULL),
64         fBGHandler(NULL),
65         fInputEvent(NULL),
66         fMCEvent(NULL),
67         fMCStack(NULL),
68         fCutFolder(NULL),
69         fESDList(NULL),
70         fBackList(NULL),
71         fMotherList(NULL),
72         fTrueList(NULL),
73         fMCList(NULL),
74         fHeaderNameList(NULL),
75         fOutputContainer(NULL),
76         fClusterCandidates(NULL),
77         fEventCutArray(NULL),
78         fEventCuts(NULL),
79         fClusterCutArray(NULL),
80         fCaloPhotonCuts(NULL),
81         fMesonCutArray(NULL),
82         fMesonCuts(NULL),
83         fHistoMotherInvMassPt(NULL),
84         fHistoMotherInvMass3ClusterPt(NULL),
85         fSparseMotherInvMassPtZM(NULL),
86         fHistoMotherBackInvMassPt(NULL),
87         fSparseMotherBackInvMassPtZM(NULL),
88         fHistoMotherInvMassEalpha(NULL),
89         fHistoMotherPi0PtY(NULL),
90         fHistoMotherEtaPtY(NULL),
91         fHistoMotherPi0PtAlpha(NULL),
92         fHistoMotherEtaPtAlpha(NULL),
93         fHistoMotherPi0PtOpenAngle(NULL),
94         fHistoMotherEtaPtOpenAngle(NULL),
95         fHistoMotherInvMassECalib(NULL),
96         fHistoMotherInvMassECalibalpha(NULL),
97         fHistoClusGammaPt(NULL),
98         fHistoMCHeaders(NULL),
99         fHistoMCAllGammaPt(NULL),
100         fHistoMCDecayGammaPi0Pt(NULL),
101         fHistoMCDecayGammaRhoPt(NULL),
102         fHistoMCDecayGammaEtaPt(NULL),
103         fHistoMCDecayGammaOmegaPt(NULL),
104         fHistoMCDecayGammaEtapPt(NULL),
105         fHistoMCDecayGammaPhiPt(NULL),
106         fHistoMCDecayGammaSigmaPt(NULL),
107         fHistoMCPi0Pt(NULL),
108         fHistoMCPi0WOWeightPt(NULL),
109         fHistoMCEtaPt(NULL),
110         fHistoMCEtaWOWeightPt(NULL),
111         fHistoMCPi0InAccPt(NULL),
112         fHistoMCEtaInAccPt(NULL),
113         fHistoMCPi0PtY(NULL),
114         fHistoMCEtaPtY(NULL),
115         fHistoMCK0sPt(NULL),
116         fHistoMCK0sWOWeightPt(NULL),
117         fHistoMCK0sPtY(NULL),
118         fHistoMCSecPi0PtvsSource(NULL),
119         fHistoMCSecPi0Source(NULL),
120         fHistoMCSecEtaPt(NULL),
121         fHistoMCSecEtaSource(NULL),
122         fHistoTruePi0InvMassPt(NULL),
123         fHistoTrueEtaInvMassPt(NULL),
124         fHistoTruePi0CaloPhotonInvMassPt(NULL),
125         fHistoTrueEtaCaloPhotonInvMassPt(NULL),
126         fHistoTruePi0CaloConvertedPhotonInvMassPt(NULL),
127         fHistoTrueEtaCaloConvertedPhotonInvMassPt(NULL),
128         fHistoTruePi0CaloMixedPhotonConvPhotonInvMassPt(NULL),
129         fHistoTrueEtaCaloMixedPhotonConvPhotonInvMassPt(NULL),
130         fHistoTruePi0CaloElectronInvMassPt(NULL),
131         fHistoTrueEtaCaloElectronInvMassPt(NULL),
132         fHistoTruePi0CaloMergedClusterInvMassPt(NULL),
133         fHistoTrueEtaCaloMergedClusterInvMassPt(NULL),
134         fHistoTruePi0CaloMergedClusterPartConvInvMassPt(NULL),
135         fHistoTrueEtaCaloMergedClusterPartConvInvMassPt(NULL),
136         fHistoTruePrimaryPi0InvMassPt(NULL),
137         fHistoTruePrimaryEtaInvMassPt(NULL),
138         fHistoTruePrimaryPi0W0WeightingInvMassPt(NULL),
139         fHistoTruePrimaryEtaW0WeightingInvMassPt(NULL),
140         fProfileTruePrimaryPi0WeightsInvMassPt(NULL),
141         fProfileTruePrimaryEtaWeightsInvMassPt(NULL),
142         fHistoTruePrimaryPi0MCPtResolPt(NULL),
143         fHistoTruePrimaryEtaMCPtResolPt(NULL),
144         fHistoTrueSecondaryPi0InvMassPt(NULL),
145         fHistoTrueSecondaryEtaInvMassPt(NULL),
146         fHistoTrueSecondaryPi0FromK0sInvMassPt(NULL),
147         fHistoTrueK0sWithPi0DaughterMCPt(NULL),
148         fHistoTrueSecondaryPi0FromEtaInvMassPt(NULL),
149         fHistoTrueEtaWithPi0DaughterMCPt(NULL),
150         fHistoTrueSecondaryPi0FromLambdaInvMassPt(NULL),
151         fHistoTrueLambdaWithPi0DaughterMCPt(NULL),
152         fHistoTrueBckGGInvMassPt(NULL),
153         fHistoTrueBckContInvMassPt(NULL),
154         fHistoTruePi0PtY(NULL),
155         fHistoTrueEtaPtY(NULL),
156         fHistoTruePi0PtAlpha(NULL),
157         fHistoTrueEtaPtAlpha(NULL),
158         fHistoTruePi0PtOpenAngle(NULL),
159         fHistoTrueEtaPtOpenAngle(NULL),
160         fHistoClusPhotonBGPt(NULL),
161         fHistoClusPhotonPlusConvBGPt(NULL),
162         fHistoTrueClusGammaPt(NULL),
163         fHistoTrueClusUnConvGammaPt(NULL),
164         fHistoTrueClusUnConvGammaMCPt(NULL),
165         fHistoTrueClusElectronPt(NULL),
166         fHistoTrueClusConvGammaPt(NULL),
167         fHistoTrueClusConvGammaMCPt(NULL),
168         fHistoTrueClusConvGammaFullyPt(NULL),
169         fHistoTrueClusMergedGammaPt(NULL),
170         fHistoTrueClusMergedPartConvGammaPt(NULL),
171         fHistoTrueClusDalitzPt(NULL),
172         fHistoTrueClusDalitzMergedPt(NULL),
173         fHistoTrueClusPhotonFromElecMotherPt(NULL),
174         fHistoTrueClusShowerPt(NULL),
175         fHistoTrueClusSubLeadingPt(NULL),
176         fHistoTrueClusNParticles(NULL),
177         fHistoTrueClusEMNonLeadingPt(NULL),
178         fHistoTrueNLabelsInClus(NULL),
179         fHistoTruePrimaryClusGammaPt(NULL),
180         fHistoTruePrimaryClusGammaESDPtMCPt(NULL),
181         fHistoTruePrimaryClusConvGammaPt(NULL),
182         fHistoTruePrimaryClusConvGammaESDPtMCPt(NULL),
183         fHistoTrueSecondaryClusGammaPt(NULL),
184         fHistoTrueSecondaryClusConvGammaPt(NULL),
185         fHistoTrueSecondaryClusGammaFromXFromK0sPt(NULL),
186         fHistoTrueSecondaryClusConvGammaFromXFromK0sPt(NULL),
187         fHistoTrueSecondaryClusGammaFromXFromLambdaPt(NULL),
188         fHistoTrueSecondaryClusConvGammaFromXFromLambdaPt(NULL),
189         fHistoTrueSecondaryClusGammaFromXFromEtasPt(NULL),
190         fHistoTrueSecondaryClusConvGammaFromXFromEtasPt(NULL),
191         fHistoNEvents(NULL),
192         fHistoNGoodESDTracks(NULL),
193         fHistoNGammaCandidates(NULL),
194         fHistoNGoodESDTracksVsNGammaCanditates(NULL),
195         fHistoNV0Tracks(NULL),
196         fProfileEtaShift(NULL),
197         fEventPlaneAngle(-100),
198         fRandom(0),
199         fnCuts(0),
200         fiCut(0),
201         fIsHeavyIon(0),
202         fDoMesonAnalysis(kTRUE),
203         fDoMesonQA(0),
204         fDoClusterQA(0),
205         fIsFromMBHeader(kTRUE),
206         fIsMC(kFALSE)
207 {
208   
209 }
210
211 //________________________________________________________________________
212 AliAnalysisTaskGammaCalo::AliAnalysisTaskGammaCalo(const char *name):
213         AliAnalysisTaskSE(name),
214         fV0Reader(NULL),
215         fBGHandler(NULL),
216         fInputEvent(NULL),
217         fMCEvent(NULL),
218         fMCStack(NULL),
219         fCutFolder(NULL),
220         fESDList(NULL),
221         fBackList(NULL),
222         fMotherList(NULL),
223         fTrueList(NULL),
224         fMCList(NULL),
225         fHeaderNameList(NULL),
226         fOutputContainer(0),
227         fClusterCandidates(NULL),
228         fEventCutArray(NULL),
229         fEventCuts(NULL),
230         fClusterCutArray(NULL),
231         fCaloPhotonCuts(NULL),
232         fMesonCutArray(NULL),
233         fMesonCuts(NULL),
234         fHistoMotherInvMassPt(NULL),
235         fHistoMotherInvMass3ClusterPt(NULL),
236         fSparseMotherInvMassPtZM(NULL),
237         fHistoMotherBackInvMassPt(NULL),
238         fSparseMotherBackInvMassPtZM(NULL),
239         fHistoMotherInvMassEalpha(NULL),
240         fHistoMotherPi0PtY(NULL),
241         fHistoMotherEtaPtY(NULL),
242         fHistoMotherPi0PtAlpha(NULL),
243         fHistoMotherEtaPtAlpha(NULL),
244         fHistoMotherPi0PtOpenAngle(NULL),
245         fHistoMotherEtaPtOpenAngle(NULL),
246         fHistoMotherInvMassECalib(NULL),
247         fHistoMotherInvMassECalibalpha(NULL),
248         fHistoClusGammaPt(NULL),
249         fHistoMCHeaders(NULL),
250         fHistoMCAllGammaPt(NULL),
251         fHistoMCDecayGammaPi0Pt(NULL),
252         fHistoMCDecayGammaRhoPt(NULL),
253         fHistoMCDecayGammaEtaPt(NULL),
254         fHistoMCDecayGammaOmegaPt(NULL),
255         fHistoMCDecayGammaEtapPt(NULL),
256         fHistoMCDecayGammaPhiPt(NULL),
257         fHistoMCDecayGammaSigmaPt(NULL),
258         fHistoMCPi0Pt(NULL),
259         fHistoMCPi0WOWeightPt(NULL),
260         fHistoMCEtaPt(NULL),
261         fHistoMCEtaWOWeightPt(NULL),
262         fHistoMCPi0InAccPt(NULL),
263         fHistoMCEtaInAccPt(NULL),
264         fHistoMCPi0PtY(NULL),
265         fHistoMCEtaPtY(NULL),
266         fHistoMCK0sPt(NULL),
267         fHistoMCK0sWOWeightPt(NULL),
268         fHistoMCK0sPtY(NULL),
269         fHistoMCSecPi0PtvsSource(NULL),
270         fHistoMCSecPi0Source(NULL),
271         fHistoMCSecEtaPt(NULL),
272         fHistoMCSecEtaSource(NULL),
273         fHistoTruePi0InvMassPt(NULL),
274         fHistoTrueEtaInvMassPt(NULL),
275         fHistoTruePi0CaloPhotonInvMassPt(NULL),
276         fHistoTrueEtaCaloPhotonInvMassPt(NULL),
277         fHistoTruePi0CaloConvertedPhotonInvMassPt(NULL),
278         fHistoTrueEtaCaloConvertedPhotonInvMassPt(NULL),
279         fHistoTruePi0CaloMixedPhotonConvPhotonInvMassPt(NULL),
280         fHistoTrueEtaCaloMixedPhotonConvPhotonInvMassPt(NULL),
281         fHistoTruePi0CaloElectronInvMassPt(NULL),
282         fHistoTrueEtaCaloElectronInvMassPt(NULL),
283         fHistoTruePi0CaloMergedClusterInvMassPt(NULL),
284         fHistoTrueEtaCaloMergedClusterInvMassPt(NULL),
285         fHistoTruePi0CaloMergedClusterPartConvInvMassPt(NULL),
286         fHistoTrueEtaCaloMergedClusterPartConvInvMassPt(NULL),
287         fHistoTruePrimaryPi0InvMassPt(NULL),
288         fHistoTruePrimaryEtaInvMassPt(NULL),
289         fHistoTruePrimaryPi0W0WeightingInvMassPt(NULL),
290         fHistoTruePrimaryEtaW0WeightingInvMassPt(NULL),
291         fProfileTruePrimaryPi0WeightsInvMassPt(NULL),
292         fProfileTruePrimaryEtaWeightsInvMassPt(NULL),
293         fHistoTruePrimaryPi0MCPtResolPt(NULL),
294         fHistoTruePrimaryEtaMCPtResolPt(NULL),
295         fHistoTrueSecondaryPi0InvMassPt(NULL),
296         fHistoTrueSecondaryEtaInvMassPt(NULL),
297         fHistoTrueSecondaryPi0FromK0sInvMassPt(NULL),
298         fHistoTrueK0sWithPi0DaughterMCPt(NULL),
299         fHistoTrueSecondaryPi0FromEtaInvMassPt(NULL),
300         fHistoTrueEtaWithPi0DaughterMCPt(NULL),
301         fHistoTrueSecondaryPi0FromLambdaInvMassPt(NULL),
302         fHistoTrueLambdaWithPi0DaughterMCPt(NULL),
303         fHistoTrueBckGGInvMassPt(NULL),
304         fHistoTrueBckContInvMassPt(NULL),
305         fHistoTruePi0PtY(NULL),
306         fHistoTrueEtaPtY(NULL),
307         fHistoTruePi0PtAlpha(NULL),
308         fHistoTrueEtaPtAlpha(NULL),
309         fHistoTruePi0PtOpenAngle(NULL),
310         fHistoTrueEtaPtOpenAngle(NULL),
311         fHistoClusPhotonBGPt(NULL),
312         fHistoClusPhotonPlusConvBGPt(NULL),
313         fHistoTrueClusGammaPt(NULL),
314         fHistoTrueClusUnConvGammaPt(NULL),
315         fHistoTrueClusUnConvGammaMCPt(NULL),
316         fHistoTrueClusElectronPt(NULL),
317         fHistoTrueClusConvGammaPt(NULL),
318         fHistoTrueClusConvGammaMCPt(NULL),
319         fHistoTrueClusConvGammaFullyPt(NULL),
320         fHistoTrueClusMergedGammaPt(NULL),
321         fHistoTrueClusMergedPartConvGammaPt(NULL),
322         fHistoTrueClusDalitzPt(NULL),
323         fHistoTrueClusDalitzMergedPt(NULL),
324         fHistoTrueClusPhotonFromElecMotherPt(NULL),
325         fHistoTrueClusShowerPt(NULL),
326         fHistoTrueClusSubLeadingPt(NULL),
327         fHistoTrueClusNParticles(NULL),
328         fHistoTrueClusEMNonLeadingPt(NULL),
329         fHistoTrueNLabelsInClus(NULL),
330         fHistoTruePrimaryClusGammaPt(NULL),
331         fHistoTruePrimaryClusGammaESDPtMCPt(NULL),
332         fHistoTruePrimaryClusConvGammaPt(NULL),
333         fHistoTruePrimaryClusConvGammaESDPtMCPt(NULL),
334         fHistoTrueSecondaryClusGammaPt(NULL),
335         fHistoTrueSecondaryClusConvGammaPt(NULL),
336         fHistoTrueSecondaryClusGammaFromXFromK0sPt(NULL),
337         fHistoTrueSecondaryClusConvGammaFromXFromK0sPt(NULL),
338         fHistoTrueSecondaryClusGammaFromXFromLambdaPt(NULL),
339         fHistoTrueSecondaryClusConvGammaFromXFromLambdaPt(NULL),
340         fHistoTrueSecondaryClusGammaFromXFromEtasPt(NULL),
341         fHistoTrueSecondaryClusConvGammaFromXFromEtasPt(NULL),
342         fHistoNEvents(NULL),
343         fHistoNGoodESDTracks(NULL),
344         fHistoNGammaCandidates(NULL),
345         fHistoNGoodESDTracksVsNGammaCanditates(NULL),
346         fHistoNV0Tracks(NULL),
347         fProfileEtaShift(NULL),
348         fEventPlaneAngle(-100),
349         fRandom(0),
350         fnCuts(0),
351         fiCut(0),
352         fIsHeavyIon(0),
353         fDoMesonAnalysis(kTRUE),
354         fDoMesonQA(0),
355         fDoClusterQA(0),
356         fIsFromMBHeader(kTRUE),
357         fIsMC(kFALSE)
358 {
359   // Define output slots here
360   DefineOutput(1, TList::Class());
361 }
362
363 AliAnalysisTaskGammaCalo::~AliAnalysisTaskGammaCalo()
364 {
365         if(fClusterCandidates){
366                 delete fClusterCandidates;
367                 fClusterCandidates = 0x0;
368         }
369         if(fBGHandler){
370                 delete[] fBGHandler;
371                 fBGHandler = 0x0;
372         }
373 }
374 //___________________________________________________________
375 void AliAnalysisTaskGammaCalo::InitBack(){
376         
377         const Int_t nDim = 4;
378         Int_t nBins[nDim] = {800,250,7,4};
379         Double_t xMin[nDim] = {0,0, 0,0};
380         Double_t xMax[nDim] = {0.8,25,7,4};
381         
382         fSparseMotherInvMassPtZM = new THnSparseF*[fnCuts];
383         fSparseMotherBackInvMassPtZM = new THnSparseF*[fnCuts];
384         
385         fBGHandler = new AliGammaConversionAODBGHandler*[fnCuts];
386
387         
388         for(Int_t iCut = 0; iCut<fnCuts;iCut++){
389                 if (((AliConversionMesonCuts*)fMesonCutArray->At(iCut))->DoBGCalculation()){
390                         TString cutstringEvent  = ((AliConvEventCuts*)fEventCutArray->At(iCut))->GetCutNumber();
391                         TString cutstringCalo   = ((AliCaloPhotonCuts*)fClusterCutArray->At(iCut))->GetCutNumber();
392                         TString cutstringMeson  = ((AliConversionMesonCuts*)fMesonCutArray->At(iCut))->GetCutNumber();
393                         
394                         Int_t collisionSystem = atoi((TString)(((AliConvEventCuts*)fEventCutArray->At(iCut))->GetCutNumber())(0,1));
395                         Int_t centMin = atoi((TString)(((AliConvEventCuts*)fEventCutArray->At(iCut))->GetCutNumber())(1,1));
396                         Int_t centMax = atoi((TString)(((AliConvEventCuts*)fEventCutArray->At(iCut))->GetCutNumber())(2,1));
397                         
398                         if(collisionSystem == 1 || collisionSystem == 2 ||
399                                 collisionSystem == 5 || collisionSystem == 8 ||
400                                 collisionSystem == 9){
401                                 centMin = centMin*10;
402                                 centMax = centMax*10;
403                                 if(centMax ==0 && centMax!=centMin) centMax=100;
404                         } else if(collisionSystem == 3 || collisionSystem == 6){
405                                 centMin = centMin*5;
406                                 centMax = centMax*5;
407                         } else if(collisionSystem == 4 || collisionSystem == 7){
408                                 centMin = ((centMin*5)+45);
409                                 centMax = ((centMax*5)+45);
410                         }
411                         
412                         fBackList[iCut] = new TList();
413                         fBackList[iCut]->SetName(Form("%s_%s_%s Back histograms",cutstringEvent.Data(),cutstringCalo.Data(), cutstringMeson.Data()));
414                         fBackList[iCut]->SetOwner(kTRUE);
415                         fCutFolder[iCut]->Add(fBackList[iCut]);
416                         
417                         fSparseMotherBackInvMassPtZM[iCut] = new THnSparseF("Back_Back_InvMass_Pt_z_m","Back_Back_InvMass_Pt_z_m",nDim,nBins,xMin,xMax);
418                         fBackList[iCut]->Add(fSparseMotherBackInvMassPtZM[iCut]);
419                         
420                         fMotherList[iCut] = new TList();
421                         fMotherList[iCut]->SetName(Form("%s_%s_%s Mother histograms",cutstringEvent.Data(),cutstringCalo.Data(),cutstringMeson.Data()));
422                         fMotherList[iCut]->SetOwner(kTRUE);
423                         fCutFolder[iCut]->Add(fMotherList[iCut]);
424                         
425                         fSparseMotherInvMassPtZM[iCut] = new THnSparseF("Back_Mother_InvMass_Pt_z_m","Back_Mother_InvMass_Pt_z_m",nDim,nBins,xMin,xMax);
426                         fMotherList[iCut]->Add(fSparseMotherInvMassPtZM[iCut]);
427                         
428                         if(((AliConversionMesonCuts*)fMesonCutArray->At(iCut))->BackgroundHandlerType() == 0){
429                                 fBGHandler[iCut] = new AliGammaConversionAODBGHandler(
430                                                                                                                                         collisionSystem,centMin,centMax,
431                                                                                                                                         ((AliConversionMesonCuts*)fMesonCutArray->At(iCut))->GetNumberOfBGEvents(),
432                                                                                                                                         ((AliConversionMesonCuts*)fMesonCutArray->At(iCut))->UseTrackMultiplicity());
433                         }
434                 }
435         }
436 }
437 //________________________________________________________________________
438 void AliAnalysisTaskGammaCalo::UserCreateOutputObjects(){
439   
440         // Create histograms
441         if(fOutputContainer != NULL){
442                 delete fOutputContainer;
443                 fOutputContainer = NULL;
444         }
445         if(fOutputContainer == NULL){
446                 fOutputContainer = new TList();
447                 fOutputContainer->SetOwner(kTRUE);
448         }
449   
450         // Array of current cut's gammas
451         fClusterCandidates = new TList();
452   
453         fCutFolder = new TList*[fnCuts];
454         fESDList = new TList*[fnCuts];
455         fBackList = new TList*[fnCuts];
456         fMotherList = new TList*[fnCuts];
457         fHistoNEvents = new TH1I*[fnCuts];
458         fHistoNGoodESDTracks = new TH1I*[fnCuts];
459         fHistoNGammaCandidates = new TH1I*[fnCuts];
460         fHistoNGoodESDTracksVsNGammaCanditates = new TH2F*[fnCuts];
461         fHistoNV0Tracks = new TH1I*[fnCuts];
462         fProfileEtaShift = new TProfile*[fnCuts];
463         
464         if(fDoMesonAnalysis){
465                 fHistoMotherInvMassPt = new TH2F*[fnCuts];
466                 fHistoMotherInvMass3ClusterPt = new TH2F*[fnCuts];
467                 fHistoMotherBackInvMassPt = new TH2F*[fnCuts];
468                 fHistoMotherInvMassEalpha = new TH2F*[fnCuts];
469                 if (fDoMesonQA > 0){
470                         fHistoMotherPi0PtY =  new TH2F*[fnCuts];
471                         fHistoMotherEtaPtY =  new TH2F*[fnCuts];
472                         fHistoMotherPi0PtAlpha =  new TH2F*[fnCuts];
473                         fHistoMotherEtaPtAlpha =  new TH2F*[fnCuts];
474                         fHistoMotherPi0PtOpenAngle =  new TH2F*[fnCuts];
475                         fHistoMotherEtaPtOpenAngle =  new TH2F*[fnCuts];
476                 }
477                 if(fDoMesonQA == 1){
478                         fHistoMotherInvMassECalib = new TH2F*[fnCuts];
479                         fHistoMotherInvMassECalibalpha = new TH2F*[fnCuts];
480                 }
481         }
482                 
483         fHistoClusGammaPt = new TH1F*[fnCuts];
484         
485         for(Int_t iCut = 0; iCut<fnCuts;iCut++){
486                 TString cutstringEvent  = ((AliConvEventCuts*)fEventCutArray->At(iCut))->GetCutNumber();
487                 TString cutstringCalo   = ((AliCaloPhotonCuts*)fClusterCutArray->At(iCut))->GetCutNumber();
488                 TString cutstringMeson  = "NoMesonCut";
489                 if(fDoMesonAnalysis)cutstringMeson = ((AliConversionMesonCuts*)fMesonCutArray->At(iCut))->GetCutNumber();
490     
491                 fCutFolder[iCut] = new TList();
492                 fCutFolder[iCut]->SetName(Form("Cut Number %s_%s_%s",cutstringEvent.Data(),cutstringCalo.Data(),cutstringMeson.Data()));
493                 fCutFolder[iCut]->SetOwner(kTRUE);
494                 fOutputContainer->Add(fCutFolder[iCut]);
495                 fESDList[iCut] = new TList();
496                 fESDList[iCut]->SetName(Form("%s_%s_%s ESD histograms",cutstringEvent.Data(),cutstringCalo.Data(),cutstringMeson.Data()));
497                 fESDList[iCut]->SetOwner(kTRUE);
498                 fCutFolder[iCut]->Add(fESDList[iCut]);
499     
500                 fHistoNEvents[iCut] = new TH1I("NEvents","NEvents",9,-0.5,8.5);
501                 fHistoNEvents[iCut]->GetXaxis()->SetBinLabel(1,"Accepted");
502                 fHistoNEvents[iCut]->GetXaxis()->SetBinLabel(2,"Centrality");
503                 fHistoNEvents[iCut]->GetXaxis()->SetBinLabel(3,"Missing MC");
504                 if (((AliConvEventCuts*)fEventCutArray->At(iCut))->IsSpecialTrigger() == 4 ){
505                         TString TriggerNames = "Not Trigger: ";
506                         TriggerNames = TriggerNames+ ( (AliConvEventCuts*)fEventCutArray->At(iCut))->GetSpecialTriggerName();
507                         fHistoNEvents[iCut]->GetXaxis()->SetBinLabel(4,TriggerNames.Data());
508                 } else {
509                         fHistoNEvents[iCut]->GetXaxis()->SetBinLabel(4,"Trigger");
510                 }
511                 fHistoNEvents[iCut]->GetXaxis()->SetBinLabel(5,"Vertex Z");
512                 fHistoNEvents[iCut]->GetXaxis()->SetBinLabel(6,"Cont. Vertex");
513                 fHistoNEvents[iCut]->GetXaxis()->SetBinLabel(7,"Pile-Up");
514                 fHistoNEvents[iCut]->GetXaxis()->SetBinLabel(8,"no SDD");
515                 fHistoNEvents[iCut]->GetXaxis()->SetBinLabel(9,"no V0AND");
516                 fESDList[iCut]->Add(fHistoNEvents[iCut]);
517                 
518                 if(fIsHeavyIon == 1) fHistoNGoodESDTracks[iCut] = new TH1I("GoodESDTracks","GoodESDTracks",4000,0,4000);
519                 else if(fIsHeavyIon == 2) fHistoNGoodESDTracks[iCut] = new TH1I("GoodESDTracks","GoodESDTracks",400,0,400);
520                 else fHistoNGoodESDTracks[iCut] = new TH1I("GoodESDTracks","GoodESDTracks",200,0,200);
521                 fESDList[iCut]->Add(fHistoNGoodESDTracks[iCut]);
522                 if(fIsHeavyIon == 1) fHistoNGammaCandidates[iCut] = new TH1I("GammaCandidates","GammaCandidates",600,0,600);
523                 else if(fIsHeavyIon == 2) fHistoNGammaCandidates[iCut] = new TH1I("GammaCandidates","GammaCandidates",400,0,400);
524                 else fHistoNGammaCandidates[iCut] = new TH1I("GammaCandidates","GammaCandidates",100,0,100);
525                 fESDList[iCut]->Add(fHistoNGammaCandidates[iCut]);
526                 if(fIsHeavyIon == 1) fHistoNGoodESDTracksVsNGammaCanditates[iCut] = new TH2F("GoodESDTracksVsGammaCandidates","GoodESDTracksVsGammaCandidates",4000,0,4000,600,0,600);
527                 else if(fIsHeavyIon == 2) fHistoNGoodESDTracksVsNGammaCanditates[iCut] = new TH2F("GoodESDTracksVsGammaCandidates","GoodESDTracksVsGammaCandidates",400,0,400,400,0,400);
528                 else fHistoNGoodESDTracksVsNGammaCanditates[iCut] = new TH2F("GoodESDTracksVsGammaCandidates","GoodESDTracksVsGammaCandidates",200,0,200,100,0,100);
529                 fESDList[iCut]->Add(fHistoNGoodESDTracksVsNGammaCanditates[iCut]);
530     
531                 
532                 if(fIsHeavyIon == 1) fHistoNV0Tracks[iCut] = new TH1I("V0 Multiplicity","V0 Multiplicity",30000,0,30000);
533                 else if(fIsHeavyIon == 2) fHistoNV0Tracks[iCut] = new TH1I("V0 Multiplicity","V0 Multiplicity",2500,0,2500);
534                 else fHistoNV0Tracks[iCut] = new TH1I("V0 Multiplicity","V0 Multiplicity",1500,0,1500);
535                 fESDList[iCut]->Add(fHistoNV0Tracks[iCut]);
536                 fProfileEtaShift[iCut] = new TProfile("Eta Shift","Eta Shift",1, -0.5,0.5);
537                 fESDList[iCut]->Add(fProfileEtaShift[iCut]);
538     
539                 fHistoClusGammaPt[iCut] = new TH1F("ClusGamma_Pt","ClusGamma_Pt",250,0,25);
540                 fESDList[iCut]->Add(fHistoClusGammaPt[iCut]);
541
542                 
543                 if(fDoMesonAnalysis){
544                         fHistoMotherInvMassPt[iCut] = new TH2F("ESD_Mother_InvMass_Pt","ESD_Mother_InvMass_Pt",800,0,0.8,250,0,25);
545                         fESDList[iCut]->Add(fHistoMotherInvMassPt[iCut]);
546                         fHistoMotherInvMass3ClusterPt[iCut] = new TH2F("ESD_Mother_InvMass3Cluster_Pt","ESD_Mother_InvMass3Cluster_Pt",800,0,0.8,250,0,25);
547                         fESDList[iCut]->Add(fHistoMotherInvMass3ClusterPt[iCut]);
548                         fHistoMotherBackInvMassPt[iCut] = new TH2F("ESD_Background_InvMass_Pt","ESD_Background_InvMass_Pt",800,0,0.8,250,0,25);
549                         fESDList[iCut]->Add(fHistoMotherBackInvMassPt[iCut]);
550                         fHistoMotherInvMassEalpha[iCut] = new TH2F("ESD_Mother_InvMass_vs_E_alpha","ESD_Mother_InvMass_vs_E_alpha",800,0,0.8,250,0,25);
551                         fESDList[iCut]->Add(fHistoMotherInvMassEalpha[iCut]);
552                         if(fDoMesonQA == 1){
553                                 fHistoMotherInvMassECalib[iCut] = new TH2F("ESD_Mother_InvMass_Pt_Calib","ESD_Mother_InvMass_Pt_Calib",800,0,0.8,250,0,25);
554                                 fESDList[iCut]->Add(fHistoMotherInvMassECalib[iCut]);
555                                 fHistoMotherInvMassECalibalpha[iCut] = new TH2F("ESD_Mother_InvMass_vs_E_Calib_alpha","ESD_Mother_InvMass_vs_E_Calib_alpha",800,0,0.8,250,0,25);
556                                 fESDList[iCut]->Add(fHistoMotherInvMassECalibalpha[iCut]);
557                         }
558
559                         if (fDoMesonQA > 0 ){
560                                 fHistoMotherPi0PtY[iCut] = new TH2F("ESD_MotherPi0_Pt_Y","ESD_MotherPi0_Pt_Y",150,0.03,15.,150,-1.5,1.5);
561                                 SetLogBinningXTH2(fHistoMotherPi0PtY[iCut]);
562                                 fESDList[iCut]->Add(fHistoMotherPi0PtY[iCut]);
563                                 fHistoMotherEtaPtY[iCut] = new TH2F("ESD_MotherEta_Pt_Y","ESD_MotherEta_Pt_Y",150,0.03,15.,150,-1.5,1.5);
564                                 SetLogBinningXTH2(fHistoMotherEtaPtY[iCut]);
565                                 fESDList[iCut]->Add(fHistoMotherEtaPtY[iCut]);
566                                 fHistoMotherPi0PtAlpha[iCut] = new TH2F("ESD_MotherPi0_Pt_Alpha","ESD_MotherPi0_Pt_Alpha",150,0.03,15.,100,0,1);
567                                 SetLogBinningXTH2(fHistoMotherPi0PtAlpha[iCut]);
568                                 fESDList[iCut]->Add(fHistoMotherPi0PtAlpha[iCut]);
569                                 fHistoMotherEtaPtAlpha[iCut] = new TH2F("ESD_MotherEta_Pt_Alpha","ESD_MotherEta_Pt_Alpha",150,0.03,15.,100,0,1);
570                                 SetLogBinningXTH2(fHistoMotherEtaPtAlpha[iCut]);
571                                 fESDList[iCut]->Add(fHistoMotherEtaPtAlpha[iCut]);
572                                 fHistoMotherPi0PtOpenAngle[iCut] = new TH2F("ESD_MotherPi0_Pt_OpenAngle","ESD_MotherPi0_Pt_OpenAngle",150,0.03,15.,200,0,2*TMath::Pi());
573                                 SetLogBinningXTH2(fHistoMotherPi0PtOpenAngle[iCut]);
574                                 fESDList[iCut]->Add(fHistoMotherPi0PtOpenAngle[iCut]);
575                                 fHistoMotherEtaPtOpenAngle[iCut] = new TH2F("ESD_MotherEta_Pt_OpenAngle","ESD_MotherEta_Pt_OpenAngle",150,0.03,15.,200,0,2*TMath::Pi());
576                                 SetLogBinningXTH2(fHistoMotherEtaPtOpenAngle[iCut]);
577                                 fESDList[iCut]->Add(fHistoMotherEtaPtOpenAngle[iCut]);
578                         }
579                 }    
580         }
581         if(fDoMesonAnalysis){
582                 InitBack(); // Init Background Handler
583         }
584   
585         if(fIsMC){
586                 // MC Histogramms
587                 fMCList         = new TList*[fnCuts];
588                 // True Histogramms
589                 fTrueList       = new TList*[fnCuts];
590                 // Selected Header List
591                 fHeaderNameList                                         = new TList*[fnCuts];
592                 fHistoMCHeaders                                         = new TH1I*[fnCuts];
593                 fHistoMCAllGammaPt                                      = new TH1F*[fnCuts];
594                 fHistoMCDecayGammaPi0Pt                         = new TH1F*[fnCuts];
595                 fHistoMCDecayGammaRhoPt                         = new TH1F*[fnCuts];
596                 fHistoMCDecayGammaEtaPt                         = new TH1F*[fnCuts];
597                 fHistoMCDecayGammaOmegaPt                       = new TH1F*[fnCuts];
598                 fHistoMCDecayGammaEtapPt                        = new TH1F*[fnCuts];
599                 fHistoMCDecayGammaPhiPt                         = new TH1F*[fnCuts];
600                 fHistoMCDecayGammaSigmaPt                       = new TH1F*[fnCuts];
601         
602                 fHistoClusPhotonBGPt                                                    = new TH2F*[fnCuts];
603                 fHistoClusPhotonPlusConvBGPt                                    = new TH2F*[fnCuts];
604         
605                 fHistoTrueClusGammaPt                                                           = new TH1F*[fnCuts];
606                 fHistoTruePrimaryClusGammaPt                                            = new TH1F*[fnCuts];
607                 fHistoTruePrimaryClusGammaESDPtMCPt                             = new TH2F*[fnCuts];
608                 fHistoTruePrimaryClusConvGammaPt                                        = new TH1F*[fnCuts];
609                 fHistoTruePrimaryClusConvGammaESDPtMCPt                         = new TH2F*[fnCuts];
610                 fHistoTrueSecondaryClusGammaPt                                          = new TH1F*[fnCuts];
611                 fHistoTrueSecondaryClusConvGammaPt                                      = new TH1F*[fnCuts];
612                 fHistoTrueSecondaryClusGammaFromXFromK0sPt                      = new TH1F*[fnCuts];
613                 fHistoTrueSecondaryClusConvGammaFromXFromK0sPt          = new TH1F*[fnCuts];
614                 fHistoTrueSecondaryClusGammaFromXFromLambdaPt           = new TH1F*[fnCuts];
615                 fHistoTrueSecondaryClusConvGammaFromXFromLambdaPt       = new TH1F*[fnCuts];
616                 fHistoTrueSecondaryClusGammaFromXFromEtasPt             = new TH1F*[fnCuts];
617                 fHistoTrueSecondaryClusConvGammaFromXFromEtasPt         = new TH1F*[fnCuts];
618     
619                 
620                 if (fDoClusterQA > 0){  
621                         fHistoTrueClusUnConvGammaPt             = new TH1F*[fnCuts];
622                         fHistoTrueClusUnConvGammaMCPt           = new TH1F*[fnCuts];
623                         fHistoTrueClusElectronPt                        = new TH1F*[fnCuts];
624                         fHistoTrueClusConvGammaPt                       = new TH1F*[fnCuts];
625                         fHistoTrueClusConvGammaMCPt                     = new TH1F*[fnCuts];
626                         fHistoTrueClusConvGammaFullyPt          = new TH1F*[fnCuts];
627                         fHistoTrueClusMergedGammaPt             = new TH1F*[fnCuts];
628                         fHistoTrueClusMergedPartConvGammaPt = new TH1F*[fnCuts];
629                         fHistoTrueClusDalitzPt                          = new TH1F*[fnCuts];
630                         fHistoTrueClusDalitzMergedPt            = new TH1F*[fnCuts];
631                         fHistoTrueClusPhotonFromElecMotherPt= new TH1F*[fnCuts];
632                         fHistoTrueClusShowerPt                          = new TH1F*[fnCuts];
633                         fHistoTrueClusSubLeadingPt                      = new TH1F*[fnCuts];
634                         fHistoTrueClusNParticles                        = new TH1I*[fnCuts];
635                         fHistoTrueClusEMNonLeadingPt            = new TH1F*[fnCuts];
636                         fHistoTrueNLabelsInClus                         = new TH1F*[fnCuts];                    
637                 }
638     
639                 if(fDoMesonAnalysis){
640                         fHistoMCPi0Pt                                   = new TH1F*[fnCuts];
641                         fHistoMCPi0WOWeightPt                   = new TH1F*[fnCuts];
642                         fHistoMCEtaPt                                   = new TH1F*[fnCuts];
643                         fHistoMCEtaWOWeightPt                   = new TH1F*[fnCuts];
644                         fHistoMCPi0InAccPt                              = new TH1F*[fnCuts];
645                         fHistoMCEtaInAccPt                              = new TH1F*[fnCuts];
646       
647                         fHistoTruePi0InvMassPt                                  = new TH2F*[fnCuts];
648                         fHistoTrueEtaInvMassPt                                  = new TH2F*[fnCuts];
649                         fHistoTruePrimaryPi0InvMassPt                   = new TH2F*[fnCuts];
650                         fHistoTruePrimaryEtaInvMassPt                   = new TH2F*[fnCuts];
651                         fHistoTruePrimaryPi0W0WeightingInvMassPt = new TH2F*[fnCuts];
652                         fHistoTruePrimaryEtaW0WeightingInvMassPt = new TH2F*[fnCuts];
653                         fProfileTruePrimaryPi0WeightsInvMassPt  = new TProfile2D*[fnCuts];
654                         fProfileTruePrimaryEtaWeightsInvMassPt  = new TProfile2D*[fnCuts];
655                         fHistoTrueSecondaryPi0InvMassPt                         = new TH2F*[fnCuts];
656                         fHistoTrueSecondaryEtaInvMassPt                         = new TH2F*[fnCuts];
657                         fHistoTrueSecondaryPi0FromK0sInvMassPt  = new TH2F*[fnCuts];
658                         fHistoTrueSecondaryPi0FromEtaInvMassPt  = new TH2F*[fnCuts];
659                         fHistoTrueSecondaryPi0FromLambdaInvMassPt = new TH2F*[fnCuts];
660                         if (fDoMesonQA > 0){
661                                 fHistoMCPi0PtY                                                          = new TH2F*[fnCuts];
662                                 fHistoMCEtaPtY                                                          = new TH2F*[fnCuts];
663                                 fHistoMCK0sPt                                                           = new TH1F*[fnCuts];
664                                 fHistoMCK0sWOWeightPt                                           = new TH1F*[fnCuts];
665                                 fHistoMCK0sPtY                                                          = new TH2F*[fnCuts];
666                                 fHistoMCSecPi0PtvsSource                                        = new TH2F*[fnCuts];
667                                 fHistoMCSecPi0Source                                            = new TH1F*[fnCuts];
668                                 fHistoMCSecEtaPt                                                        = new TH1F*[fnCuts];
669                                 fHistoMCSecEtaSource                                            = new TH1F*[fnCuts];
670                                 fHistoTruePi0CaloPhotonInvMassPt                        = new TH2F*[fnCuts];
671                                 fHistoTrueEtaCaloPhotonInvMassPt                        = new TH2F*[fnCuts];
672                                 fHistoTruePi0CaloMixedPhotonConvPhotonInvMassPt= new TH2F*[fnCuts];
673                                 fHistoTrueEtaCaloMixedPhotonConvPhotonInvMassPt= new TH2F*[fnCuts];
674                                 fHistoTruePi0CaloConvertedPhotonInvMassPt       = new TH2F*[fnCuts];
675                                 fHistoTrueEtaCaloConvertedPhotonInvMassPt       = new TH2F*[fnCuts];
676                                 fHistoTruePi0CaloElectronInvMassPt              = new TH2F*[fnCuts];
677                                 fHistoTrueEtaCaloElectronInvMassPt              = new TH2F*[fnCuts];
678                                 fHistoTruePi0CaloMergedClusterInvMassPt         = new TH2F*[fnCuts];
679                                 fHistoTrueEtaCaloMergedClusterInvMassPt         = new TH2F*[fnCuts];
680                                 fHistoTruePi0CaloMergedClusterPartConvInvMassPt         = new TH2F*[fnCuts];
681                                 fHistoTrueEtaCaloMergedClusterPartConvInvMassPt         = new TH2F*[fnCuts];
682                                 fHistoTruePrimaryPi0MCPtResolPt                         = new TH2F*[fnCuts];
683                                 fHistoTruePrimaryEtaMCPtResolPt                         = new TH2F*[fnCuts];
684                                 fHistoTrueK0sWithPi0DaughterMCPt                        = new TH1F*[fnCuts];
685                                 fHistoTrueEtaWithPi0DaughterMCPt                        = new TH1F*[fnCuts];
686                                 fHistoTrueLambdaWithPi0DaughterMCPt             = new TH1F*[fnCuts];
687                                 fHistoTrueBckGGInvMassPt                                        = new TH2F*[fnCuts];
688                                 fHistoTrueBckContInvMassPt                                      = new TH2F*[fnCuts];
689                                 fHistoTruePi0PtY                                                        = new TH2F*[fnCuts];
690                                 fHistoTrueEtaPtY                                                        = new TH2F*[fnCuts];
691                                 fHistoTruePi0PtAlpha                                            = new TH2F*[fnCuts];
692                                 fHistoTrueEtaPtAlpha                                            = new TH2F*[fnCuts];
693                                 fHistoTruePi0PtOpenAngle                                        = new TH2F*[fnCuts];
694                                 fHistoTrueEtaPtOpenAngle                                        = new TH2F*[fnCuts];
695                         }
696                 }
697     
698     
699     
700                 for(Int_t iCut = 0; iCut<fnCuts;iCut++){
701                         TString cutstringEvent  = ((AliConvEventCuts*)fEventCutArray->At(iCut))->GetCutNumber();
702                         TString cutstringCalo   = ((AliCaloPhotonCuts*)fClusterCutArray->At(iCut))->GetCutNumber();
703                         TString cutstringMeson  = "NoMesonCut";
704                         if(fDoMesonAnalysis)cutstringMeson = ((AliConversionMesonCuts*)fMesonCutArray->At(iCut))->GetCutNumber();
705       
706                         fMCList[iCut] = new TList();
707                         fMCList[iCut]->SetName(Form("%s_%s_%s MC histograms",cutstringEvent.Data(),cutstringCalo.Data(),cutstringMeson.Data()));
708                         fMCList[iCut]->SetOwner(kTRUE);
709                         fCutFolder[iCut]->Add(fMCList[iCut]);
710                         fHistoMCHeaders[iCut] = new TH1I("MC_Headers","MC_Headers",20,0,20);
711                         fMCList[iCut]->Add(fHistoMCHeaders[iCut]);
712                         fHistoMCAllGammaPt[iCut] = new TH1F("MC_AllGamma_Pt","MC_AllGamma_Pt",250,0,25);
713                         fMCList[iCut]->Add(fHistoMCAllGammaPt[iCut]);
714                         fHistoMCDecayGammaPi0Pt[iCut] = new TH1F("MC_DecayGammaPi0_Pt","MC_DecayGammaPi0_Pt",250,0,25);
715                         fMCList[iCut]->Add(fHistoMCDecayGammaPi0Pt[iCut]);
716                         fHistoMCDecayGammaRhoPt[iCut] = new TH1F("MC_DecayGammaRho_Pt","MC_DecayGammaRho_Pt",250,0,25);
717                         fMCList[iCut]->Add(fHistoMCDecayGammaRhoPt[iCut]);
718                         fHistoMCDecayGammaEtaPt[iCut] = new TH1F("MC_DecayGammaEta_Pt","MC_DecayGammaEta_Pt",250,0,25);
719                         fMCList[iCut]->Add(fHistoMCDecayGammaEtaPt[iCut]);
720                         fHistoMCDecayGammaOmegaPt[iCut] = new TH1F("MC_DecayGammaOmega_Pt","MC_DecayGammaOmmega_Pt",250,0,25);
721                         fMCList[iCut]->Add(fHistoMCDecayGammaOmegaPt[iCut]);
722                         fHistoMCDecayGammaEtapPt[iCut] = new TH1F("MC_DecayGammaEtap_Pt","MC_DecayGammaEtap_Pt",250,0,25);
723                         fMCList[iCut]->Add(fHistoMCDecayGammaEtapPt[iCut]);
724                         fHistoMCDecayGammaPhiPt[iCut] = new TH1F("MC_DecayGammaPhi_Pt","MC_DecayGammaPhi_Pt",250,0,25);
725                         fMCList[iCut]->Add(fHistoMCDecayGammaPhiPt[iCut]);
726                         fHistoMCDecayGammaSigmaPt[iCut] = new TH1F("MC_DecayGammaSigma_Pt","MC_DecayGammaSigma_Pt",250,0,25);
727                         fMCList[iCut]->Add(fHistoMCDecayGammaSigmaPt[iCut]);
728                         
729                         if(fDoMesonAnalysis){
730                                 fHistoMCPi0Pt[iCut] = new TH1F("MC_Pi0_Pt","MC_Pi0_Pt",250,0,25);
731                                 fHistoMCPi0Pt[iCut]->Sumw2();
732                                 fMCList[iCut]->Add(fHistoMCPi0Pt[iCut]);
733                                 fHistoMCPi0WOWeightPt[iCut] = new TH1F("MC_Pi0_WOWeights_Pt","MC_Pi0_WOWeights_Pt",250,0,25);
734                                 fHistoMCPi0WOWeightPt[iCut]->Sumw2();
735                                 fMCList[iCut]->Add(fHistoMCPi0WOWeightPt[iCut]);
736                                 
737                                 fHistoMCEtaPt[iCut] = new TH1F("MC_Eta_Pt","MC_Eta_Pt",250,0,25);
738                                 fHistoMCEtaPt[iCut]->Sumw2();
739                                 fMCList[iCut]->Add(fHistoMCEtaPt[iCut]);
740                                 fHistoMCEtaWOWeightPt[iCut] = new TH1F("MC_Eta_WOWeights_Pt","MC_Eta_WOWeights_Pt",250,0,25);
741                                 fHistoMCEtaWOWeightPt[iCut]->Sumw2();
742                                 fMCList[iCut]->Add(fHistoMCEtaWOWeightPt[iCut]);
743                                 fHistoMCPi0InAccPt[iCut] = new TH1F("MC_Pi0InAcc_Pt","MC_Pi0InAcc_Pt",250,0,25);
744                                 fHistoMCPi0InAccPt[iCut]->Sumw2();
745                                 fMCList[iCut]->Add(fHistoMCPi0InAccPt[iCut]);
746                                 fHistoMCEtaInAccPt[iCut] = new TH1F("MC_EtaInAcc_Pt","MC_EtaInAcc_Pt",250,0,25);
747                                 fHistoMCEtaInAccPt[iCut]->Sumw2();
748                                 fMCList[iCut]->Add(fHistoMCEtaInAccPt[iCut]);
749                                 if (fDoMesonQA > 0){
750                                         fHistoMCPi0PtY[iCut] = new TH2F("MC_Pi0_Pt_Y","MC_Pi0_Pt_Y",150,0.03,15.,150,-1.5,1.5);
751                                         fHistoMCPi0PtY[iCut]->Sumw2();
752                                         SetLogBinningXTH2(fHistoMCPi0PtY[iCut]);
753                                         fMCList[iCut]->Add(fHistoMCPi0PtY[iCut]);
754                                         fHistoMCEtaPtY[iCut] = new TH2F("MC_Eta_Pt_Y","MC_Eta_Pt_Y",150,0.03,15.,150,-1.5,1.5);
755                                         fHistoMCEtaPtY[iCut]->Sumw2();
756                                         SetLogBinningXTH2(fHistoMCEtaPtY[iCut]);
757                                         fMCList[iCut]->Add(fHistoMCEtaPtY[iCut]);
758                                         fHistoMCK0sPt[iCut] = new TH1F("MC_K0s_Pt","MC_K0s_Pt",150,0,15);
759                                         fHistoMCK0sPt[iCut]->Sumw2();
760                                         fMCList[iCut]->Add(fHistoMCK0sPt[iCut]);
761                                         fHistoMCK0sWOWeightPt[iCut] = new TH1F("MC_K0s_WOWeights_Pt","MC_K0s_WOWeights_Pt",150,0,15);
762                                         fHistoMCK0sWOWeightPt[iCut]->Sumw2();
763                                         fMCList[iCut]->Add(fHistoMCK0sWOWeightPt[iCut]);
764                                         fHistoMCK0sPtY[iCut] = new TH2F("MC_K0s_Pt_Y","MC_K0s_Pt_Y",150,0.03,15.,150,-1.5,1.5);
765                                         fHistoMCK0sPtY[iCut]->Sumw2();
766                                         SetLogBinningXTH2(fHistoMCK0sPtY[iCut]);
767                                         fMCList[iCut]->Add(fHistoMCK0sPtY[iCut]);
768                                         
769                                         fHistoMCSecPi0Source[iCut] = new TH1F("MC_SecPi0_Source","MC_SecPi0_Source",5000,0.,5000);
770                                         fMCList[iCut]->Add(fHistoMCSecPi0Source[iCut]);
771                                         fHistoMCSecEtaSource[iCut] = new TH1F("MC_SecEta_Source","MC_SecEta_Source",5000,0,5000);
772                                         fMCList[iCut]->Add(fHistoMCSecEtaSource[iCut]);
773                                         fHistoMCSecPi0PtvsSource[iCut] = new TH2F("MC_SecPi0_Pt_Source","MC_SecPi0_Pt_Source",250,0.0,25.,16,-0.5,15.5);
774                                         fHistoMCSecPi0PtvsSource[iCut]->Sumw2();
775                                         fMCList[iCut]->Add(fHistoMCSecPi0PtvsSource[iCut]);
776                                         fHistoMCSecEtaPt[iCut] = new TH1F("MC_SecEta_Pt","MC_SecEta_Pt",250,0,25);
777                                         fHistoMCSecEtaPt[iCut]->Sumw2();
778                                         fMCList[iCut]->Add(fHistoMCSecEtaPt[iCut]);
779                                 }
780         
781                         }
782                         fTrueList[iCut] = new TList();
783                         fTrueList[iCut]->SetName(Form("%s_%s_%s True histograms",cutstringEvent.Data(),cutstringCalo.Data(),cutstringMeson.Data()));
784                         fTrueList[iCut]->SetOwner(kTRUE);
785                         fCutFolder[iCut]->Add(fTrueList[iCut]);
786                   
787                         fHistoClusPhotonBGPt[iCut] = new TH2F("ESD_TrueClusPhotonBG_Pt","ESD_TrueClusPhotonBG_Pt",250,0,25,9,-0.5,8.5);
788                         fHistoClusPhotonBGPt[iCut]->GetYaxis()->SetBinLabel( 1,"Elec");
789                         fHistoClusPhotonBGPt[iCut]->GetYaxis()->SetBinLabel( 2,"Pion");
790                         fHistoClusPhotonBGPt[iCut]->GetYaxis()->SetBinLabel( 3,"Proton");
791                         fHistoClusPhotonBGPt[iCut]->GetYaxis()->SetBinLabel( 4,"Kaon");
792                         fHistoClusPhotonBGPt[iCut]->GetYaxis()->SetBinLabel( 5,"Neutron");
793                         fHistoClusPhotonBGPt[iCut]->GetYaxis()->SetBinLabel( 6,"K0s");
794                         fHistoClusPhotonBGPt[iCut]->GetYaxis()->SetBinLabel( 7,"Lambda");
795                         fHistoClusPhotonBGPt[iCut]->GetYaxis()->SetBinLabel( 8,"Muon");
796                         fHistoClusPhotonBGPt[iCut]->GetYaxis()->SetBinLabel(9,"Rest");
797                         fTrueList[iCut]->Add(fHistoClusPhotonBGPt[iCut]);
798                         fHistoClusPhotonPlusConvBGPt[iCut] = new TH2F("ESD_TrueClusPhotonPlusConvBG_Pt","ESD_TrueClusPhotonPlusConvBG_Pt",250,0,25,9,-0.5,8.5);
799                         fHistoClusPhotonPlusConvBGPt[iCut]->GetYaxis()->SetBinLabel( 1,"Elec");
800                         fHistoClusPhotonPlusConvBGPt[iCut]->GetYaxis()->SetBinLabel( 2,"Pion");
801                         fHistoClusPhotonPlusConvBGPt[iCut]->GetYaxis()->SetBinLabel( 3,"Proton");
802                         fHistoClusPhotonPlusConvBGPt[iCut]->GetYaxis()->SetBinLabel( 4,"Kaon");
803                         fHistoClusPhotonPlusConvBGPt[iCut]->GetYaxis()->SetBinLabel( 5,"Neutron");
804                         fHistoClusPhotonPlusConvBGPt[iCut]->GetYaxis()->SetBinLabel( 6,"K0s");
805                         fHistoClusPhotonPlusConvBGPt[iCut]->GetYaxis()->SetBinLabel( 7,"Lambda");
806                         fHistoClusPhotonPlusConvBGPt[iCut]->GetYaxis()->SetBinLabel( 8,"Muon");
807                         fHistoClusPhotonPlusConvBGPt[iCut]->GetYaxis()->SetBinLabel(9,"Rest");
808                         fTrueList[iCut]->Add(fHistoClusPhotonPlusConvBGPt[iCut]);
809                 
810                         fHistoTrueClusGammaPt[iCut] = new TH1F("TrueClusGamma_Pt","ESD_TrueClusGamma_Pt",250,0,25);
811                         fTrueList[iCut]->Add(fHistoTrueClusGammaPt[iCut]);
812                         fHistoTruePrimaryClusGammaPt[iCut] = new TH1F("TruePrimaryClusGamma_Pt","ESD_TruePrimaryClusGamma_Pt",250,0,25);
813                         fTrueList[iCut]->Add(fHistoTruePrimaryClusGammaPt[iCut]);
814                         fHistoTruePrimaryClusGammaESDPtMCPt[iCut] = new TH2F("TruePrimaryClusGamma_Pt_MCPt","ESD_TruePrimaryClusGamma_MCPt",250,0,25,250,0,25);
815                         fTrueList[iCut]->Add(fHistoTruePrimaryClusGammaESDPtMCPt[iCut]);
816                         fHistoTruePrimaryClusConvGammaPt[iCut] = new TH1F("TruePrimaryClusConvGamma_Pt","ESD_TruePrimaryClusConvGamma_Pt",250,0,25);
817                         fTrueList[iCut]->Add(fHistoTruePrimaryClusConvGammaPt[iCut]);
818                         fHistoTruePrimaryClusConvGammaESDPtMCPt[iCut] = new TH2F("TruePrimaryClusConvGamma_Pt_MCPt","ESD_TruePrimaryClusConvGamma_MCPt",250,0,25,250,0,25);
819                         fTrueList[iCut]->Add(fHistoTruePrimaryClusConvGammaESDPtMCPt[iCut]);
820                         fHistoTrueSecondaryClusGammaPt[iCut] = new TH1F("ESD_TrueSecondaryClusGamma_Pt","ESD_TrueSecondaryClusGamma_Pt",250,0,25);
821                         fTrueList[iCut]->Add(fHistoTrueSecondaryClusGammaPt[iCut]);      
822                         fHistoTrueSecondaryClusConvGammaPt[iCut] = new TH1F("ESD_TrueSecondaryClusConvGamma_Pt","ESD_TrueSecondaryClusConvGamma_Pt",250,0,25);
823                         fTrueList[iCut]->Add(fHistoTrueSecondaryClusConvGammaPt[iCut]);      
824
825                         fHistoTrueSecondaryClusGammaFromXFromK0sPt[iCut] = new TH1F("ESD_TrueSecondaryClusGammaFromXFromK0s_Pt", "ESD_TrueSecondaryClusGammaFromXFromK0s_Pt",250,0,25);
826                         fTrueList[iCut]->Add(fHistoTrueSecondaryClusGammaFromXFromK0sPt[iCut]);
827                         fHistoTrueSecondaryClusConvGammaFromXFromK0sPt[iCut] = new TH1F("ESD_TrueSecondaryClusConvGammaFromXFromK0s_Pt", "ESD_TrueSecondaryClusConvGammaFromXFromK0s_Pt",250,0,25);
828                         fTrueList[iCut]->Add(fHistoTrueSecondaryClusConvGammaFromXFromK0sPt[iCut]);
829                         fHistoTrueSecondaryClusGammaFromXFromLambdaPt[iCut] = new TH1F("ESD_TrueSecondaryClusGammaFromXFromLambda_Pt", "ESD_TrueSecondaryClusGammaFromXFromLambda_Pt",250,0,25);
830                         fTrueList[iCut]->Add(fHistoTrueSecondaryClusGammaFromXFromLambdaPt[iCut]);
831                         fHistoTrueSecondaryClusConvGammaFromXFromLambdaPt[iCut] = new TH1F("ESD_TrueSecondaryClusConvGammaFromXFromLambda_Pt", "ESD_TrueSecondaryClusConvGammaFromXFromLambda_Pt",250,0,25);
832                         fTrueList[iCut]->Add(fHistoTrueSecondaryClusConvGammaFromXFromLambdaPt[iCut]);
833                         fHistoTrueSecondaryClusGammaFromXFromEtasPt[iCut] = new TH1F("ESD_TrueSecondaryClusGammaFromXFromEtas_Pt", "ESD_TrueSecondaryClusGammaFromXFromEtas_Pt",250,0,25);
834                         fTrueList[iCut]->Add(fHistoTrueSecondaryClusGammaFromXFromEtasPt[iCut]);
835                         fHistoTrueSecondaryClusConvGammaFromXFromEtasPt[iCut] = new TH1F("ESD_TrueSecondaryClusConvGammaFromXFromEtas_Pt", "ESD_TrueSecondaryClusConvGammaFromXFromEtas_Pt",250,0,25);
836                         fTrueList[iCut]->Add(fHistoTrueSecondaryClusConvGammaFromXFromEtasPt[iCut]);
837
838                         
839                         if (fDoClusterQA > 0){  
840                                 fHistoTrueClusUnConvGammaPt[iCut] = new TH1F("TrueClusUnConvGamma_Pt","TrueClusUnConvGamma_Pt",250,0,25);
841                                 fTrueList[iCut]->Add(fHistoTrueClusUnConvGammaPt[iCut]);
842                                 fHistoTrueClusUnConvGammaMCPt[iCut] = new TH1F("TrueClusUnConvGamma_MCPt","TrueClusUnConvGamma_MCPt",250,0,25);
843                                 fTrueList[iCut]->Add(fHistoTrueClusUnConvGammaMCPt[iCut]);
844                                 fHistoTrueClusElectronPt[iCut] = new TH1F("TrueClusElectron_Pt","TrueElectronGamma_Pt",250,0,25);
845                                 fTrueList[iCut]->Add(fHistoTrueClusElectronPt[iCut]);
846                                 fHistoTrueClusConvGammaPt[iCut] = new TH1F("TrueClusConvGamma_Pt","TrueClusConvGamma_Pt",250,0,25);
847                                 fTrueList[iCut]->Add(fHistoTrueClusConvGammaPt[iCut]);
848                                 fHistoTrueClusConvGammaMCPt[iCut] = new TH1F("TrueClusConvGamma_MCPt","TrueClusConvGamma_MCPt",250,0,25);
849                                 fTrueList[iCut]->Add(fHistoTrueClusConvGammaMCPt[iCut]);
850                                 fHistoTrueClusConvGammaFullyPt[iCut] = new TH1F("TrueClusConvGammaFullyContained_Pt","TrueClusConvGammaFullyContained_Pt",250,0,25);
851                                 fTrueList[iCut]->Add(fHistoTrueClusConvGammaFullyPt[iCut]);
852                                 fHistoTrueClusMergedGammaPt[iCut] = new TH1F("TrueClusMergedGamma_Pt","TrueClusMergedGamma_Pt",250,0,25);
853                                 fTrueList[iCut]->Add(fHistoTrueClusMergedGammaPt[iCut]);
854                                 fHistoTrueClusMergedPartConvGammaPt[iCut] = new TH1F("TrueClusMergedPartConvGamma_Pt","TrueClusMergedPartConvGamma_Pt",250,0,25);
855                                 fTrueList[iCut]->Add(fHistoTrueClusMergedPartConvGammaPt[iCut]);
856                                 fHistoTrueClusDalitzPt[iCut] = new TH1F("TrueClusDalitz_Pt","TrueClusDalitz_Pt",250,0,25);
857                                 fTrueList[iCut]->Add(fHistoTrueClusDalitzPt[iCut]);
858                                 fHistoTrueClusDalitzMergedPt[iCut] = new TH1F("TrueClusDalitzMerged_Pt","TrueClusDalitzMerged_Pt",250,0,25);
859                                 fTrueList[iCut]->Add(fHistoTrueClusDalitzMergedPt[iCut]);
860                                 fHistoTrueClusPhotonFromElecMotherPt[iCut] = new TH1F("TrueClusPhotonFromElecMother_Pt","TrueClusPhotonFromElecMother_Pt",250,0,25);
861                                 fTrueList[iCut]->Add(fHistoTrueClusPhotonFromElecMotherPt[iCut]);
862                                 fHistoTrueClusShowerPt[iCut] = new TH1F("TrueClusShower_Pt","TrueClusShower_Pt",250,0,25);
863                                 fTrueList[iCut]->Add(fHistoTrueClusShowerPt[iCut]);
864                                 fHistoTrueClusSubLeadingPt[iCut] = new TH1F("TrueClusSubleading_Pt","TrueClusSubleading_Pt",250,0,25);
865                                 fTrueList[iCut]->Add(fHistoTrueClusSubLeadingPt[iCut]);
866                                 fHistoTrueClusNParticles[iCut] = new TH1I("TrueClusNParticles","TrueClusNParticles",20,0,20);
867                                 fTrueList[iCut]->Add(fHistoTrueClusNParticles[iCut]);
868                                 fHistoTrueClusEMNonLeadingPt[iCut] = new TH1F("TrueClusEMNonLeading_Pt","TrueClusEMNonLeading_Pt",250,0,25);
869                                 fTrueList[iCut]->Add(fHistoTrueClusEMNonLeadingPt[iCut]);
870                                 fHistoTrueNLabelsInClus[iCut] = new TH1F("TrueNLabelsInClus","TrueNLabelsInClus",100,-0.5,99.5);
871                                 fTrueList[iCut]->Add(fHistoTrueNLabelsInClus[iCut]);    
872                         }       
873
874                         if(fDoMesonAnalysis){
875                                 fHistoTruePi0InvMassPt[iCut] = new TH2F("ESD_TruePi0_InvMass_Pt","ESD_TruePi0_InvMass_Pt",800,0,0.8,250,0,25);
876                                 fTrueList[iCut]->Add(fHistoTruePi0InvMassPt[iCut]);
877                                 fHistoTrueEtaInvMassPt[iCut] = new TH2F("ESD_TrueEta_InvMass_Pt","ESD_TrueEta_InvMass_Pt",800,0,0.8,250,0,25);
878                                 fTrueList[iCut]->Add(fHistoTrueEtaInvMassPt[iCut]);
879
880                                 fHistoTruePrimaryPi0InvMassPt[iCut] = new TH2F("ESD_TruePrimaryPi0_InvMass_Pt", "ESD_TruePrimaryPi0_InvMass_Pt", 800,0,0.8,250,0,25);
881                                 fHistoTruePrimaryPi0InvMassPt[iCut]->Sumw2();
882                                 fTrueList[iCut]->Add(fHistoTruePrimaryPi0InvMassPt[iCut]);
883                                 fHistoTruePrimaryEtaInvMassPt[iCut] = new TH2F("ESD_TruePrimaryEta_InvMass_Pt", "ESD_TruePrimaryEta_InvMass_Pt", 800,0,0.8,250,0,25);
884                                 fHistoTruePrimaryEtaInvMassPt[iCut]->Sumw2();
885                                 fTrueList[iCut]->Add(fHistoTruePrimaryEtaInvMassPt[iCut]);
886
887                                 fHistoTruePrimaryPi0W0WeightingInvMassPt[iCut] = new TH2F("ESD_TruePrimaryPi0W0Weights_InvMass_Pt", "ESD_TruePrimaryPi0W0Weights_InvMass_Pt", 800,0,0.8,250,0,25);
888                                 fHistoTruePrimaryPi0W0WeightingInvMassPt[iCut]->Sumw2();
889                                 fTrueList[iCut]->Add(fHistoTruePrimaryPi0W0WeightingInvMassPt[iCut]);
890                                 fHistoTruePrimaryEtaW0WeightingInvMassPt[iCut] = new TH2F("ESD_TruePrimaryEtaW0Weights_InvMass_Pt", "ESD_TruePrimaryEtaW0Weights_InvMass_Pt", 800,0,0.8,250,0,25);
891                                 fHistoTruePrimaryEtaW0WeightingInvMassPt[iCut]->Sumw2();
892                                 fTrueList[iCut]->Add(fHistoTruePrimaryEtaW0WeightingInvMassPt[iCut]);
893
894                                 fProfileTruePrimaryPi0WeightsInvMassPt[iCut] = new TProfile2D("ESD_TruePrimaryPi0Weights_InvMass_Pt", "ESD_TruePrimaryPi0Weights_InvMass_Pt", 800,0,0.8,250,0,25);
895                                 fProfileTruePrimaryPi0WeightsInvMassPt[iCut]->Sumw2();
896                                 fTrueList[iCut]->Add(fProfileTruePrimaryPi0WeightsInvMassPt[iCut]);
897                                 fProfileTruePrimaryEtaWeightsInvMassPt[iCut] = new TProfile2D("ESD_TruePrimaryEtaWeights_InvMass_Pt", "ESD_TruePrimaryEtaWeights_InvMass_Pt", 800,0,0.8,250,0,25);
898                                 fProfileTruePrimaryEtaWeightsInvMassPt[iCut]->Sumw2();
899                                 fTrueList[iCut]->Add(fProfileTruePrimaryEtaWeightsInvMassPt[iCut]);
900
901                                 fHistoTrueSecondaryPi0InvMassPt[iCut] = new TH2F("ESD_TrueSecondaryPi0_InvMass_Pt", "ESD_TrueSecondaryPi0_InvMass_Pt", 800,0,0.8,250,0,25);
902                                 fHistoTrueSecondaryPi0InvMassPt[iCut]->Sumw2();
903                                 fTrueList[iCut]->Add(fHistoTrueSecondaryPi0InvMassPt[iCut]);
904                                 fHistoTrueSecondaryEtaInvMassPt[iCut] = new TH2F("ESD_TrueSecondaryEta_InvMass_Pt", "ESD_TrueSecondaryEta_InvMass_Pt", 800,0,0.8,250,0,25);
905                                 fHistoTrueSecondaryEtaInvMassPt[iCut]->Sumw2();
906                                 fTrueList[iCut]->Add(fHistoTrueSecondaryEtaInvMassPt[iCut]);
907
908                                 fHistoTrueSecondaryPi0FromK0sInvMassPt[iCut] = new TH2F("ESD_TrueSecondaryPi0FromK0s_InvMass_Pt","ESD_TrueSecondaryPi0FromK0s_InvMass_Pt",800,0,0.8,250,0,25);
909                                 fHistoTrueSecondaryPi0FromK0sInvMassPt[iCut]->Sumw2();
910                                 fTrueList[iCut]->Add(fHistoTrueSecondaryPi0FromK0sInvMassPt[iCut]);
911                                 fHistoTrueSecondaryPi0FromEtaInvMassPt[iCut] = new TH2F("ESD_TrueSecondaryPi0FromEta_InvMass_Pt","ESD_TrueSecondaryPi0FromEta_InvMass_Pt",800,0,0.8,250,0,25);
912                                 fTrueList[iCut]->Add(fHistoTrueSecondaryPi0FromEtaInvMassPt[iCut]);
913                                 fHistoTrueSecondaryPi0FromLambdaInvMassPt[iCut] = new TH2F("ESD_TrueSecondaryPi0FromLambda_InvMass_Pt","ESD_TrueSecondaryPi0FromLambda_InvMass_Pt",800,0,0.8,250,0,25);
914                                 fTrueList[iCut]->Add(fHistoTrueSecondaryPi0FromLambdaInvMassPt[iCut]);
915                                 if (fDoMesonQA > 0){
916                                         fHistoTruePi0CaloPhotonInvMassPt[iCut] = new TH2F("ESD_TruePi0CaloPhoton_InvMass_Pt","ESD_TruePi0CaloPhoton_InvMass_Pt",800,0,0.8,250,0,25);
917                                         fTrueList[iCut]->Add(fHistoTruePi0CaloPhotonInvMassPt[iCut]);
918                                         fHistoTrueEtaCaloPhotonInvMassPt[iCut] = new TH2F("ESD_TrueEtaCaloPhoton_InvMass_Pt","ESD_TrueEtaCaloPhoton_InvMass_Pt",800,0,0.8,250,0,25);
919                                         fTrueList[iCut]->Add(fHistoTrueEtaCaloPhotonInvMassPt[iCut]);
920
921                                         fHistoTruePi0CaloMixedPhotonConvPhotonInvMassPt[iCut] = new TH2F("ESD_TruePi0CaloMixedPhotonConvertedPhoton_InvMass_Pt","ESD_TruePi0CaloMixedPhotonConvertedPhoton_InvMass_Pt",800,0,0.8,250,0,25);
922                                         fTrueList[iCut]->Add(fHistoTruePi0CaloMixedPhotonConvPhotonInvMassPt[iCut]);
923                                         fHistoTrueEtaCaloMixedPhotonConvPhotonInvMassPt[iCut] = new TH2F("ESD_TrueEtaCaloMixedPhotonConvertedPhoton_InvMass_Pt","ESD_TrueEtaCaloMixedPhotonConvertedPhoton_InvMass_Pt",800,0,0.8,250,0,25);
924                                         fTrueList[iCut]->Add(fHistoTrueEtaCaloMixedPhotonConvPhotonInvMassPt[iCut]);
925
926                                         fHistoTruePi0CaloConvertedPhotonInvMassPt[iCut] = new TH2F("ESD_TruePi0CaloConvertedPhoton_InvMass_Pt","ESD_TruePi0CaloConvertedPhoton_InvMass_Pt",800,0,0.8,250,0,25);
927                                         fTrueList[iCut]->Add(fHistoTruePi0CaloConvertedPhotonInvMassPt[iCut]);
928                                         fHistoTrueEtaCaloConvertedPhotonInvMassPt[iCut] = new TH2F("ESD_TrueEtaCaloConvertedPhoton_InvMass_Pt","ESD_TrueEtaCaloConvertedPhoton_InvMass_Pt",800,0,0.8,250,0,25);
929                                         fTrueList[iCut]->Add(fHistoTrueEtaCaloConvertedPhotonInvMassPt[iCut]);
930
931                                         fHistoTruePi0CaloElectronInvMassPt[iCut] = new TH2F("ESD_TruePi0CaloElectron_InvMass_Pt","ESD_TruePi0CaloElectron_InvMass_Pt",800,0,0.8,250,0,25);
932                                         fTrueList[iCut]->Add(fHistoTruePi0CaloElectronInvMassPt[iCut]);
933                                         fHistoTrueEtaCaloElectronInvMassPt[iCut] = new TH2F("ESD_TrueEtaCaloElectron_InvMass_Pt","ESD_TrueEtaCaloElectron_InvMass_Pt",800,0,0.8,250,0,25);
934                                         fTrueList[iCut]->Add(fHistoTrueEtaCaloElectronInvMassPt[iCut]);
935
936                                         fHistoTruePi0CaloMergedClusterInvMassPt[iCut] = new TH2F("ESD_TruePi0CaloMergedCluster_InvMass_Pt","ESD_TruePi0CaloMergedCluster_InvMass_Pt",800,0,0.8,250,0,25);
937                                         fTrueList[iCut]->Add(fHistoTruePi0CaloMergedClusterInvMassPt[iCut]);
938                                         fHistoTrueEtaCaloMergedClusterInvMassPt[iCut] = new TH2F("ESD_TrueEtaCaloMergedCluster_InvMass_Pt","ESD_TrueEtaCaloMergedCluster_InvMass_Pt",800,0,0.8,250,0,25);
939                                         fTrueList[iCut]->Add(fHistoTrueEtaCaloMergedClusterInvMassPt[iCut]);
940
941                                         fHistoTruePi0CaloMergedClusterPartConvInvMassPt[iCut] = new TH2F("ESD_TruePi0CaloMergedClusterPartConv_InvMass_Pt","ESD_TruePi0CaloMergedClusterPartConv_InvMass_Pt",800,0,0.8,250,0,25);
942                                         fTrueList[iCut]->Add(fHistoTruePi0CaloMergedClusterPartConvInvMassPt[iCut]);
943                                         fHistoTrueEtaCaloMergedClusterPartConvInvMassPt[iCut] = new TH2F("ESD_TrueEtaCaloMergedClusterPartConv_InvMass_Pt","ESD_TrueEtaCaloMergedClusterPartConv_InvMass_Pt",800,0,0.8,250,0,25);
944                                         fTrueList[iCut]->Add(fHistoTrueEtaCaloMergedClusterPartConvInvMassPt[iCut]);
945
946                                         
947                                         fHistoTruePrimaryPi0MCPtResolPt[iCut] = new TH2F("ESD_TruePrimaryPi0_MCPt_ResolPt","ESD_TruePrimaryPi0_ResolPt_MCPt",500,0.03,25,1000,-1.,1.);
948                                         fHistoTruePrimaryPi0MCPtResolPt[iCut]->Sumw2();
949                                         SetLogBinningXTH2(fHistoTruePrimaryPi0MCPtResolPt[iCut]);
950                                         fTrueList[iCut]->Add(fHistoTruePrimaryPi0MCPtResolPt[iCut]);
951                                         fHistoTruePrimaryEtaMCPtResolPt[iCut]  = new TH2F("ESD_TruePrimaryEta_MCPt_ResolPt","ESD_TruePrimaryEta_ResolPt_MCPt",500,0.03,25,1000,-1.,1.);
952                                         fHistoTruePrimaryEtaMCPtResolPt[iCut]->Sumw2();
953                                         SetLogBinningXTH2(fHistoTruePrimaryEtaMCPtResolPt[iCut]);
954                                         fTrueList[iCut]->Add(fHistoTruePrimaryEtaMCPtResolPt[iCut]);
955                                         fHistoTrueBckGGInvMassPt[iCut] = new TH2F("ESD_TrueBckGG_InvMass_Pt","ESD_TrueBckGG_InvMass_Pt",800,0,0.8,250,0,25);
956                                         fTrueList[iCut]->Add(fHistoTrueBckGGInvMassPt[iCut]);
957                                         fHistoTrueBckContInvMassPt[iCut] = new TH2F("ESD_TrueBckCont_InvMass_Pt","ESD_TrueBckCont_InvMass_Pt",800,0,0.8,250,0,25);
958                                         fTrueList[iCut]->Add(fHistoTrueBckContInvMassPt[iCut]);
959                                         fHistoTrueK0sWithPi0DaughterMCPt[iCut] = new TH1F("ESD_TrueK0sWithPi0Daughter_MCPt","ESD_TrueK0sWithPi0Daughter_MCPt",250,0,25);
960                                         fTrueList[iCut]->Add(fHistoTrueK0sWithPi0DaughterMCPt[iCut]);
961                                         fHistoTrueEtaWithPi0DaughterMCPt[iCut] = new TH1F("ESD_TrueEtaWithPi0Daughter_MCPt","ESD_TrueEtaWithPi0Daughter_MCPt",250,0,25);
962                                         fTrueList[iCut]->Add(fHistoTrueEtaWithPi0DaughterMCPt[iCut]);
963                                         fHistoTrueLambdaWithPi0DaughterMCPt[iCut] = new TH1F("ESD_TrueLambdaWithPi0Daughter_MCPt","ESD_TrueLambdaWithPi0Daughter_MCPt",250,0,25);
964                                         fTrueList[iCut]->Add(fHistoTrueLambdaWithPi0DaughterMCPt[iCut]);
965           
966                                         fHistoTruePi0PtY[iCut] = new TH2F("ESD_TruePi0_Pt_Y","ESD_TruePi0_Pt_Y",150,0.03,15.,150,-1.5,1.5);
967                                         SetLogBinningXTH2(fHistoTruePi0PtY[iCut]);
968                                         fTrueList[iCut]->Add(fHistoTruePi0PtY[iCut]);
969                                         fHistoTrueEtaPtY[iCut] = new TH2F("ESD_TrueEta_Pt_Y","ESD_TrueEta_Pt_Y",150,0.03,15.,150,-1.5,1.5);
970                                         SetLogBinningXTH2(fHistoTrueEtaPtY[iCut]);
971                                         fTrueList[iCut]->Add(fHistoTrueEtaPtY[iCut]);
972                                         fHistoTruePi0PtAlpha[iCut] = new TH2F("ESD_TruePi0_Pt_Alpha","ESD_TruePi0_Pt_Alpha",150,0.03,15.,100,0,1);
973                                         SetLogBinningXTH2(fHistoTruePi0PtAlpha[iCut]);
974                                         fTrueList[iCut]->Add(fHistoTruePi0PtAlpha[iCut]);
975                                         fHistoTrueEtaPtAlpha[iCut] = new TH2F("ESD_TrueEta_Pt_Alpha","ESD_TrueEta_Pt_Alpha",150,0.03,15.,100,0,1);
976                                         SetLogBinningXTH2(fHistoTrueEtaPtAlpha[iCut]);
977                                         fTrueList[iCut]->Add(fHistoTrueEtaPtAlpha[iCut]);
978                                         
979                                         fHistoTruePi0PtOpenAngle[iCut] = new TH2F("ESD_TruePi0_Pt_OpenAngle","ESD_TruePi0_Pt_OpenAngle",150,0.03,15.,200,0,2*TMath::Pi());
980                                         SetLogBinningXTH2(fHistoTruePi0PtOpenAngle[iCut]);
981                                         fTrueList[iCut]->Add(fHistoTruePi0PtOpenAngle[iCut]);
982                                         fHistoTrueEtaPtOpenAngle[iCut] = new TH2F("ESD_TrueEta_Pt_OpenAngle","ESD_TrueEta_Pt_OpenAngle",150,0.03,15.,200,0,2*TMath::Pi());
983                                         SetLogBinningXTH2(fHistoTrueEtaPtOpenAngle[iCut]);
984                                         fTrueList[iCut]->Add(fHistoTrueEtaPtOpenAngle[iCut]);
985                                 }
986                         }
987                 }
988         }  
989     
990         fV0Reader=(AliV0ReaderV1*)AliAnalysisManager::GetAnalysisManager()->GetTask("V0ReaderV1");
991         if(!fV0Reader){printf("Error: No V0 Reader");return;} // GetV0Reader
992   
993         if(fV0Reader)
994                 if((AliConversionPhotonCuts*)fV0Reader->GetConversionCuts())
995                         if(((AliConversionPhotonCuts*)fV0Reader->GetConversionCuts())->GetCutHistograms())
996                                 fOutputContainer->Add(((AliConversionPhotonCuts*)fV0Reader->GetConversionCuts())->GetCutHistograms());
997         if(fV0Reader)
998                 if((AliConvEventCuts*)fV0Reader->GetEventCuts())
999                         if(((AliConvEventCuts*)fV0Reader->GetEventCuts())->GetCutHistograms())
1000                                 fOutputContainer->Add(((AliConvEventCuts*)fV0Reader->GetEventCuts())->GetCutHistograms());
1001
1002                         
1003         for(Int_t iCut = 0; iCut<fnCuts;iCut++){
1004                 if(!((AliConvEventCuts*)fEventCutArray->At(iCut))) continue;
1005                 if(((AliConvEventCuts*)fEventCutArray->At(iCut))->GetCutHistograms()){
1006                         fCutFolder[iCut]->Add(((AliConvEventCuts*)fEventCutArray->At(iCut))->GetCutHistograms());
1007                 }
1008                 if(!((AliCaloPhotonCuts*)fClusterCutArray->At(iCut))) continue;
1009                 if(((AliCaloPhotonCuts*)fClusterCutArray->At(iCut))->GetCutHistograms()){
1010                         fCutFolder[iCut]->Add(((AliCaloPhotonCuts*)fClusterCutArray->At(iCut))->GetCutHistograms());
1011                 }
1012                 if(fDoMesonAnalysis){
1013                         if(!((AliConversionMesonCuts*)fMesonCutArray->At(iCut))) continue;
1014                         if(((AliConversionMesonCuts*)fMesonCutArray->At(iCut))->GetCutHistograms()){
1015                                 fCutFolder[iCut]->Add(((AliConversionMesonCuts*)fMesonCutArray->At(iCut))->GetCutHistograms());
1016                         }
1017                 }
1018         }
1019         PostData(1, fOutputContainer);
1020 }
1021 //_____________________________________________________________________________
1022 Bool_t AliAnalysisTaskGammaCalo::Notify()
1023 {
1024         for(Int_t iCut = 0; iCut<fnCuts;iCut++){
1025                 if(!((AliConvEventCuts*)fEventCutArray->At(iCut))->GetDoEtaShift()){
1026                         fProfileEtaShift[iCut]->Fill(0.,(((AliConvEventCuts*)fEventCutArray->At(iCut))->GetEtaShift()));
1027                         continue; // No Eta Shift requested, continue
1028                 }
1029                 if(((AliConvEventCuts*)fEventCutArray->At(iCut))->GetEtaShift() == 0.0){ // Eta Shift requested but not set, get shift automatically
1030                         ((AliConvEventCuts*)fEventCutArray->At(iCut))->GetCorrectEtaShiftFromPeriod(fV0Reader->GetPeriodName());
1031                         fProfileEtaShift[iCut]->Fill(0.,(((AliConvEventCuts*)fEventCutArray->At(iCut))->GetEtaShift()));
1032                         ((AliConvEventCuts*)fEventCutArray->At(iCut))->DoEtaShift(kFALSE); // Eta Shift Set, make sure that it is called only once
1033                         continue;
1034                 }
1035                 else{
1036                         printf(" Gamma Conversion Task %s :: Eta Shift Manually Set to %f \n\n",
1037                                         (((AliConvEventCuts*)fEventCutArray->At(iCut))->GetCutNumber()).Data(),((AliConvEventCuts*)fEventCutArray->At(iCut))->GetEtaShift());
1038                         fProfileEtaShift[iCut]->Fill(0.,(((AliConvEventCuts*)fEventCutArray->At(iCut))->GetEtaShift()));
1039                         ((AliConvEventCuts*)fEventCutArray->At(iCut))->DoEtaShift(kFALSE); // Eta Shift Set, make sure that it is called only once
1040                 }
1041         }
1042         
1043         return kTRUE;
1044 }
1045 //_____________________________________________________________________________
1046 void AliAnalysisTaskGammaCalo::UserExec(Option_t *)
1047 {
1048         //
1049         // Called for each event
1050         //
1051         Int_t eventQuality = ((AliConvEventCuts*)fV0Reader->GetEventCuts())->GetEventQuality();
1052         if(eventQuality == 2 || eventQuality == 3){// Event Not Accepted due to MC event missing or wrong trigger for V0ReaderV1
1053                 for(Int_t iCut = 0; iCut<fnCuts; iCut++){
1054                 fHistoNEvents[iCut]->Fill(eventQuality);
1055                 }
1056                 return;
1057         }
1058         
1059         if(fIsMC) fMCEvent = MCEvent();
1060         if(fMCEvent == NULL) fIsMC = kFALSE;
1061         
1062         fInputEvent = InputEvent();
1063         
1064         if(fIsMC && fInputEvent->IsA()==AliESDEvent::Class()){
1065                 fMCStack = fMCEvent->Stack();
1066                 if(fMCStack == NULL) fIsMC = kFALSE;
1067         }
1068         
1069         // ------------------- BeginEvent ----------------------------
1070         
1071         AliEventplane *EventPlane = fInputEvent->GetEventplane();
1072         if(fIsHeavyIon ==1)fEventPlaneAngle = EventPlane->GetEventplane("V0",fInputEvent,2);
1073         else fEventPlaneAngle=0.0;
1074         
1075         for(Int_t iCut = 0; iCut<fnCuts; iCut++){
1076                 
1077                 fiCut = iCut;
1078                 Int_t eventNotAccepted = ((AliConvEventCuts*)fEventCutArray->At(iCut))->IsEventAcceptedByCut(fV0Reader->GetEventCuts(),fInputEvent,fMCEvent,fIsHeavyIon);
1079                 
1080                 if(eventNotAccepted){
1081                 // cout << "event rejected due to wrong trigger: " <<eventNotAccepted << endl;
1082                         fHistoNEvents[iCut]->Fill(eventNotAccepted); // Check Centrality, PileUp, SDD and V0AND --> Not Accepted => eventQuality = 1
1083                         continue;
1084                 }
1085
1086                 if(eventQuality != 0){// Event Not Accepted
1087                         //cout << "event rejected due to: " <<eventQuality << endl;
1088                         fHistoNEvents[iCut]->Fill(eventQuality);
1089                         continue;
1090                 }
1091
1092                 fHistoNEvents[iCut]->Fill(eventQuality); // Should be 0 here
1093                 fHistoNGoodESDTracks[iCut]->Fill(fV0Reader->GetNumberOfPrimaryTracks());
1094                 if(((AliConvEventCuts*)fEventCutArray->At(iCut))->IsHeavyIon() == 2)    fHistoNV0Tracks[iCut]->Fill(fInputEvent->GetVZEROData()->GetMTotV0A());
1095                         else fHistoNV0Tracks[iCut]->Fill(fInputEvent->GetVZEROData()->GetMTotV0A()+fInputEvent->GetVZEROData()->GetMTotV0C());
1096
1097                 if(fIsMC){
1098                         // Process MC Particle
1099                         if(((AliConvEventCuts*)fEventCutArray->At(iCut))->GetSignalRejection() != 0){
1100                                 if(fInputEvent->IsA()==AliESDEvent::Class()){
1101                                 ((AliConvEventCuts*)fEventCutArray->At(iCut))->GetNotRejectedParticles(((AliConvEventCuts*)fEventCutArray->At(iCut))->GetSignalRejection(),
1102                                                                                                                                                                         ((AliConvEventCuts*)fEventCutArray->At(iCut))->GetAcceptedHeader(),
1103                                                                                                                                                                         fMCEvent);
1104                                 }
1105                                 else if(fInputEvent->IsA()==AliAODEvent::Class()){
1106                                 ((AliConvEventCuts*)fEventCutArray->At(iCut))->GetNotRejectedParticles(((AliConvEventCuts*)fEventCutArray->At(iCut))->GetSignalRejection(),
1107                                                                                                                                                                         ((AliConvEventCuts*)fEventCutArray->At(iCut))->GetAcceptedHeader(),
1108                                                                                                                                                                         fInputEvent);
1109                                 }
1110
1111                                 if(((AliConvEventCuts*)fEventCutArray->At(iCut))->GetAcceptedHeader()){
1112                                         for(Int_t i = 0;i<(((AliConvEventCuts*)fEventCutArray->At(iCut))->GetAcceptedHeader())->GetEntries();i++){
1113                                                 TString nameBin= fHistoMCHeaders[iCut]->GetXaxis()->GetBinLabel(i+1);
1114                                                 if (nameBin.CompareTo("")== 0){
1115                                                         TString nameHeader = ((TObjString*)((TList*)((AliConvEventCuts*)fEventCutArray->At(iCut))
1116                                                                                                                                 ->GetAcceptedHeader())->At(i))->GetString();
1117                                                         fHistoMCHeaders[iCut]->GetXaxis()->SetBinLabel(i+1,nameHeader.Data());
1118                                                 }
1119                                         }
1120                                 }
1121                         }
1122                 }
1123                 if(fIsMC){
1124                 if(fInputEvent->IsA()==AliESDEvent::Class())
1125                         ProcessMCParticles();
1126                 if(fInputEvent->IsA()==AliAODEvent::Class())
1127                         ProcessAODMCParticles();
1128                 }
1129                 
1130                 // 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)
1131                 ProcessClusters();                                      // process calo clusters
1132
1133                 fHistoNGammaCandidates[iCut]->Fill(fClusterCandidates->GetEntries());
1134                 fHistoNGoodESDTracksVsNGammaCanditates[iCut]->Fill(fV0Reader->GetNumberOfPrimaryTracks(),fClusterCandidates->GetEntries());
1135                 if(fDoMesonAnalysis){ // Meson Analysis
1136
1137                         CalculatePi0Candidates(); // Combine Gammas from conversion and from calo
1138
1139                         if(((AliConversionMesonCuts*)fMesonCutArray->At(iCut))->DoBGCalculation()){
1140                                 if(((AliConversionMesonCuts*)fMesonCutArray->At(iCut))->BackgroundHandlerType() == 0){
1141                                         CalculateBackground(); // Combinatorial Background
1142                                         UpdateEventByEventData(); // Store Event for mixed Events
1143                                 }
1144                                 
1145                         }
1146                 }
1147
1148                 fClusterCandidates->Clear(); // delete cluster candidates
1149         }
1150         
1151         PostData(1, fOutputContainer);
1152 }
1153
1154 //________________________________________________________________________
1155 void AliAnalysisTaskGammaCalo::ProcessClusters()
1156 {
1157         
1158         Int_t nclus = 0;
1159         nclus = fInputEvent->GetNumberOfCaloClusters();
1160         
1161 //      cout << nclus << endl;
1162         
1163         if(nclus == 0)  return;
1164         
1165         // vertex
1166         Double_t vertex[3] = {0};
1167         InputEvent()->GetPrimaryVertex()->GetXYZ(vertex);
1168         
1169         // Loop over EMCal clusters
1170         for(Long_t i = 0; i < nclus; i++){
1171                 
1172                 AliVCluster* clus = NULL;
1173                 clus = fInputEvent->GetCaloCluster(i);          
1174                 if (!clus) continue;
1175                 if(!((AliCaloPhotonCuts*)fClusterCutArray->At(fiCut))->ClusterIsSelected(clus,fInputEvent,fIsMC)) continue;
1176                 // TLorentzvector with cluster
1177                 TLorentzVector clusterVector;
1178                 clus->GetMomentum(clusterVector,vertex);
1179                 
1180                 TLorentzVector* tmpvec = new TLorentzVector();
1181                 tmpvec->SetPxPyPzE(clusterVector.Px(),clusterVector.Py(),clusterVector.Pz(),clusterVector.E());
1182                 
1183                 // convert to AODConversionPhoton
1184                 AliAODConversionPhoton *PhotonCandidate=new AliAODConversionPhoton(tmpvec);
1185                 if(!PhotonCandidate) continue;
1186                 
1187                 // Flag Photon as CaloPhoton
1188                 PhotonCandidate->SetIsCaloPhoton();
1189                 PhotonCandidate->SetCaloClusterRef(i);
1190                 // get MC label
1191                 if(fIsMC){
1192                         Int_t* mclabelsCluster = clus->GetLabels();
1193                         PhotonCandidate->SetNCaloPhotonMCLabels(clus->GetNLabels());
1194 //                      cout << clus->GetNLabels() << endl;
1195                         if (clus->GetNLabels()>0){
1196                                 for (Int_t k =0; k< (Int_t)clus->GetNLabels(); k++){
1197                                         if (k< 20)PhotonCandidate->SetCaloPhotonMCLabel(k,mclabelsCluster[k]);
1198 //                                      Int_t pdgCode = fMCStack->Particle(mclabelsCluster[k])->GetPdgCode();
1199 //                                      cout << "label " << k << "\t" << mclabelsCluster[k] << " pdg code: " << pdgCode << endl;
1200                                 }       
1201                         }
1202                 }
1203                 
1204                 fIsFromMBHeader = kTRUE; 
1205                 // test whether largest contribution to cluster orginates in added signals
1206                 if (fIsMC && ((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsParticleFromBGEvent(PhotonCandidate->GetCaloPhotonMCLabel(0), fMCStack, fInputEvent) == 0) fIsFromMBHeader = kFALSE;
1207                 
1208                 if (fIsFromMBHeader)fHistoClusGammaPt[fiCut]->Fill(PhotonCandidate->Pt());
1209                 fClusterCandidates->Add(PhotonCandidate); // if no second loop is required add to events good gammas
1210                 
1211                 if(fIsMC){
1212                         if(fInputEvent->IsA()==AliESDEvent::Class()){
1213                                 ProcessTrueClusterCandidates(PhotonCandidate);
1214                         } else {
1215                                 ProcessTrueClusterCandidatesAOD(PhotonCandidate);
1216                         }       
1217                 }
1218                 
1219                 delete tmpvec;
1220         }
1221         
1222 }
1223
1224 //________________________________________________________________________
1225 void AliAnalysisTaskGammaCalo::ProcessTrueClusterCandidates(AliAODConversionPhoton *TruePhotonCandidate)
1226 {
1227                 
1228         TParticle *Photon = NULL;
1229         if (!TruePhotonCandidate->GetIsCaloPhoton()) AliFatal("CaloPhotonFlag has not been set task will abort");
1230         fHistoTrueNLabelsInClus[fiCut]->Fill(TruePhotonCandidate->GetNCaloPhotonMCLabels());
1231         
1232         if (TruePhotonCandidate->GetNCaloPhotonMCLabels()>0)Photon = fMCStack->Particle(TruePhotonCandidate->GetCaloPhotonMCLabel(0));
1233                 else return;
1234                 
1235         if(Photon == NULL){
1236         //    cout << "no photon" << endl;
1237                 return;
1238         }
1239
1240         Int_t pdgCodeParticle = Photon->GetPdgCode();
1241         TruePhotonCandidate->SetCaloPhotonMCFlags(fMCStack);
1242         
1243         // True Photon
1244         if(fIsFromMBHeader){
1245                 if (TruePhotonCandidate->IsLargestComponentPhoton() || TruePhotonCandidate->IsLargestComponentElectron() )fHistoTrueClusGammaPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1246                         else fHistoTrueClusEMNonLeadingPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1247                 if (fDoClusterQA > 0){
1248                         if (TruePhotonCandidate->IsLargestComponentPhoton()){ 
1249                                 fHistoTrueClusUnConvGammaPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1250                                 fHistoTrueClusUnConvGammaMCPt[fiCut]->Fill(Photon->Pt());
1251                         }       
1252                         if (TruePhotonCandidate->IsLargestComponentElectron()) 
1253                                 fHistoTrueClusElectronPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1254                         if (TruePhotonCandidate->IsLargestComponentElectron() && TruePhotonCandidate->IsConversion()){
1255                                 fHistoTrueClusConvGammaPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1256                                 fHistoTrueClusConvGammaMCPt[fiCut]->Fill(((TParticle*)fMCStack->Particle(Photon->GetMother(0)))->Pt());
1257                         }       
1258                         if (TruePhotonCandidate->IsLargestComponentElectron() && TruePhotonCandidate->IsConversion() && TruePhotonCandidate->IsConversionFullyContained()) 
1259                                 fHistoTrueClusConvGammaFullyPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1260                         if (TruePhotonCandidate->IsMerged() || TruePhotonCandidate->IsMergedPartConv() || TruePhotonCandidate->IsDalitzMerged())
1261                                 fHistoTrueClusMergedGammaPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1262                         if (TruePhotonCandidate->IsMergedPartConv())
1263                                 fHistoTrueClusMergedPartConvGammaPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1264                         if (TruePhotonCandidate->IsDalitz()) 
1265                                 fHistoTrueClusDalitzPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1266                         if (TruePhotonCandidate->IsDalitzMerged()) 
1267                                 fHistoTrueClusDalitzMergedPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1268                         if (TruePhotonCandidate->IsPhotonWithElecMother()) 
1269                                 fHistoTrueClusPhotonFromElecMotherPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1270                         if (TruePhotonCandidate->IsShower()) 
1271                                 fHistoTrueClusShowerPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1272                         if (TruePhotonCandidate->IsSubLeadingEM())
1273                                 fHistoTrueClusSubLeadingPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1274                         fHistoTrueClusNParticles[fiCut]->Fill(TruePhotonCandidate->GetNCaloPhotonMotherMCLabels());
1275                         if (!TruePhotonCandidate->IsLargestComponentPhoton())
1276                                 FillPhotonBackgroundHist(TruePhotonCandidate,pdgCodeParticle);
1277                         if (!(TruePhotonCandidate->IsLargestComponentPhoton() || (TruePhotonCandidate->IsLargestComponentElectron() && TruePhotonCandidate->IsConversion())) )
1278                                 FillPhotonPlusConversionBackgroundHist(TruePhotonCandidate,pdgCodeParticle);
1279                 }
1280                 
1281                 if(Photon->GetMother(0) <= fMCStack->GetNprimary()){
1282                         if (TruePhotonCandidate->IsLargestComponentPhoton()){
1283                                 fHistoTruePrimaryClusGammaPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1284                                 fHistoTruePrimaryClusGammaESDPtMCPt[fiCut]->Fill(TruePhotonCandidate->Pt(),Photon->Pt()); // Allways Filled
1285                         }
1286                         if (TruePhotonCandidate->IsLargestComponentElectron() && TruePhotonCandidate->IsConversion()){
1287                                 fHistoTruePrimaryClusConvGammaPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1288                                 fHistoTruePrimaryClusConvGammaESDPtMCPt[fiCut]->Fill(TruePhotonCandidate->Pt(),Photon->Pt()); // Allways Filled
1289                         }
1290                         
1291                 } else {
1292                         if (TruePhotonCandidate->IsLargestComponentPhoton())
1293                                 fHistoTrueSecondaryClusGammaPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1294                         if (TruePhotonCandidate->IsLargestComponentElectron() && TruePhotonCandidate->IsConversion())
1295                                 fHistoTrueSecondaryClusConvGammaPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1296                         if(fMCStack->Particle(Photon->GetMother(0))->GetMother(0) > -1){
1297                                 if(fMCStack->Particle(fMCStack->Particle(Photon->GetMother(0))->GetMother(0))->GetPdgCode() == 3122){
1298                                         if (TruePhotonCandidate->IsLargestComponentPhoton())
1299                                                 fHistoTrueSecondaryClusGammaFromXFromLambdaPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1300                                         if (TruePhotonCandidate->IsLargestComponentElectron() && TruePhotonCandidate->IsConversion())
1301                                                 fHistoTrueSecondaryClusConvGammaFromXFromLambdaPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1302                                 }
1303                                 if(fMCStack->Particle(fMCStack->Particle(Photon->GetMother(0))->GetMother(0))->GetPdgCode() == 310){
1304                                         if (TruePhotonCandidate->IsLargestComponentPhoton())
1305                                                 fHistoTrueSecondaryClusGammaFromXFromK0sPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1306                                         if (TruePhotonCandidate->IsLargestComponentElectron() && TruePhotonCandidate->IsConversion())
1307                                                 fHistoTrueSecondaryClusConvGammaFromXFromK0sPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1308                                 }
1309                                 if(fMCStack->Particle(fMCStack->Particle(Photon->GetMother(0))->GetMother(0))->GetPdgCode() == 221){
1310                                         if (TruePhotonCandidate->IsLargestComponentPhoton())
1311                                                 fHistoTrueSecondaryClusGammaFromXFromEtasPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1312                                         if (TruePhotonCandidate->IsLargestComponentElectron() && TruePhotonCandidate->IsConversion())
1313                                                 fHistoTrueSecondaryClusConvGammaFromXFromEtasPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1314                                                 
1315                                 }       
1316                         }
1317                 }
1318         }
1319         return;
1320 }
1321
1322
1323 //________________________________________________________________________
1324 void AliAnalysisTaskGammaCalo::ProcessTrueClusterCandidatesAOD(AliAODConversionPhoton *TruePhotonCandidate)
1325 {
1326         AliAODMCParticle *Photon = NULL;
1327         TClonesArray *AODMCTrackArray = dynamic_cast<TClonesArray*>(fInputEvent->FindListObject(AliAODMCParticle::StdBranchName()));
1328         if (AODMCTrackArray){
1329                 if (!TruePhotonCandidate->GetIsCaloPhoton()) AliFatal("CaloPhotonFlag has not been set task will abort");
1330                 if (TruePhotonCandidate->GetNCaloPhotonMCLabels()>0) Photon = (AliAODMCParticle*) AODMCTrackArray->At(TruePhotonCandidate->GetCaloPhotonMCLabel(0));
1331                         else return;
1332         } else {
1333                 AliInfo("AODMCTrackArray could not be loaded");
1334                 return;
1335         }
1336
1337         if(Photon == NULL){
1338         //      cout << "no photon" << endl;
1339                 return;
1340         }
1341         Int_t pdgCodeParticle = Photon->GetPdgCode();
1342         TruePhotonCandidate->SetCaloPhotonMCFlagsAOD(fInputEvent);
1343         
1344         // True Photon
1345         if(fIsFromMBHeader){
1346                 if (TruePhotonCandidate->IsLargestComponentPhoton() || TruePhotonCandidate->IsLargestComponentElectron() )fHistoTrueClusGammaPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1347                         else fHistoTrueClusEMNonLeadingPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1348                 if (fDoClusterQA > 0){
1349                         if (TruePhotonCandidate->IsLargestComponentPhoton()) {
1350                                 fHistoTrueClusUnConvGammaPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1351                                 fHistoTrueClusUnConvGammaMCPt[fiCut]->Fill(Photon->Pt());
1352                         }       
1353                         if (TruePhotonCandidate->IsLargestComponentElectron()) 
1354                                 fHistoTrueClusElectronPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1355                         if (TruePhotonCandidate->IsLargestComponentElectron() && TruePhotonCandidate->IsConversion()) {
1356                                 fHistoTrueClusConvGammaPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1357                                 AliAODMCParticle *Mother = (AliAODMCParticle*) AODMCTrackArray->At(Photon->GetMother());
1358                                 fHistoTrueClusConvGammaMCPt[fiCut]->Fill(Mother->Pt());
1359                         }       
1360                         if (TruePhotonCandidate->IsLargestComponentElectron() && TruePhotonCandidate->IsConversion() && TruePhotonCandidate->IsConversionFullyContained()) 
1361                                 fHistoTrueClusConvGammaFullyPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1362                         if (TruePhotonCandidate->IsMerged() || TruePhotonCandidate->IsMergedPartConv() || TruePhotonCandidate->IsDalitzMerged())
1363                                 fHistoTrueClusMergedGammaPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1364                         if (TruePhotonCandidate->IsMergedPartConv())
1365                                 fHistoTrueClusMergedPartConvGammaPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1366                         if (TruePhotonCandidate->IsDalitz()) 
1367                                 fHistoTrueClusDalitzPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1368                         if (TruePhotonCandidate->IsDalitzMerged()) 
1369                                 fHistoTrueClusDalitzMergedPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1370                         if (TruePhotonCandidate->IsPhotonWithElecMother()) 
1371                                 fHistoTrueClusPhotonFromElecMotherPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1372                         if (TruePhotonCandidate->IsShower()) 
1373                                 fHistoTrueClusShowerPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1374                         if (TruePhotonCandidate->IsSubLeadingEM())
1375                                 fHistoTrueClusSubLeadingPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1376                         fHistoTrueClusNParticles[fiCut]->Fill(TruePhotonCandidate->GetNCaloPhotonMotherMCLabels());
1377                         
1378                         if (!TruePhotonCandidate->IsLargestComponentPhoton())
1379                                 FillPhotonBackgroundHist(TruePhotonCandidate,pdgCodeParticle);
1380                         if (!(TruePhotonCandidate->IsLargestComponentPhoton() || (TruePhotonCandidate->IsLargestComponentElectron() && TruePhotonCandidate->IsConversion())) )
1381                                 FillPhotonPlusConversionBackgroundHist(TruePhotonCandidate,pdgCodeParticle);
1382                 }
1383                 
1384                 if(Photon->IsPrimary()){
1385                         if (TruePhotonCandidate->IsLargestComponentPhoton()){
1386                                 fHistoTruePrimaryClusGammaPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1387                                 fHistoTruePrimaryClusGammaESDPtMCPt[fiCut]->Fill(TruePhotonCandidate->Pt(),Photon->Pt()); // Allways Filled
1388                         }
1389                         if (TruePhotonCandidate->IsLargestComponentElectron() && TruePhotonCandidate->IsConversion()){
1390                                 fHistoTruePrimaryClusConvGammaPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1391                                 fHistoTruePrimaryClusConvGammaESDPtMCPt[fiCut]->Fill(TruePhotonCandidate->Pt(),Photon->Pt()); // Allways Filled
1392                         }
1393                         
1394                 } else {
1395                         if (TruePhotonCandidate->IsLargestComponentPhoton())
1396                                 fHistoTrueSecondaryClusGammaPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1397                         if (TruePhotonCandidate->IsLargestComponentElectron() && TruePhotonCandidate->IsConversion())
1398                                 fHistoTrueSecondaryClusConvGammaPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1399                         if(((AliAODMCParticle*)AODMCTrackArray->At(Photon->GetMother()))->GetMother() > -1){
1400                                 if(((AliAODMCParticle*)AODMCTrackArray->At(((AliAODMCParticle*)AODMCTrackArray->At(Photon->GetMother()))->GetMother()))->GetPdgCode() == 3122){
1401                                         if (TruePhotonCandidate->IsLargestComponentPhoton())
1402                                                 fHistoTrueSecondaryClusGammaFromXFromLambdaPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1403                                         if (TruePhotonCandidate->IsLargestComponentElectron() && TruePhotonCandidate->IsConversion())
1404                                                 fHistoTrueSecondaryClusConvGammaFromXFromLambdaPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1405                                 }
1406                                 if(((AliAODMCParticle*)AODMCTrackArray->At(((AliAODMCParticle*)AODMCTrackArray->At(Photon->GetMother()))->GetMother()))->GetPdgCode() == 310){
1407                                         if (TruePhotonCandidate->IsLargestComponentPhoton())
1408                                                 fHistoTrueSecondaryClusGammaFromXFromK0sPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1409                                         if (TruePhotonCandidate->IsLargestComponentElectron() && TruePhotonCandidate->IsConversion())
1410                                                 fHistoTrueSecondaryClusConvGammaFromXFromK0sPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1411                                 }
1412                                 if(((AliAODMCParticle*)AODMCTrackArray->At(((AliAODMCParticle*)AODMCTrackArray->At(Photon->GetMother()))->GetMother()))->GetPdgCode() == 221){
1413                                         if (TruePhotonCandidate->IsLargestComponentPhoton())
1414                                                 fHistoTrueSecondaryClusGammaFromXFromEtasPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1415                                         if (TruePhotonCandidate->IsLargestComponentElectron() && TruePhotonCandidate->IsConversion())
1416                                                 fHistoTrueSecondaryClusConvGammaFromXFromEtasPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1417                                                 
1418                                 }       
1419                         }
1420                 }       
1421         }
1422 }
1423
1424 //________________________________________________________________________
1425 void AliAnalysisTaskGammaCalo::ProcessAODMCParticles()
1426 {
1427         
1428         TClonesArray *AODMCTrackArray = dynamic_cast<TClonesArray*>(fInputEvent->FindListObject(AliAODMCParticle::StdBranchName()));
1429         
1430         // Loop over all primary MC particle
1431         for(Int_t i = 0; i < AODMCTrackArray->GetEntriesFast(); i++) {
1432                 
1433                 AliAODMCParticle* particle = static_cast<AliAODMCParticle*>(AODMCTrackArray->At(i));
1434                 if (!particle) continue;
1435                 if (!particle->IsPrimary()) continue;
1436                 
1437                 Int_t isMCFromMBHeader = -1;
1438                 if(((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetSignalRejection() != 0){
1439                         isMCFromMBHeader = ((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsParticleFromBGEvent(i, fMCStack, fInputEvent);
1440                         if(isMCFromMBHeader == 0 && ((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetSignalRejection() != 3) continue;
1441                 }
1442                 
1443                 if(((AliCaloPhotonCuts*)fClusterCutArray->At(fiCut))->ClusterIsSelectedAODMC(particle,AODMCTrackArray)){
1444                         fHistoMCAllGammaPt[fiCut]->Fill(particle->Pt()); // All MC Gamma
1445                         if(particle->GetMother() >-1){ // Meson Decay Gamma
1446                                 switch((static_cast<AliAODMCParticle*>(AODMCTrackArray->At(particle->GetMother())))->GetPdgCode()){
1447                                 case 111: // Pi0
1448                                         fHistoMCDecayGammaPi0Pt[fiCut]->Fill(particle->Pt());
1449                                         break;
1450                                 case 113: // Rho0
1451                                         fHistoMCDecayGammaRhoPt[fiCut]->Fill(particle->Pt());
1452                                         break;
1453                                 case 221: // Eta
1454                                         fHistoMCDecayGammaEtaPt[fiCut]->Fill(particle->Pt());
1455                                         break;
1456                                 case 223: // Omega
1457                                         fHistoMCDecayGammaOmegaPt[fiCut]->Fill(particle->Pt());
1458                                         break;
1459                                 case 331: // Eta'
1460                                         fHistoMCDecayGammaEtapPt[fiCut]->Fill(particle->Pt());
1461                                         break;
1462                                 case 333: // Phi
1463                                         fHistoMCDecayGammaPhiPt[fiCut]->Fill(particle->Pt());
1464                                         break;
1465                                 case 3212: // Sigma
1466                                         fHistoMCDecayGammaSigmaPt[fiCut]->Fill(particle->Pt());
1467                                         break;
1468                                 }
1469                         }
1470                 }
1471                         // Converted MC Gamma
1472                 if(fDoMesonAnalysis){
1473                         if(particle->GetPdgCode() == 310 && fDoMesonQA > 0){
1474                                 Double_t mesonY = 10.;
1475                                 if(particle->E() - particle->Pz() == 0 || particle->E() + particle->Pz() == 0){
1476                                         mesonY=10.-((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetEtaShift();
1477                                 } else {
1478                                         mesonY = 0.5*(TMath::Log((particle->E()+particle->Pz()) / (particle->E()-particle->Pz())))-((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetEtaShift();
1479                                 }
1480                                 Float_t weightedK0s= 1;
1481                                 if(((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsParticleFromBGEvent(i, fMCStack, fInputEvent)){
1482                                         if (particle->Pt()>0.005){
1483                                                 weightedK0s= ((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetWeightForMeson(fV0Reader->GetPeriodName(),i, 0x0, fInputEvent);
1484                                                 //cout << "MC input \t"<<i << "\t" <<  particle->Pt()<<"\t"<<weighted << endl;
1485                                         }
1486                                 }
1487                                 fHistoMCK0sPt[fiCut]->Fill(particle->Pt(),weightedK0s);
1488                                 fHistoMCK0sWOWeightPt[fiCut]->Fill(particle->Pt());
1489                                 fHistoMCK0sPtY[fiCut]->Fill(particle->Pt(),mesonY,weightedK0s);
1490                         }
1491                         if(((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))
1492                                 ->MesonIsSelectedAODMC(particle,AODMCTrackArray,((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetEtaShift())){
1493                                 AliAODMCParticle* daughter0 = static_cast<AliAODMCParticle*>(AODMCTrackArray->At(particle->GetDaughter(0)));
1494                                 AliAODMCParticle* daughter1 = static_cast<AliAODMCParticle*>(AODMCTrackArray->At(particle->GetDaughter(1)));
1495                                 Float_t weighted= 1;
1496                                 if(((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsParticleFromBGEvent(i, fMCStack, fInputEvent)){
1497                                         if (particle->Pt()>0.005){
1498                                                 weighted= ((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetWeightForMeson(fV0Reader->GetPeriodName(),i, 0x0, fInputEvent);
1499                                                 //                   if(particle->GetPdgCode() == 221){
1500                                                 //                      cout << "MC input \t"<<i << "\t" <<  particle->Pt()<<"\t"<<weighted << endl;
1501                                                 //                   }
1502                                         }
1503                                 }
1504                                 Double_t mesonY = 10.;
1505                                 if(particle->E() - particle->Pz() == 0 || particle->E() + particle->Pz() == 0){
1506                                         mesonY=10.-((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetEtaShift();
1507                                 } else{
1508                                         mesonY = 0.5*(TMath::Log((particle->E()+particle->Pz()) / (particle->E()-particle->Pz())))-((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetEtaShift();
1509                                 }
1510                                 
1511                                 if(particle->GetPdgCode() == 111){
1512                                         fHistoMCPi0Pt[fiCut]->Fill(particle->Pt(),weighted); // All MC Pi0
1513                                         fHistoMCPi0WOWeightPt[fiCut]->Fill(particle->Pt());
1514                                         if (fDoMesonQA > 0) fHistoMCPi0PtY[fiCut]->Fill(particle->Pt(),mesonY,weighted); // All MC Pi0
1515                                 } else if(particle->GetPdgCode() == 221){
1516                                         fHistoMCEtaPt[fiCut]->Fill(particle->Pt(),weighted); // All MC Eta
1517                                         fHistoMCEtaWOWeightPt[fiCut]->Fill(particle->Pt());
1518                                         if (fDoMesonQA > 0) fHistoMCEtaPtY[fiCut]->Fill(particle->Pt(),mesonY,weighted); // All MC Pi0
1519                                 }
1520                                 
1521                                 // Check the acceptance for both gammas
1522                                 if(((AliCaloPhotonCuts*)fClusterCutArray->At(fiCut))->ClusterIsSelectedAODMC(daughter0,AODMCTrackArray) &&
1523                                 ((AliCaloPhotonCuts*)fClusterCutArray->At(fiCut))->ClusterIsSelectedAODMC(daughter1,AODMCTrackArray) ){
1524                                         
1525                                         if(particle->GetPdgCode() == 111){
1526                                                 fHistoMCPi0InAccPt[fiCut]->Fill(particle->Pt(),weighted); // MC Pi0 with gamma in acc
1527                                         } else if(particle->GetPdgCode() == 221){
1528                                                 fHistoMCEtaInAccPt[fiCut]->Fill(particle->Pt(),weighted); // MC Eta with gamma in acc
1529                                         }
1530                                 }
1531                         }
1532                 }
1533         }
1534         
1535 }
1536 //________________________________________________________________________
1537 void AliAnalysisTaskGammaCalo::ProcessMCParticles()
1538 {
1539         // Loop over all primary MC particle
1540         for(Int_t i = 0; i < fMCStack->GetNprimary(); i++) {
1541                 TParticle* particle = (TParticle *)fMCStack->Particle(i);
1542                 if (!particle) continue;
1543                 
1544                 Int_t isMCFromMBHeader = -1;
1545                 if(((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetSignalRejection() != 0){
1546                         isMCFromMBHeader = ((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsParticleFromBGEvent(i, fMCStack, fInputEvent);
1547                         if(isMCFromMBHeader == 0 && ((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetSignalRejection() != 3) continue;
1548                 }
1549                 
1550                 if(((AliCaloPhotonCuts*)fClusterCutArray->At(fiCut))->ClusterIsSelectedMC(particle,fMCStack)){
1551                         fHistoMCAllGammaPt[fiCut]->Fill(particle->Pt()); // All MC Gamma                        
1552                         if(particle->GetMother(0) >-1){ // Meson Decay Gamma
1553                                 switch(fMCStack->Particle(particle->GetMother(0))->GetPdgCode()){
1554                                 case 111: // Pi0
1555                                         fHistoMCDecayGammaPi0Pt[fiCut]->Fill(particle->Pt());
1556                                         break;
1557                                 case 113: // Rho0
1558                                         fHistoMCDecayGammaRhoPt[fiCut]->Fill(particle->Pt());
1559                                         break;
1560                                 case 221: // Eta
1561                                         fHistoMCDecayGammaEtaPt[fiCut]->Fill(particle->Pt());
1562                                         break;
1563                                 case 223: // Omega
1564                                         fHistoMCDecayGammaOmegaPt[fiCut]->Fill(particle->Pt());
1565                                         break;
1566                                 case 331: // Eta'
1567                                         fHistoMCDecayGammaEtapPt[fiCut]->Fill(particle->Pt());
1568                                         break;
1569                                 case 333: // Phi
1570                                         fHistoMCDecayGammaPhiPt[fiCut]->Fill(particle->Pt());
1571                                         break;
1572                                 case 3212: // Sigma
1573                                         fHistoMCDecayGammaSigmaPt[fiCut]->Fill(particle->Pt());
1574                                         break;
1575                                 }
1576                         }
1577                 }
1578                 if(fDoMesonAnalysis){
1579                         if(particle->GetPdgCode() == 310 && fDoMesonQA > 0){
1580                                 Double_t mesonY = 10.;
1581                                 if(particle->Energy() - particle->Pz() == 0 || particle->Energy() + particle->Pz() == 0){
1582                                         mesonY=10.-((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetEtaShift();
1583                                 } else{
1584                                         mesonY = 0.5*(TMath::Log((particle->Energy()+particle->Pz()) / (particle->Energy()-particle->Pz())))-((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetEtaShift();
1585                                 }
1586                                 Float_t weightedK0s= 1;
1587                                 if(((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsParticleFromBGEvent(i, fMCStack, fInputEvent)){
1588                                         if (particle->Pt()>0.005){
1589                                                 weightedK0s= ((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetWeightForMeson(fV0Reader->GetPeriodName(),i, fMCStack, fInputEvent);
1590                                                 //cout << "MC input \t"<<i << "\t" <<  particle->Pt()<<"\t"<<weighted << endl;
1591                                         }
1592                                 }
1593                                 if (fMCStack->IsPhysicalPrimary(i)){
1594                                         fHistoMCK0sPt[fiCut]->Fill(particle->Pt(),weightedK0s);
1595                                         fHistoMCK0sWOWeightPt[fiCut]->Fill(particle->Pt());
1596                                         fHistoMCK0sPtY[fiCut]->Fill(particle->Pt(),mesonY,weightedK0s);
1597                                 }
1598                         }
1599                         if(((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))
1600                                 ->MesonIsSelectedMC(particle,fMCStack,((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetEtaShift())){
1601                                 TParticle* daughter0 = (TParticle*)fMCStack->Particle(particle->GetFirstDaughter());
1602                                 TParticle* daughter1 = (TParticle*)fMCStack->Particle(particle->GetLastDaughter());
1603                                 
1604                                 Float_t weighted= 1;
1605                                 if(((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsParticleFromBGEvent(i, fMCStack, fInputEvent)){
1606                                         if (particle->Pt()>0.005){
1607                                                 weighted= ((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetWeightForMeson(fV0Reader->GetPeriodName(),i, fMCStack, fInputEvent);
1608                                                 //                   if(particle->GetPdgCode() == 221){
1609                                                 //                      cout << "MC input \t"<<i << "\t" <<  particle->Pt()<<"\t"<<weighted << endl;
1610                                                 //                   }
1611                                         }
1612                                 }
1613                                 Double_t mesonY = 10.;
1614                                 if(particle->Energy() - particle->Pz() == 0 || particle->Energy() + particle->Pz() == 0){
1615                                         mesonY=10.-((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetEtaShift();
1616                                 } else{
1617                                         mesonY = 0.5*(TMath::Log((particle->Energy()+particle->Pz()) / (particle->Energy()-particle->Pz())))-((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetEtaShift();
1618                                 }
1619                                 
1620                                 if(particle->GetPdgCode() == 111){
1621                                         fHistoMCPi0Pt[fiCut]->Fill(particle->Pt(),weighted); // All MC Pi0
1622                                         fHistoMCPi0WOWeightPt[fiCut]->Fill(particle->Pt());
1623                                         if (fDoMesonQA > 0) fHistoMCPi0PtY[fiCut]->Fill(particle->Pt(),mesonY,weighted); // All MC Pi0
1624                                 } else if(particle->GetPdgCode() == 221){
1625                                         fHistoMCEtaPt[fiCut]->Fill(particle->Pt(),weighted); // All MC Eta
1626                                         fHistoMCEtaWOWeightPt[fiCut]->Fill(particle->Pt());
1627                                         if (fDoMesonQA > 0) fHistoMCEtaPtY[fiCut]->Fill(particle->Pt(),mesonY,weighted); // All MC Pi0
1628                                 }
1629                                 
1630                                 // Check the acceptance for both gammas
1631                                 if(((AliCaloPhotonCuts*)fClusterCutArray->At(fiCut))->ClusterIsSelectedMC(daughter0,fMCStack) &&
1632                                 ((AliCaloPhotonCuts*)fClusterCutArray->At(fiCut))->ClusterIsSelectedMC(daughter1,fMCStack) ){                                   
1633                                         if(particle->GetPdgCode() == 111){
1634                                                 fHistoMCPi0InAccPt[fiCut]->Fill(particle->Pt(),weighted); // MC Pi0 with gamma in acc
1635                                         } else if(particle->GetPdgCode() == 221){
1636                                                 fHistoMCEtaInAccPt[fiCut]->Fill(particle->Pt(),weighted); // MC Eta with gamma in acc
1637                                         }
1638                                 }
1639                         }
1640                 }
1641         }
1642   
1643         if (fDoMesonQA){
1644                 for(Int_t i = fMCStack->GetNprimary(); i < fMCStack->GetNtrack(); i++) {
1645                         TParticle* particle = (TParticle *)fMCStack->Particle(i);
1646                         if (!particle) continue;
1647       
1648                         Int_t isMCFromMBHeader = -1;
1649                         if(((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetSignalRejection() != 0){
1650                                 isMCFromMBHeader = ((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsParticleFromBGEvent(i, fMCStack, fInputEvent);
1651                                 if(isMCFromMBHeader == 0 && ((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetSignalRejection() != 3) continue;
1652                         }
1653       
1654                         if(fDoMesonAnalysis){
1655                                 if(((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->MesonIsSelectedMC(particle,fMCStack,((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetEtaShift())){
1656                                         Float_t weighted= 1;
1657                                         if(((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsParticleFromBGEvent(i, fMCStack, fInputEvent)){
1658                                                 if (particle->Pt()>0.005){
1659                                                         weighted= ((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetWeightForMeson(fV0Reader->GetPeriodName(),i, fMCStack, fInputEvent);
1660               //                   if(particle->GetPdgCode() == 221){
1661               //                      cout << "MC input \t"<<i << "\t" <<  particle->Pt()<<"\t"<<weighted << endl;
1662               //                   }
1663                                                 }
1664                                         }
1665                                         
1666                                         if(particle->GetPdgCode() == 111){
1667                                                 Int_t pdgCode = ((TParticle*)fMCStack->Particle( particle->GetFirstMother() ))->GetPdgCode();
1668                                                 Int_t source = GetSourceClassification(111,pdgCode);
1669                                                 fHistoMCSecPi0PtvsSource[fiCut]->Fill(particle->Pt(),source,weighted); // All MC Pi0
1670                                                 fHistoMCSecPi0Source[fiCut]->Fill(pdgCode);
1671                                         } else if(particle->GetPdgCode() == 221){
1672                                                 Int_t pdgCode = ((TParticle*)fMCStack->Particle( particle->GetFirstMother() ))->GetPdgCode();
1673                                                 fHistoMCSecEtaPt[fiCut]->Fill(particle->Pt(),weighted); // All MC Pi0
1674                                                 fHistoMCSecEtaSource[fiCut]->Fill(pdgCode);
1675                                         }
1676                                 }
1677                         }
1678                 }
1679         }
1680 }
1681
1682 //________________________________________________________________________
1683 void AliAnalysisTaskGammaCalo::CalculatePi0Candidates(){
1684         
1685         // Conversion Gammas
1686         if(fClusterCandidates->GetEntries()>0){
1687
1688                 // vertex
1689                 Double_t vertex[3] = {0};
1690                 InputEvent()->GetPrimaryVertex()->GetXYZ(vertex);
1691
1692                 for(Int_t firstGammaIndex=0;firstGammaIndex<fClusterCandidates->GetEntries();firstGammaIndex++){
1693                         AliAODConversionPhoton *gamma0=dynamic_cast<AliAODConversionPhoton*>(fClusterCandidates->At(firstGammaIndex));
1694                         if (gamma0==NULL) continue;
1695                         
1696                         for(Int_t secondGammaIndex=0;secondGammaIndex<fClusterCandidates->GetEntries();secondGammaIndex++){
1697                                 if (firstGammaIndex == secondGammaIndex) continue;
1698                                 AliAODConversionPhoton *gamma1=dynamic_cast<AliAODConversionPhoton*>(fClusterCandidates->At(secondGammaIndex));
1699                                 if (gamma1==NULL) continue;
1700                                 
1701                                 AliAODConversionMother *pi0cand = new AliAODConversionMother(gamma0,gamma1);
1702                                 pi0cand->SetLabels(firstGammaIndex,secondGammaIndex);
1703                                 
1704                                 
1705                                 
1706                                 if((((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->MesonIsSelected(pi0cand,kTRUE,((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetEtaShift()))){
1707                                         fHistoMotherInvMassPt[fiCut]->Fill(pi0cand->M(),pi0cand->Pt());
1708                                         // fill new histograms
1709                                         if(pi0cand->GetAlpha()<0.1)
1710                                                 fHistoMotherInvMassEalpha[fiCut]->Fill(pi0cand->M(),pi0cand->E());
1711                                         
1712                                         if (fDoMesonQA > 0){
1713                                                 if ( pi0cand->M() > 0.05 && pi0cand->M() < 0.17){
1714                                                         fHistoMotherPi0PtY[fiCut]->Fill(pi0cand->Pt(),pi0cand->Rapidity()-((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetEtaShift());
1715                                                         fHistoMotherPi0PtAlpha[fiCut]->Fill(pi0cand->Pt(),pi0cand->GetAlpha());
1716                                                         fHistoMotherPi0PtOpenAngle[fiCut]->Fill(pi0cand->Pt(),pi0cand->GetOpeningAngle());      
1717                                                 }
1718                                                 if ( pi0cand->M() > 0.45 && pi0cand->M() < 0.65){
1719                                                         fHistoMotherEtaPtY[fiCut]->Fill(pi0cand->Pt(),pi0cand->Rapidity()-((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetEtaShift());
1720                                                         fHistoMotherEtaPtAlpha[fiCut]->Fill(pi0cand->Pt(),pi0cand->GetAlpha());
1721                                                         fHistoMotherEtaPtOpenAngle[fiCut]->Fill(pi0cand->Pt(),pi0cand->GetOpeningAngle());
1722                                                 }
1723                                         }
1724                                         if(((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->DoBGCalculation()){
1725                                                 Int_t zbin = 0;
1726                                                 Int_t mbin = 0;
1727                                                 
1728                                                 if(((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->BackgroundHandlerType() == 0){
1729                                                         zbin = fBGHandler[fiCut]->GetZBinIndex(fInputEvent->GetPrimaryVertex()->GetZ());
1730                                                         if(((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->UseTrackMultiplicity()){
1731                                                                 mbin = fBGHandler[fiCut]->GetMultiplicityBinIndex(fV0Reader->GetNumberOfPrimaryTracks());
1732                                                         } else {
1733                                                                 mbin = fBGHandler[fiCut]->GetMultiplicityBinIndex(fClusterCandidates->GetEntries());
1734                                                         }
1735                                                 } 
1736                                                 Double_t sparesFill[4] = {pi0cand->M(),pi0cand->Pt(),(Double_t)zbin,(Double_t)mbin};
1737                                                 fSparseMotherInvMassPtZM[fiCut]->Fill(sparesFill,1);
1738                                         }
1739                                 
1740                                         if(fIsMC){
1741                                                 if(fInputEvent->IsA()==AliESDEvent::Class())
1742                                                         ProcessTrueMesonCandidates(pi0cand,gamma0,gamma1);
1743                                                 if(fInputEvent->IsA()==AliAODEvent::Class())
1744                                                         ProcessTrueMesonCandidatesAOD(pi0cand,gamma0,gamma1);
1745                                         }
1746                                         
1747                                         if (fDoMesonQA == 1){
1748                                                 fHistoMotherInvMassECalib[fiCut]->Fill(pi0cand->M(),gamma1->E());
1749                                                 if(pi0cand->GetAlpha()<0.1)
1750                                                 fHistoMotherInvMassECalibalpha[fiCut]->Fill(pi0cand->M(),gamma1->E());            
1751                                         }
1752                                         
1753                                 }
1754                                 for (Int_t thirdGammaIndex=0;thirdGammaIndex<fClusterCandidates->GetEntries();thirdGammaIndex++){
1755                                         if (firstGammaIndex == thirdGammaIndex || secondGammaIndex == thirdGammaIndex ) continue;
1756                                         AliAODConversionPhoton *gamma2=dynamic_cast<AliAODConversionPhoton*>(fClusterCandidates->At(thirdGammaIndex));
1757                                         if (gamma2==NULL) continue;
1758                                         
1759                                         AliAODConversionMother *pi0cand2 = new AliAODConversionMother(pi0cand,gamma2);
1760                                         if((((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->MesonIsSelected(pi0cand2,kTRUE,((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetEtaShift()))){
1761                                                 fHistoMotherInvMass3ClusterPt[fiCut]->Fill(pi0cand2->M(),pi0cand2->Pt());
1762                                         }
1763                                         delete pi0cand2;
1764                                         pi0cand2=0x0;
1765                                         
1766                                 }
1767                                 
1768                                 delete pi0cand;
1769                                 pi0cand=0x0;
1770                         }
1771                 }
1772         }
1773 }
1774 //______________________________________________________________________
1775 void AliAnalysisTaskGammaCalo::ProcessTrueMesonCandidates(AliAODConversionMother *Pi0Candidate, AliAODConversionPhoton *TrueGammaCandidate0, AliAODConversionPhoton *TrueGammaCandidate1)
1776 {
1777         // Process True Mesons
1778         AliStack *MCStack = fMCEvent->Stack();
1779         
1780         Bool_t isTruePi0 = kFALSE;
1781         Bool_t isTrueEta = kFALSE;
1782         Int_t gamma0MCLabel = TrueGammaCandidate0->GetCaloPhotonMCLabel(0);     // get most probable MC label
1783         Int_t gamma0MotherLabel = -1;
1784
1785         if(gamma0MCLabel != -1){ // Gamma is Combinatorial; MC Particles don't belong to the same Mother
1786                 TParticle * gammaMC0 = (TParticle*)MCStack->Particle(gamma0MCLabel);
1787                 if (TrueGammaCandidate0->IsLargestComponentPhoton() || TrueGammaCandidate0->IsLargestComponentElectron()){              // largest component is electro magnetic
1788                         // get mother of interest (pi0 or eta)
1789                         if (TrueGammaCandidate0->IsLargestComponentPhoton()){                                                                                                           // for photons its the direct mother 
1790                                 gamma0MotherLabel=gammaMC0->GetMother(0);
1791                         } else if (TrueGammaCandidate0->IsLargestComponentElectron()){                                                                                          // for electrons its either the direct mother or for conversions the grandmother
1792                                 if (TrueGammaCandidate0->IsConversion()) gamma0MotherLabel=MCStack->Particle(gammaMC0->GetMother(0))->GetMother(0);
1793                                 else gamma0MotherLabel=gammaMC0->GetMother(0); 
1794                         }
1795                 }       
1796         }
1797         if (!TrueGammaCandidate1->GetIsCaloPhoton()) AliFatal("CaloPhotonFlag has not been set. Aborting");
1798         
1799         Int_t gamma1MCLabel = TrueGammaCandidate1->GetCaloPhotonMCLabel(0);     // get most probable MC label
1800         Int_t gamma1MotherLabel = -1;
1801         // check if 
1802
1803         if(gamma1MCLabel != -1){ // Gamma is Combinatorial; MC Particles don't belong to the same Mother
1804                 // Daughters Gamma 1
1805                 TParticle * gammaMC1 = (TParticle*)MCStack->Particle(gamma1MCLabel);
1806                 if (TrueGammaCandidate1->IsLargestComponentPhoton() || TrueGammaCandidate1->IsLargestComponentElectron()){              // largest component is electro magnetic
1807                         // get mother of interest (pi0 or eta)
1808                         if (TrueGammaCandidate1->IsLargestComponentPhoton()){                                                                                                           // for photons its the direct mother 
1809                                 gamma1MotherLabel=gammaMC1->GetMother(0);
1810                         } else if (TrueGammaCandidate1->IsLargestComponentElectron()){                                                                                          // for electrons its either the direct mother or for conversions the grandmother
1811                                 if (TrueGammaCandidate1->IsConversion()) gamma1MotherLabel=MCStack->Particle(gammaMC1->GetMother(0))->GetMother(0);
1812                                 else gamma1MotherLabel=gammaMC1->GetMother(0); 
1813                         }
1814                 }       
1815         }
1816                         
1817         if(gamma0MotherLabel>=0 && gamma0MotherLabel==gamma1MotherLabel){
1818                 if(((TParticle*)MCStack->Particle(gamma1MotherLabel))->GetPdgCode() == 111){
1819                         isTruePi0=kTRUE;
1820                 }
1821                 if(((TParticle*)MCStack->Particle(gamma1MotherLabel))->GetPdgCode() == 221){
1822                         isTrueEta=kTRUE;
1823                 }
1824         }
1825         
1826         if(isTruePi0 || isTrueEta){// True Pion or Eta
1827                 if (isTruePi0)  fHistoTruePi0InvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt());
1828                 if (isTrueEta)  fHistoTrueEtaInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt());
1829                 if (fDoMesonQA > 0){
1830                         // both gammas are real gammas
1831                         if (TrueGammaCandidate0->IsLargestComponentPhoton() && TrueGammaCandidate1->IsLargestComponentPhoton()) {
1832                                 if (isTruePi0) fHistoTruePi0CaloPhotonInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt());
1833                                 if (isTrueEta) fHistoTrueEtaCaloPhotonInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt());
1834                         }       
1835                         // both particles are electrons
1836                         if (TrueGammaCandidate0->IsLargestComponentElectron() && TrueGammaCandidate1->IsLargestComponentElectron() ) {
1837                                 if (isTruePi0) fHistoTruePi0CaloElectronInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt());
1838                                 if (isTrueEta) fHistoTrueEtaCaloElectronInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt());
1839                         }       
1840                         // both particles are converted electrons
1841                         if ((TrueGammaCandidate0->IsLargestComponentElectron() && TrueGammaCandidate0->IsConversion()) && (TrueGammaCandidate1->IsLargestComponentElectron() && TrueGammaCandidate1->IsConversion()) ){
1842                                 if (isTruePi0 )fHistoTruePi0CaloConvertedPhotonInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt());
1843                                 if (isTrueEta )fHistoTrueEtaCaloConvertedPhotonInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt());
1844                         }
1845                         // 1 gamma is converted the other one is normals
1846                         if ( (TrueGammaCandidate0->IsLargestComponentPhoton() && TrueGammaCandidate1->IsLargestComponentElectron() && TrueGammaCandidate1->IsConversion()) ||
1847                                  (TrueGammaCandidate1->IsLargestComponentPhoton() && TrueGammaCandidate0->IsLargestComponentElectron() && TrueGammaCandidate0->IsConversion())
1848                         ) {
1849                                 if (isTruePi0) fHistoTruePi0CaloMixedPhotonConvPhotonInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt());
1850                                 if (isTrueEta) fHistoTrueEtaCaloMixedPhotonConvPhotonInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt());
1851                         }
1852                         // at least one of the photon is merged
1853                         if (TrueGammaCandidate0->IsMerged() || TrueGammaCandidate0->IsMergedPartConv() || TrueGammaCandidate0->IsDalitzMerged() || TrueGammaCandidate1->IsMerged() || TrueGammaCandidate1->IsMergedPartConv() || TrueGammaCandidate1->IsDalitzMerged() ){
1854                                 if (isTruePi0) fHistoTruePi0CaloMergedClusterInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt());
1855                                 if (isTrueEta) fHistoTruePi0CaloMergedClusterInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt());
1856                         }       
1857                         // at least one of the photon is merged and part conv
1858                         if (TrueGammaCandidate1->IsMergedPartConv() || TrueGammaCandidate0->IsMergedPartConv()) {       
1859                                 if (isTruePi0) fHistoTruePi0CaloMergedClusterPartConvInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt());
1860                                 if (isTrueEta) fHistoTrueEtaCaloMergedClusterPartConvInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt());
1861                         }       
1862                 }
1863         
1864                 if (fDoMesonQA > 0){
1865                         if (isTruePi0){
1866                                 if ( Pi0Candidate->M() > 0.05 && Pi0Candidate->M() < 0.17){
1867                                         fHistoTruePi0PtY[fiCut]->Fill(Pi0Candidate->Pt(),Pi0Candidate->Rapidity()-((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetEtaShift());
1868                                         fHistoTruePi0PtAlpha[fiCut]->Fill(Pi0Candidate->Pt(),Pi0Candidate->GetAlpha());
1869                                         fHistoTruePi0PtOpenAngle[fiCut]->Fill(Pi0Candidate->Pt(),Pi0Candidate->GetOpeningAngle());
1870                                 }
1871                         } else if (isTrueEta){
1872                                 if ( Pi0Candidate->M() > 0.45 && Pi0Candidate->M() < 0.65){
1873                                         fHistoTrueEtaPtY[fiCut]->Fill(Pi0Candidate->Pt(),Pi0Candidate->Rapidity()-((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetEtaShift());
1874                                         fHistoTrueEtaPtAlpha[fiCut]->Fill(Pi0Candidate->Pt(),Pi0Candidate->GetAlpha());
1875                                         fHistoTrueEtaPtOpenAngle[fiCut]->Fill(Pi0Candidate->Pt(),Pi0Candidate->GetOpeningAngle());
1876                                 }
1877                         }
1878                 }
1879                 
1880                 if(gamma0MotherLabel >= MCStack->GetNprimary()){ // Secondary Meson
1881                         Int_t secMotherLabel = ((TParticle*)MCStack->Particle(gamma0MotherLabel))->GetMother(0);
1882                         Float_t weightedSec= 1;
1883                         if(((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsParticleFromBGEvent(secMotherLabel, fMCStack, fInputEvent) && MCStack->Particle(secMotherLabel)->GetPdgCode()==310){
1884                                 weightedSec= ((AliConvEventCuts*)fEventCutArray->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
1885                                 //cout << "MC input \t"<<i << "\t" <<  particle->Pt()<<"\t"<<weighted << endl;
1886                         }
1887                         if (isTruePi0) fHistoTrueSecondaryPi0InvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt(),weightedSec);
1888                         if (isTrueEta) fHistoTrueSecondaryEtaInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt(),weightedSec);
1889                         if (secMotherLabel >-1){
1890                                 if(MCStack->Particle(secMotherLabel)->GetPdgCode()==310 && isTruePi0){
1891                                         fHistoTrueSecondaryPi0FromK0sInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt(),weightedSec);
1892                                         if (fDoMesonQA > 0)fHistoTrueK0sWithPi0DaughterMCPt[fiCut]->Fill(MCStack->Particle(secMotherLabel)->Pt());
1893                                 }
1894                                 if(MCStack->Particle(secMotherLabel)->GetPdgCode()==221 && isTruePi0){
1895                                         fHistoTrueSecondaryPi0FromEtaInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt(),weightedSec);
1896                                         if (fDoMesonQA > 0)fHistoTrueEtaWithPi0DaughterMCPt[fiCut]->Fill(MCStack->Particle(secMotherLabel)->Pt());
1897                                 }
1898                                 if(MCStack->Particle(secMotherLabel)->GetPdgCode()==3122 && isTruePi0){
1899                                         fHistoTrueSecondaryPi0FromLambdaInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt(),weightedSec);
1900                                         if (fDoMesonQA > 0)fHistoTrueLambdaWithPi0DaughterMCPt[fiCut]->Fill(MCStack->Particle(secMotherLabel)->Pt());
1901                                 }
1902                         }
1903                 } else { // Only primary pi0 for efficiency calculation
1904                         Float_t weighted= 1;
1905                         if(((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsParticleFromBGEvent(gamma1MotherLabel, fMCStack, fInputEvent)){
1906                                 if (((TParticle*)MCStack->Particle(gamma1MotherLabel))->Pt()>0.005){
1907                                         weighted= ((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetWeightForMeson(fV0Reader->GetPeriodName(),gamma1MotherLabel, fMCStack, fInputEvent);
1908                                         //                      cout << "rec \t " <<gamma1MotherLabel << "\t" <<  weighted << endl;
1909                                 }
1910                         }
1911                         if (isTruePi0){
1912                                 fHistoTruePrimaryPi0InvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt(),weighted);
1913                                 fHistoTruePrimaryPi0W0WeightingInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt());
1914                                 fProfileTruePrimaryPi0WeightsInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt(),weighted);
1915                         } else if (isTrueEta){
1916                                 fHistoTruePrimaryEtaInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt(),weighted);
1917                                 fHistoTruePrimaryEtaW0WeightingInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt());
1918                                 fProfileTruePrimaryEtaWeightsInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt(),weighted);
1919                         }       
1920                                 
1921                         if (fDoMesonQA > 0){
1922                                 if(isTruePi0){ // Only primary pi0 for resolution
1923                                         fHistoTruePrimaryPi0MCPtResolPt[fiCut]->Fill(((TParticle*)MCStack->Particle(gamma1MotherLabel))->Pt(),(Pi0Candidate->Pt()-((TParticle*)MCStack->Particle(gamma1MotherLabel))->Pt())/((TParticle*)MCStack->Particle(gamma1MotherLabel))->Pt(),weighted);
1924                                 }
1925                                 if (isTrueEta){ // Only primary eta for resolution
1926                                         fHistoTruePrimaryEtaMCPtResolPt[fiCut]->Fill(((TParticle*)MCStack->Particle(gamma1MotherLabel))->Pt(),(Pi0Candidate->Pt()-((TParticle*)MCStack->Particle(gamma1MotherLabel))->Pt())/((TParticle*)MCStack->Particle(gamma1MotherLabel))->Pt(),weighted);
1927                                 }
1928                         }
1929                 }       
1930         } else if(!isTruePi0 && !isTrueEta){ // Background
1931                 if (fDoMesonQA > 0){
1932                         if(gamma0MotherLabel>-1 && gamma1MotherLabel>-1){ // Both Tracks are Photons and have a mother but not Pi0 or Eta
1933                                 fHistoTrueBckGGInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt());
1934                         } else { // No photon or without mother
1935                                 fHistoTrueBckContInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt());
1936                         }
1937                 }
1938         }
1939
1940 }
1941 //______________________________________________________________________
1942 void AliAnalysisTaskGammaCalo::ProcessTrueMesonCandidatesAOD(AliAODConversionMother *Pi0Candidate, AliAODConversionPhoton *TrueGammaCandidate0, AliAODConversionPhoton *TrueGammaCandidate1)
1943 {
1944         
1945         // Process True Mesons
1946         TClonesArray *AODMCTrackArray = dynamic_cast<TClonesArray*>(fInputEvent->FindListObject(AliAODMCParticle::StdBranchName()));
1947         Bool_t isTruePi0 = kFALSE;
1948         Bool_t isTrueEta = kFALSE;
1949                 
1950         Int_t gamma0MCLabel = TrueGammaCandidate1->GetCaloPhotonMCLabel(0);     // get most probable MC label
1951         Int_t gamma0MotherLabel = -1;
1952                 // check if 
1953
1954         if(gamma0MCLabel != -1){ // Gamma is Combinatorial; MC Particles don't belong to the same Mother
1955                 // Daughters Gamma 0
1956                 AliAODMCParticle * gammaMC0 = static_cast<AliAODMCParticle*>(AODMCTrackArray->At(gamma0MCLabel));
1957                 if (TrueGammaCandidate0->IsLargestComponentPhoton() || TrueGammaCandidate0->IsLargestComponentElectron()){              // largest component is electro magnetic
1958                         // get mother of interest (pi0 or eta)
1959                         if (TrueGammaCandidate0->IsLargestComponentPhoton()){                                                                                                           // for photons its the direct mother 
1960                                 gamma0MotherLabel=gammaMC0->GetMother();
1961                         } else if (TrueGammaCandidate0->IsLargestComponentElectron()){                                                                                          // for electrons its either the direct mother or for conversions the grandmother
1962                                 if (TrueGammaCandidate0->IsConversion()){
1963                                         AliAODMCParticle * gammaGrandMotherMC0 =  static_cast<AliAODMCParticle*>(AODMCTrackArray->At(gammaMC0->GetMother()));
1964                                         gamma0MotherLabel=gammaGrandMotherMC0->GetMother();
1965                                 } else gamma0MotherLabel=gammaMC0->GetMother(); 
1966                         }
1967                 }       
1968         }
1969
1970         Int_t gamma1MCLabel = TrueGammaCandidate1->GetCaloPhotonMCLabel(0);     // get most probable MC label
1971         Int_t gamma1MotherLabel = -1;
1972                 // check if 
1973
1974         if(gamma1MCLabel != -1){ // Gamma is Combinatorial; MC Particles don't belong to the same Mother
1975                 // Daughters Gamma 1
1976                 AliAODMCParticle * gammaMC1 = static_cast<AliAODMCParticle*>(AODMCTrackArray->At(gamma1MCLabel));
1977                 if (TrueGammaCandidate1->IsLargestComponentPhoton() || TrueGammaCandidate1->IsLargestComponentElectron()){              // largest component is electro magnetic
1978                         // get mother of interest (pi0 or eta)
1979                         if (TrueGammaCandidate1->IsLargestComponentPhoton()){                                                                                                           // for photons its the direct mother 
1980                                 gamma1MotherLabel=gammaMC1->GetMother();
1981                         } else if (TrueGammaCandidate1->IsLargestComponentElectron()){                                                                                          // for electrons its either the direct mother or for conversions the grandmother
1982                                 if (TrueGammaCandidate1->IsConversion()){
1983                                         AliAODMCParticle * gammaGrandMotherMC1 =  static_cast<AliAODMCParticle*>(AODMCTrackArray->At(gammaMC1->GetMother()));
1984                                         gamma1MotherLabel=gammaGrandMotherMC1->GetMother();
1985                                 } else gamma1MotherLabel=gammaMC1->GetMother(); 
1986                         }
1987                 }       
1988         }
1989                         
1990         if(gamma0MotherLabel>=0 && gamma0MotherLabel==gamma1MotherLabel){
1991                 if(((AliAODMCParticle*)AODMCTrackArray->At(gamma1MotherLabel))->GetPdgCode() == 111){
1992                         isTruePi0=kTRUE;
1993                 }
1994                 if(((AliAODMCParticle*)AODMCTrackArray->At(gamma1MotherLabel))->GetPdgCode() == 221){
1995                         isTrueEta=kTRUE;
1996                 }
1997         }
1998         
1999         if(isTruePi0 || isTrueEta){// True Pion or Eta
2000                 if (isTruePi0)fHistoTruePi0InvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt());
2001                 if (isTrueEta)fHistoTrueEtaInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt());
2002                 if (fDoMesonQA > 0){
2003                         // both gammas are real gammas
2004                         if (TrueGammaCandidate0->IsLargestComponentPhoton() && TrueGammaCandidate1->IsLargestComponentPhoton()) {
2005                                 if (isTruePi0)fHistoTruePi0CaloPhotonInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt());
2006                                 if (isTrueEta)fHistoTrueEtaCaloPhotonInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt());
2007                         }       
2008                         // both particles are electrons
2009                         if (TrueGammaCandidate0->IsLargestComponentElectron() && TrueGammaCandidate1->IsLargestComponentElectron() ) {
2010                                 if (isTruePi0) fHistoTruePi0CaloElectronInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt());
2011                                 if (isTrueEta) fHistoTrueEtaCaloElectronInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt());
2012                         }       
2013                         // both particles are converted electrons
2014                         if ((TrueGammaCandidate0->IsLargestComponentElectron() && TrueGammaCandidate0->IsConversion()) && (TrueGammaCandidate1->IsLargestComponentElectron() && TrueGammaCandidate1->IsConversion()) ){
2015                                 if (isTruePi0 )fHistoTruePi0CaloConvertedPhotonInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt());
2016                                 if (isTrueEta )fHistoTrueEtaCaloConvertedPhotonInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt());
2017                         }
2018                         // 1 gamma is converted the other one is normals
2019                         if ( (TrueGammaCandidate0->IsLargestComponentPhoton() && TrueGammaCandidate1->IsLargestComponentElectron() && TrueGammaCandidate1->IsConversion()) ||
2020                                  (TrueGammaCandidate1->IsLargestComponentPhoton() && TrueGammaCandidate0->IsLargestComponentElectron() && TrueGammaCandidate0->IsConversion())
2021                         ) {
2022                                 if (isTruePi0) fHistoTruePi0CaloMixedPhotonConvPhotonInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt());
2023                                 if (isTrueEta) fHistoTrueEtaCaloMixedPhotonConvPhotonInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt());
2024                         }
2025                         // at least one of the photon is merged
2026                         if (TrueGammaCandidate0->IsMerged() || TrueGammaCandidate0->IsMergedPartConv() || TrueGammaCandidate0->IsDalitzMerged() || TrueGammaCandidate1->IsMerged() || TrueGammaCandidate1->IsMergedPartConv() || TrueGammaCandidate1->IsDalitzMerged() ){
2027                                 if (isTruePi0) fHistoTruePi0CaloMergedClusterInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt());
2028                                 if (isTrueEta) fHistoTrueEtaCaloMergedClusterInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt());
2029                         }
2030                         // at least one of the photon is merged and part conv
2031                         if (TrueGammaCandidate1->IsMergedPartConv() || TrueGammaCandidate0->IsMergedPartConv()) {
2032                                 if (isTruePi0)fHistoTruePi0CaloMergedClusterPartConvInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt());
2033                                 if (isTrueEta)fHistoTrueEtaCaloMergedClusterPartConvInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt());
2034                         }       
2035                 }
2036
2037                 if (fDoMesonQA > 0){
2038                         if (isTruePi0){
2039                                 if ( Pi0Candidate->M() > 0.05 && Pi0Candidate->M() < 0.17){
2040                                 fHistoTruePi0PtY[fiCut]->Fill(Pi0Candidate->Pt(),Pi0Candidate->Rapidity()-((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetEtaShift());
2041                                 fHistoTruePi0PtAlpha[fiCut]->Fill(Pi0Candidate->Pt(),Pi0Candidate->GetAlpha());
2042                                 fHistoTruePi0PtOpenAngle[fiCut]->Fill(Pi0Candidate->Pt(),Pi0Candidate->GetOpeningAngle());
2043                                 }
2044                         } else if (isTrueEta){
2045                                 if ( Pi0Candidate->M() > 0.45 && Pi0Candidate->M() < 0.65){
2046                                 fHistoTrueEtaPtY[fiCut]->Fill(Pi0Candidate->Pt(),Pi0Candidate->Rapidity()-((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetEtaShift());
2047                                 fHistoTrueEtaPtAlpha[fiCut]->Fill(Pi0Candidate->Pt(),Pi0Candidate->GetAlpha());
2048                                 fHistoTrueEtaPtOpenAngle[fiCut]->Fill(Pi0Candidate->Pt(),Pi0Candidate->GetOpeningAngle());
2049                                 }
2050                         }
2051                 }
2052                 if(!(static_cast<AliAODMCParticle*>(AODMCTrackArray->At(gamma0MotherLabel))->IsPrimary())){ // Secondary Meson
2053                         Int_t secMotherLabel = static_cast<AliAODMCParticle*>(AODMCTrackArray->At(gamma1MotherLabel))->GetMother();
2054                         Float_t weightedSec= 1;
2055                         if(((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsParticleFromBGEvent(secMotherLabel, 0x0, fInputEvent) && static_cast<AliAODMCParticle*>(AODMCTrackArray->At(secMotherLabel))->GetPdgCode()==310){
2056                                 weightedSec= ((AliConvEventCuts*)fEventCutArray->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
2057                                 //cout << "MC input \t"<<i << "\t" <<  particle->Pt()<<"\t"<<weighted << endl;
2058                         }
2059                         if (isTruePi0) fHistoTrueSecondaryPi0InvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt(),weightedSec);
2060                         if (isTrueEta) fHistoTrueSecondaryEtaInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt(),weightedSec);
2061                         if (secMotherLabel >-1){
2062                                 if(static_cast<AliAODMCParticle*>(AODMCTrackArray->At(secMotherLabel))->GetPdgCode()==310 && isTruePi0 ){
2063                                         fHistoTrueSecondaryPi0FromK0sInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt(),weightedSec);
2064                                         if (fDoMesonQA > 0)fHistoTrueK0sWithPi0DaughterMCPt[fiCut]->Fill(static_cast<AliAODMCParticle*>(AODMCTrackArray->At(secMotherLabel))->Pt());
2065                                 }
2066                                 if(static_cast<AliAODMCParticle*>(AODMCTrackArray->At(secMotherLabel))->GetPdgCode()==221 && isTruePi0){
2067                                         fHistoTrueSecondaryPi0FromEtaInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt(),weightedSec);
2068                                         if (fDoMesonQA > 0)fHistoTrueEtaWithPi0DaughterMCPt[fiCut]->Fill(static_cast<AliAODMCParticle*>(AODMCTrackArray->At(secMotherLabel))->Pt());
2069                                 }
2070                                 if(static_cast<AliAODMCParticle*>(AODMCTrackArray->At(secMotherLabel))->GetPdgCode()==3122 && isTruePi0){
2071                                         fHistoTrueSecondaryPi0FromLambdaInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt(),weightedSec);
2072                                         if (fDoMesonQA > 0)fHistoTrueLambdaWithPi0DaughterMCPt[fiCut]->Fill(static_cast<AliAODMCParticle*>(AODMCTrackArray->At(secMotherLabel))->Pt());
2073                                 }
2074                         }       
2075                 } else{ // Only primary pi0 for efficiency calculation
2076                         Float_t weighted= 1;
2077                         if(((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsParticleFromBGEvent(gamma1MotherLabel, 0x0, fInputEvent)){
2078                                 if (static_cast<AliAODMCParticle*>(AODMCTrackArray->At(gamma1MotherLabel))->Pt()>0.005){
2079                                 weighted= ((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetWeightForMeson(fV0Reader->GetPeriodName(),gamma1MotherLabel, 0x0, fInputEvent);
2080                                 //                      cout << "rec \t " <<gamma1MotherLabel << "\t" <<  weighted << endl;
2081                                 }
2082                         }
2083                         if (isTruePi0){
2084                                 fHistoTruePrimaryPi0InvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt(),weighted);
2085                                 fHistoTruePrimaryPi0W0WeightingInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt());
2086                                 fProfileTruePrimaryPi0WeightsInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt(),weighted);
2087                         } else if (isTrueEta){
2088                                 fHistoTruePrimaryEtaInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt(),weighted);
2089                                 fHistoTruePrimaryEtaW0WeightingInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt());
2090                                 fProfileTruePrimaryEtaWeightsInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt(),weighted);     
2091                         }       
2092                         if (fDoMesonQA > 0){
2093                                 if(isTruePi0){ // Only primary pi0 for resolution
2094                                         fHistoTruePrimaryPi0MCPtResolPt[fiCut]->Fill(static_cast<AliAODMCParticle*>(AODMCTrackArray->At(gamma1MotherLabel))->Pt(),
2095                                                                                                                         (Pi0Candidate->Pt()-static_cast<AliAODMCParticle*>(AODMCTrackArray->At(gamma1MotherLabel))->Pt())/static_cast<AliAODMCParticle*>(AODMCTrackArray->At(gamma1MotherLabel))->Pt(),weighted);
2096                                 
2097                                 }
2098                                 if (isTrueEta){ // Only primary eta for resolution
2099                                         fHistoTruePrimaryEtaMCPtResolPt[fiCut]->Fill(static_cast<AliAODMCParticle*>(AODMCTrackArray->At(gamma1MotherLabel))->Pt(),
2100                                                                                                                         (Pi0Candidate->Pt()-static_cast<AliAODMCParticle*>(AODMCTrackArray->At(gamma1MotherLabel))->Pt())/static_cast<AliAODMCParticle*>(AODMCTrackArray->At(gamma1MotherLabel))->Pt(),weighted);
2101                                 }
2102                         }
2103                 }
2104         } else if(!isTruePi0 && !isTrueEta) { // Background
2105                 if (fDoMesonQA > 0){
2106                         if(gamma0MotherLabel>-1 && gamma1MotherLabel>-1){ // Both Tracks are Photons and have a mother but not Pi0 or Eta
2107                                 fHistoTrueBckGGInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt());
2108                         } else { // No photon or without mother
2109                                 fHistoTrueBckContInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt());
2110                         }
2111                 }
2112         }
2113 }
2114
2115 //________________________________________________________________________
2116 void AliAnalysisTaskGammaCalo::CalculateBackground(){
2117         
2118         Int_t zbin= fBGHandler[fiCut]->GetZBinIndex(fInputEvent->GetPrimaryVertex()->GetZ());
2119         Int_t mbin = 0;
2120         
2121         if(((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->UseTrackMultiplicity()){
2122                 mbin = fBGHandler[fiCut]->GetMultiplicityBinIndex(fV0Reader->GetNumberOfPrimaryTracks());
2123         } else {
2124                 mbin = fBGHandler[fiCut]->GetMultiplicityBinIndex(fClusterCandidates->GetEntries());
2125         }
2126         
2127         if(((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->UseTrackMultiplicity()){
2128                 for(Int_t nEventsInBG=0;nEventsInBG<fBGHandler[fiCut]->GetNBGEvents();nEventsInBG++){
2129                         AliGammaConversionAODVector *previousEventV0s = fBGHandler[fiCut]->GetBGGoodV0s(zbin,mbin,nEventsInBG);
2130                         
2131                         for(Int_t iCurrent=0;iCurrent<fClusterCandidates->GetEntries();iCurrent++){
2132                                 AliAODConversionPhoton currentEventGoodV0 = *(AliAODConversionPhoton*)(fClusterCandidates->At(iCurrent));
2133                                 for(UInt_t iPrevious=0;iPrevious<previousEventV0s->size();iPrevious++){
2134                                         AliAODConversionPhoton previousGoodV0 = (AliAODConversionPhoton)(*(previousEventV0s->at(iPrevious)));
2135                                         AliAODConversionMother *backgroundCandidate = new AliAODConversionMother(&currentEventGoodV0,&previousGoodV0);
2136                                         backgroundCandidate->CalculateDistanceOfClossetApproachToPrimVtx(fInputEvent->GetPrimaryVertex());
2137                                         if((((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))
2138                                                 ->MesonIsSelected(backgroundCandidate,kFALSE,((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetEtaShift()))){
2139                                                 fHistoMotherBackInvMassPt[fiCut]->Fill(backgroundCandidate->M(),backgroundCandidate->Pt());
2140                                                 Double_t sparesFill[4] = {backgroundCandidate->M(),backgroundCandidate->Pt(),(Double_t)zbin,(Double_t)mbin};
2141                                                 fSparseMotherBackInvMassPtZM[fiCut]->Fill(sparesFill,1);
2142                                         }
2143                                         delete backgroundCandidate;
2144                                         backgroundCandidate = 0x0;
2145                                 }
2146                         }
2147                 }
2148         } else {
2149                 for(Int_t nEventsInBG=0;nEventsInBG <fBGHandler[fiCut]->GetNBGEvents();nEventsInBG++){
2150                         AliGammaConversionAODVector *previousEventV0s = fBGHandler[fiCut]->GetBGGoodV0s(zbin,mbin,nEventsInBG);
2151                         if(previousEventV0s){
2152                                 for(Int_t iCurrent=0;iCurrent<fClusterCandidates->GetEntries();iCurrent++){
2153                                         AliAODConversionPhoton currentEventGoodV0 = *(AliAODConversionPhoton*)(fClusterCandidates->At(iCurrent));
2154                                         for(UInt_t iPrevious=0;iPrevious<previousEventV0s->size();iPrevious++){
2155                                 
2156                                                 AliAODConversionPhoton previousGoodV0 = (AliAODConversionPhoton)(*(previousEventV0s->at(iPrevious)));
2157                                                 AliAODConversionMother *backgroundCandidate = new AliAODConversionMother(&currentEventGoodV0,&previousGoodV0);
2158                                                 backgroundCandidate->CalculateDistanceOfClossetApproachToPrimVtx(fInputEvent->GetPrimaryVertex());
2159                                                 if((((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->MesonIsSelected(backgroundCandidate,kFALSE,((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetEtaShift()))){
2160                                                         fHistoMotherBackInvMassPt[fiCut]->Fill(backgroundCandidate->M(),backgroundCandidate->Pt());
2161                                                         Double_t sparesFill[4] = {backgroundCandidate->M(),backgroundCandidate->Pt(),(Double_t)zbin,(Double_t)mbin};
2162                                                         fSparseMotherBackInvMassPtZM[fiCut]->Fill(sparesFill,1);
2163                                                 }
2164                                                 delete backgroundCandidate;
2165                                                 backgroundCandidate = 0x0;
2166                                         }
2167                                 }
2168                         }
2169                 }
2170         }
2171 }
2172
2173 //________________________________________________________________________
2174 void AliAnalysisTaskGammaCalo::RotateParticle(AliAODConversionPhoton *gamma){
2175         Int_t fNDegreesPMBackground= ((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->NDegreesRotation();
2176         Double_t nRadiansPM = fNDegreesPMBackground*TMath::Pi()/180;
2177         Double_t rotationValue = fRandom.Rndm()*2*nRadiansPM + TMath::Pi()-nRadiansPM;
2178         gamma->RotateZ(rotationValue);
2179 }
2180
2181 //________________________________________________________________________
2182 void AliAnalysisTaskGammaCalo::UpdateEventByEventData(){
2183         //see header file for documentation
2184         if(fClusterCandidates->GetEntries() >0 ){
2185                 if(((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->UseTrackMultiplicity()){
2186                         fBGHandler[fiCut]->AddEvent(fClusterCandidates,fInputEvent->GetPrimaryVertex()->GetX(),fInputEvent->GetPrimaryVertex()->GetY(),fInputEvent->GetPrimaryVertex()->GetZ(),fV0Reader->GetNumberOfPrimaryTracks(),fEventPlaneAngle);
2187                 } else { // means we use #V0s for multiplicity
2188                         fBGHandler[fiCut]->AddEvent(fClusterCandidates,fInputEvent->GetPrimaryVertex()->GetX(),fInputEvent->GetPrimaryVertex()->GetY(),fInputEvent->GetPrimaryVertex()->GetZ(),fClusterCandidates->GetEntries(),fEventPlaneAngle);
2189                 }
2190         }
2191 }
2192
2193
2194 //________________________________________________________________________
2195 void AliAnalysisTaskGammaCalo::FillPhotonBackgroundHist(AliAODConversionPhoton *TruePhotonCandidate, Int_t pdgCode)
2196 {
2197         // Bck = 0 e+-, 1 pi+-, 2 p+-, 3 K+-, 4 n, 5 K0s, 6 Lambda, 7 mu+-, 8 rest
2198         if(fIsFromMBHeader){
2199                 if(abs(pdgCode) == 11)                  fHistoClusPhotonBGPt[fiCut]->Fill(TruePhotonCandidate->Pt(),0);
2200                 else if( abs(pdgCode)==211)     fHistoClusPhotonBGPt[fiCut]->Fill(TruePhotonCandidate->Pt(),1);
2201                 else if( abs(pdgCode)==2212)    fHistoClusPhotonBGPt[fiCut]->Fill(TruePhotonCandidate->Pt(),2);
2202                 else if( abs(pdgCode)==321)     fHistoClusPhotonBGPt[fiCut]->Fill(TruePhotonCandidate->Pt(),3);
2203                 else if( abs(pdgCode)==2112)    fHistoClusPhotonBGPt[fiCut]->Fill(TruePhotonCandidate->Pt(),4);
2204                 else if( abs(pdgCode)==310)     fHistoClusPhotonBGPt[fiCut]->Fill(TruePhotonCandidate->Pt(),5);   
2205                 else if( abs(pdgCode)==3122)    fHistoClusPhotonBGPt[fiCut]->Fill(TruePhotonCandidate->Pt(),6);
2206                 else if( abs(pdgCode)==13)              fHistoClusPhotonBGPt[fiCut]->Fill(TruePhotonCandidate->Pt(),7);
2207                 else                                                    fHistoClusPhotonBGPt[fiCut]->Fill(TruePhotonCandidate->Pt(),8);
2208         }       
2209 }
2210
2211 //________________________________________________________________________
2212 void AliAnalysisTaskGammaCalo::FillPhotonPlusConversionBackgroundHist(AliAODConversionPhoton *TruePhotonCandidate, Int_t pdgCode)
2213 {
2214         // Bck = 0 e+-, 1 pi+-, 2 p+-, 3 K+-, 4 n, 5 K0s, 6 Lambda, 7 mu+-, 8 rest
2215         if(fIsFromMBHeader){
2216                 if(abs(pdgCode) == 11)                  fHistoClusPhotonPlusConvBGPt[fiCut]->Fill(TruePhotonCandidate->Pt(),0);
2217                 else if( abs(pdgCode)==211)     fHistoClusPhotonPlusConvBGPt[fiCut]->Fill(TruePhotonCandidate->Pt(),1);
2218                 else if( abs(pdgCode)==2212)    fHistoClusPhotonPlusConvBGPt[fiCut]->Fill(TruePhotonCandidate->Pt(),2);
2219                 else if( abs(pdgCode)==321)     fHistoClusPhotonPlusConvBGPt[fiCut]->Fill(TruePhotonCandidate->Pt(),3);
2220                 else if( abs(pdgCode)==2112)    fHistoClusPhotonPlusConvBGPt[fiCut]->Fill(TruePhotonCandidate->Pt(),4);
2221                 else if( abs(pdgCode)==310)     fHistoClusPhotonPlusConvBGPt[fiCut]->Fill(TruePhotonCandidate->Pt(),5);   
2222                 else if( abs(pdgCode)==3122)    fHistoClusPhotonPlusConvBGPt[fiCut]->Fill(TruePhotonCandidate->Pt(),6);
2223                 else if( abs(pdgCode)==13)              fHistoClusPhotonPlusConvBGPt[fiCut]->Fill(TruePhotonCandidate->Pt(),7);
2224                 else                                                    fHistoClusPhotonPlusConvBGPt[fiCut]->Fill(TruePhotonCandidate->Pt(),8);
2225         }       
2226 }
2227
2228
2229 void AliAnalysisTaskGammaCalo::SetLogBinningXTH2(TH2* histoRebin){
2230         TAxis *axisafter = histoRebin->GetXaxis();
2231         Int_t bins = axisafter->GetNbins();
2232         Double_t from = axisafter->GetXmin();
2233         Double_t to = axisafter->GetXmax();
2234         Double_t *newbins = new Double_t[bins+1];
2235         newbins[0] = from;
2236         Double_t factor = TMath::Power(to/from, 1./bins);
2237         for(Int_t i=1; i<=bins; ++i) newbins[i] = factor * newbins[i-1];
2238         axisafter->Set(bins, newbins);
2239         delete [] newbins;
2240 }
2241
2242 //________________________________________________________________________
2243 void AliAnalysisTaskGammaCalo::Terminate(const Option_t *)
2244 {
2245   
2246   //fOutputContainer->Print(); // Will crash on GRID
2247 }
2248
2249 //________________________________________________________________________
2250 Int_t AliAnalysisTaskGammaCalo::GetSourceClassification(Int_t daughter, Int_t pdgCode){
2251   
2252         if (daughter == 111) {
2253                 if (abs(pdgCode) == 310) return 1; // k0s
2254                 else if (abs(pdgCode) == 3122) return 2; // Lambda
2255                 else if (abs(pdgCode) == 130) return 3; // K0L
2256                 else if (abs(pdgCode) == 2212) return 4; // proton
2257                 else if (abs(pdgCode) == 2112) return 5; // neutron
2258                 else if (abs(pdgCode) == 211) return 6; // pion
2259                 else if (abs(pdgCode) == 321) return 7; // kaon
2260                 else if (abs(pdgCode) == 113 || abs(pdgCode) == 213 ) return 8; // rho 0,+,-
2261                 else if (abs(pdgCode) == 3222 || abs(pdgCode) == 3212 || abs(pdgCode) == 3112  ) return 9; // Sigma
2262                 else if (abs(pdgCode) == 2224 || abs(pdgCode) == 2214 || abs(pdgCode) == 2114 || abs(pdgCode) == 1114  ) return 10; // Delta
2263                 else if (abs(pdgCode) == 313 || abs(pdgCode) == 323   ) return 11; // K*
2264                 else return 15;
2265         }
2266         return 15;
2267   
2268 }