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