- changes by Lucia in order to include PbPb centrality trigger according to trigger...
[u/mrichter/AliRoot.git] / PWGGA / GammaConv / AliAnalysisTaskGammaConvV1.cxx
1 /**************************************************************************
2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3  *                                                                        *
4  * Author: Martin Wilde, Daniel Lohner, Friederike Bock                   *
5  * Version 1.0                                                            *
6  *                                                                        *
7  * based on: on older version (see aliroot up to v5-04-42-AN)             *
8  *           AliAnalysisTaskGammaConversion.cxx                           *
9  *           Authors: Kathrin Koch, Kenneth Aamodt, Ana Marin             *
10  *                                                                        *
11  * Permission to use, copy, modify and distribute this software and its   *
12  * documentation strictly for non-commercial purposes is hereby granted   *
13  * without fee, provided that the above copyright notice appears in all   *
14  * copies and that both the copyright notice and this permission notice   *
15  * appear in the supporting documentation. The authors make no claims     *
16  * about the suitability of this software for any purpose. It is          *
17  * provided "as is" without express or implied warranty.                  *
18  **************************************************************************/
19
20 ////////////////////////////////////////////////
21 //---------------------------------------------
22
23 // Class used to do analysis on conversion pairs
24 //---------------------------------------------
25 ///////////////////////////////////////////////
26 #include "TChain.h"
27 #include "TTree.h"
28 #include "TBranch.h"
29 #include "TFile.h"
30 #include "TH1F.h"
31 #include "TH2F.h"
32 #include "TH3F.h"
33 #include "THnSparse.h"
34 #include "TCanvas.h"
35 #include "TNtuple.h"
36 #include "AliAnalysisTask.h"
37 #include "AliAnalysisManager.h"
38 #include "AliESDEvent.h"
39 #include "AliESDInputHandler.h"
40 #include "AliMCEventHandler.h"
41 #include "AliMCEvent.h"
42 #include "AliMCParticle.h"
43 #include "AliCentrality.h"
44 #include "AliESDVZERO.h"
45 #include "AliESDpid.h"
46 #include "AliAnalysisTaskGammaConvV1.h"
47 #include "AliVParticle.h"
48 #include "AliESDtrack.h"
49 #include "AliESDtrackCuts.h"
50 #include "AliKFVertex.h"
51 #include "AliV0ReaderV1.h"
52 #include "AliGenCocktailEventHeader.h"
53 #include "AliConversionAODBGHandlerRP.h"
54 #include "AliAODMCParticle.h"
55 #include "AliAODMCHeader.h"
56 #include "AliEventplane.h"
57
58 ClassImp(AliAnalysisTaskGammaConvV1)
59
60 //________________________________________________________________________
61 AliAnalysisTaskGammaConvV1::AliAnalysisTaskGammaConvV1(): AliAnalysisTaskSE(),
62         fV0Reader(NULL),
63         fBGHandler(NULL),
64         fBGHandlerRP(NULL),
65         fInputEvent(NULL),
66         fMCEvent(NULL),
67         fMCStack(NULL),
68         fCutFolder(NULL),
69         fESDList(NULL),
70         fBackList(NULL),
71         fMotherList(NULL),
72         fPhotonDCAList(NULL),
73         fMesonDCAList(NULL),        
74         fTrueList(NULL),
75         fMCList(NULL),
76         fHeaderNameList(NULL),
77         fOutputContainer(0),
78         fReaderGammas(NULL),
79         fGammaCandidates(NULL),
80         fEventCutArray(NULL),
81         fEventCuts(NULL),
82         fCutArray(NULL),
83         fConversionCuts(NULL),
84         fMesonCutArray(NULL),
85         fMesonCuts(NULL),
86         hESDConvGammaPt(NULL),
87         hESDConvGammaR(NULL),
88         hESDConvGammaEta(NULL),
89         tESDConvGammaPtDcazCat(NULL),
90         fPtGamma(0),
91         fDCAzPhoton(0),
92         fRConvPhoton(0),
93         fEtaPhoton(0),
94         iCatPhoton(0),
95         iPhotonMCInfo(0),
96         hESDMotherInvMassPt(NULL),
97         sESDMotherInvMassPtZM(NULL),
98         hESDMotherBackInvMassPt(NULL),
99         sESDMotherBackInvMassPtZM(NULL),
100         hESDMotherInvMassEalpha(NULL),
101         hESDMotherPi0PtY(NULL),
102         hESDMotherEtaPtY(NULL),
103         hESDMotherPi0PtAlpha(NULL),
104         hESDMotherEtaPtAlpha(NULL),
105         hESDMotherPi0PtOpenAngle(NULL),
106         hESDMotherEtaPtOpenAngle(NULL),
107         hMCHeaders(NULL),
108         hMCAllGammaPt(NULL),
109         hMCDecayGammaPi0Pt(NULL),
110         hMCDecayGammaRhoPt(NULL),
111         hMCDecayGammaEtaPt(NULL),
112         hMCDecayGammaOmegaPt(NULL),
113         hMCDecayGammaEtapPt(NULL),
114         hMCDecayGammaPhiPt(NULL),
115         hMCDecayGammaSigmaPt(NULL),
116         hMCConvGammaPt(NULL),
117         hMCConvGammaR(NULL),
118         hMCConvGammaEta(NULL),
119         hMCPi0Pt(NULL),
120         hMCPi0WOWeightPt(NULL),
121         hMCEtaPt(NULL),
122         hMCEtaWOWeightPt(NULL),
123         hMCPi0InAccPt(NULL),
124         hMCEtaInAccPt(NULL),
125         hMCPi0PtY(NULL),
126         hMCEtaPtY(NULL),
127         hMCPi0PtAlpha(NULL),
128         hMCEtaPtAlpha(NULL),
129         hMCK0sPt(NULL),
130         hMCK0sWOWeightPt(NULL),
131         hMCK0sPtY(NULL),
132         hMCSecPi0PtvsSource(NULL),
133         hMCSecPi0Source(NULL),
134         hMCSecEtaPt(NULL),
135         hMCSecEtaSource(NULL),
136         hESDTrueMotherInvMassPt(NULL),
137         hESDTruePrimaryMotherInvMassPt(NULL),
138         hESDTruePrimaryMotherW0WeightingInvMassPt(NULL),
139         pESDTruePrimaryMotherWeightsInvMassPt(NULL),
140         hESDTruePrimaryPi0MCPtResolPt(NULL),
141         hESDTruePrimaryEtaMCPtResolPt(NULL),
142         hESDTrueSecondaryMotherInvMassPt(NULL),
143         hESDTrueSecondaryMotherFromK0sInvMassPt(NULL),
144         hESDTrueK0sWithPi0DaughterMCPt(NULL),
145         hESDTrueSecondaryMotherFromEtaInvMassPt(NULL),
146         hESDTrueEtaWithPi0DaughterMCPt(NULL),
147         hESDTrueSecondaryMotherFromLambdaInvMassPt(NULL),
148         hESDTrueLambdaWithPi0DaughterMCPt(NULL),
149         hESDTrueBckGGInvMassPt(NULL),
150         hESDTrueBckContInvMassPt(NULL),
151         hESDTruePi0PtY(NULL),
152         hESDTrueEtaPtY(NULL),
153         hESDTruePi0PtAlpha(NULL),
154         hESDTrueEtaPtAlpha(NULL),
155         hESDTruePi0PtOpenAngle(NULL),
156         hESDTrueEtaPtOpenAngle(NULL),
157         hESDTrueMotherDalitzInvMassPt(NULL),
158         hESDTrueConvGammaPt(NULL),
159         hESDTrueConvGammaEta(NULL),
160         hESDCombinatorialPt(NULL),
161         hESDTruePrimaryConvGammaPt(NULL),
162         hESDTruePrimaryConvGammaESDPtMCPt(NULL),
163         hESDTrueSecondaryConvGammaPt(NULL),
164         hESDTrueSecondaryConvGammaFromXFromK0sPt(NULL),
165         hESDTrueSecondaryConvGammaFromXFromLambdaPt(NULL),
166         hESDTrueDalitzPsiPairDeltaPhi(NULL),
167         hESDTrueGammaPsiPairDeltaPhi(NULL),
168         hNEvents(NULL),
169         hNGoodESDTracks(NULL),
170         hNGammaCandidates(NULL),
171         hNGoodESDTracksVsNGammaCanditates(NULL),
172         hNV0Tracks(NULL),
173         hEtaShift(NULL),
174         tESDMesonsInvMassPtDcazMinDcazMaxFlag(NULL),
175         fInvMass(0),
176         fPt(0),
177         fDCAzGammaMin(0),
178         fDCAzGammaMax(0),
179         iFlag(0),
180         iMesonMCInfo(0),
181         fEventPlaneAngle(-100),
182         fRandom(0),
183         fnGammaCandidates(0),
184         fUnsmearedPx(NULL),
185         fUnsmearedPy(NULL),
186         fUnsmearedPz(NULL),
187         fUnsmearedE(NULL),
188         fMCStackPos(NULL),
189         fMCStackNeg(NULL),
190         fESDArrayPos(NULL),
191         fESDArrayNeg(NULL),
192         fnCuts(0),
193         fiCut(0),
194         fMoveParticleAccordingToVertex(kTRUE),
195         fIsHeavyIon(0),
196         fDoMesonAnalysis(kTRUE),
197         fDoMesonQA(0),
198         fDoPhotonQA(0),
199         fIsFromMBHeader(kTRUE),
200         fIsMC(kFALSE)
201 {
202
203 }
204
205 //________________________________________________________________________
206 AliAnalysisTaskGammaConvV1::AliAnalysisTaskGammaConvV1(const char *name):
207         AliAnalysisTaskSE(name),
208         fV0Reader(NULL),
209         fBGHandler(NULL),
210         fBGHandlerRP(NULL),
211         fInputEvent(NULL),
212         fMCEvent(NULL),
213         fMCStack(NULL),
214         fCutFolder(NULL),
215         fESDList(NULL),
216         fBackList(NULL),
217         fMotherList(NULL),
218         fPhotonDCAList(NULL),
219         fMesonDCAList(NULL),
220         fTrueList(NULL),
221         fMCList(NULL),
222         fHeaderNameList(NULL),
223         fOutputContainer(0),
224         fReaderGammas(NULL),
225         fGammaCandidates(NULL),
226         fEventCutArray(NULL),
227         fEventCuts(NULL),
228         fCutArray(NULL),
229         fConversionCuts(NULL),
230         fMesonCutArray(NULL),
231         fMesonCuts(NULL),
232         hESDConvGammaPt(NULL),
233         hESDConvGammaR(NULL),
234         hESDConvGammaEta(NULL),
235         tESDConvGammaPtDcazCat(NULL),
236         fPtGamma(0),
237         fDCAzPhoton(0),
238         fRConvPhoton(0),
239         fEtaPhoton(0),
240         iCatPhoton(0),
241         iPhotonMCInfo(0),
242         hESDMotherInvMassPt(NULL),
243         sESDMotherInvMassPtZM(NULL),
244         hESDMotherBackInvMassPt(NULL),
245         sESDMotherBackInvMassPtZM(NULL),
246         hESDMotherInvMassEalpha(NULL),
247         hESDMotherPi0PtY(NULL),
248         hESDMotherEtaPtY(NULL),
249         hESDMotherPi0PtAlpha(NULL),
250         hESDMotherEtaPtAlpha(NULL),
251         hESDMotherPi0PtOpenAngle(NULL),
252         hESDMotherEtaPtOpenAngle(NULL),
253         hMCHeaders(NULL),
254         hMCAllGammaPt(NULL),
255         hMCDecayGammaPi0Pt(NULL),
256         hMCDecayGammaRhoPt(NULL),
257         hMCDecayGammaEtaPt(NULL),
258         hMCDecayGammaOmegaPt(NULL),
259         hMCDecayGammaEtapPt(NULL),
260         hMCDecayGammaPhiPt(NULL),
261         hMCDecayGammaSigmaPt(NULL),
262         hMCConvGammaPt(NULL),
263         hMCConvGammaR(NULL),
264         hMCConvGammaEta(NULL),
265         hMCPi0Pt(NULL),
266         hMCPi0WOWeightPt(NULL),
267         hMCEtaPt(NULL),
268         hMCEtaWOWeightPt(NULL),
269         hMCPi0InAccPt(NULL),
270         hMCEtaInAccPt(NULL),
271         hMCPi0PtY(NULL),
272         hMCEtaPtY(NULL),
273         hMCPi0PtAlpha(NULL),
274         hMCEtaPtAlpha(NULL),
275         hMCK0sPt(NULL),
276         hMCK0sWOWeightPt(NULL),
277         hMCK0sPtY(NULL),
278         hMCSecPi0PtvsSource(NULL),
279         hMCSecPi0Source(NULL),
280         hMCSecEtaPt(NULL),
281         hMCSecEtaSource(NULL),
282         hESDTrueMotherInvMassPt(NULL),
283         hESDTruePrimaryMotherInvMassPt(NULL),
284         hESDTruePrimaryMotherW0WeightingInvMassPt(NULL),
285         pESDTruePrimaryMotherWeightsInvMassPt(NULL),
286         hESDTruePrimaryPi0MCPtResolPt(NULL),
287         hESDTruePrimaryEtaMCPtResolPt(NULL),
288         hESDTrueSecondaryMotherInvMassPt(NULL),
289         hESDTrueSecondaryMotherFromK0sInvMassPt(NULL),
290         hESDTrueK0sWithPi0DaughterMCPt(NULL),
291         hESDTrueSecondaryMotherFromEtaInvMassPt(NULL),
292         hESDTrueEtaWithPi0DaughterMCPt(NULL),
293         hESDTrueSecondaryMotherFromLambdaInvMassPt(NULL),
294         hESDTrueLambdaWithPi0DaughterMCPt(NULL),
295         hESDTrueBckGGInvMassPt(NULL),
296         hESDTrueBckContInvMassPt(NULL),
297         hESDTruePi0PtY(NULL),
298         hESDTrueEtaPtY(NULL),
299         hESDTruePi0PtAlpha(NULL),
300         hESDTrueEtaPtAlpha(NULL),
301         hESDTruePi0PtOpenAngle(NULL),
302         hESDTrueEtaPtOpenAngle(NULL),
303         hESDTrueMotherDalitzInvMassPt(NULL),
304         hESDTrueConvGammaPt(NULL),
305         hESDTrueConvGammaEta(NULL),
306         hESDCombinatorialPt(NULL),
307         hESDTruePrimaryConvGammaPt(NULL),
308         hESDTruePrimaryConvGammaESDPtMCPt(NULL),
309         hESDTrueSecondaryConvGammaPt(NULL),
310         hESDTrueSecondaryConvGammaFromXFromK0sPt(NULL),
311         hESDTrueSecondaryConvGammaFromXFromLambdaPt(NULL),
312         hESDTrueDalitzPsiPairDeltaPhi(NULL),
313         hESDTrueGammaPsiPairDeltaPhi(NULL),
314         hNEvents(NULL),
315         hNGoodESDTracks(NULL),
316         hNGammaCandidates(NULL),
317         hNGoodESDTracksVsNGammaCanditates(NULL),
318         hNV0Tracks(NULL),
319         hEtaShift(NULL),
320         tESDMesonsInvMassPtDcazMinDcazMaxFlag(NULL),
321         fInvMass(0),
322         fPt(0),
323         fDCAzGammaMin(0),
324         fDCAzGammaMax(0),
325         iFlag(0),
326         iMesonMCInfo(0),
327         fEventPlaneAngle(-100),
328         fRandom(0),
329         fnGammaCandidates(0),
330         fUnsmearedPx(NULL),
331         fUnsmearedPy(NULL),
332         fUnsmearedPz(NULL),
333         fUnsmearedE(NULL),
334         fMCStackPos(NULL),
335         fMCStackNeg(NULL),
336         fESDArrayPos(NULL),
337         fESDArrayNeg(NULL),
338         fnCuts(0),
339         fiCut(0),
340         fMoveParticleAccordingToVertex(kTRUE),
341         fIsHeavyIon(0),
342         fDoMesonAnalysis(kTRUE),
343         fDoMesonQA(0),
344         fDoPhotonQA(0),
345         fIsFromMBHeader(kTRUE),
346         fIsMC(kFALSE)
347 {
348    // Define output slots here
349         DefineOutput(1, TList::Class());
350 }
351
352 AliAnalysisTaskGammaConvV1::~AliAnalysisTaskGammaConvV1()
353 {
354         if(fGammaCandidates){
355                 delete fGammaCandidates;
356                 fGammaCandidates = 0x0;
357         }
358         if(fBGHandler){
359                 delete[] fBGHandler;
360                 fBGHandler = 0x0;
361         }
362         if(fBGHandlerRP){
363                 delete[] fBGHandlerRP;
364                 fBGHandlerRP = 0x0;
365         }
366 }
367 //___________________________________________________________
368 void AliAnalysisTaskGammaConvV1::InitBack(){
369
370         const Int_t nDim = 4;
371         Int_t nBins[nDim] = {800,250,7,4};
372         Double_t xMin[nDim] = {0,0, 0,0};
373         Double_t xMax[nDim] = {0.8,25,7,4};
374         
375         sESDMotherInvMassPtZM = new THnSparseF*[fnCuts];
376         sESDMotherBackInvMassPtZM = new THnSparseF*[fnCuts];
377
378         fBGHandler = new AliGammaConversionAODBGHandler*[fnCuts];
379         fBGHandlerRP = new AliConversionAODBGHandlerRP*[fnCuts];
380         for(Int_t iCut = 0; iCut<fnCuts;iCut++){
381                 if (((AliConversionMesonCuts*)fMesonCutArray->At(iCut))->DoBGCalculation()){
382                         TString cutstringEvent  = ((AliConvEventCuts*)fEventCutArray->At(iCut))->GetCutNumber();
383                         TString cutstringPhoton         = ((AliConversionPhotonCuts*)fCutArray->At(iCut))->GetCutNumber();
384                         TString cutstringMeson  = ((AliConversionMesonCuts*)fMesonCutArray->At(iCut))->GetCutNumber();
385                         
386                         Int_t collisionSystem = atoi((TString)(((AliConvEventCuts*)fEventCutArray->At(iCut))->GetCutNumber())(0,1));
387                         Int_t centMin = atoi((TString)(((AliConvEventCuts*)fEventCutArray->At(iCut))->GetCutNumber())(1,1));
388                         Int_t centMax = atoi((TString)(((AliConvEventCuts*)fEventCutArray->At(iCut))->GetCutNumber())(2,1));
389                         
390                         if(collisionSystem == 1 || collisionSystem == 2 ||
391                                 collisionSystem == 5 || collisionSystem == 8 ||
392                                 collisionSystem == 9){
393                                 centMin = centMin*10;
394                                 centMax = centMax*10; 
395                                 if(centMax ==0 && centMax!=centMin) centMax=100;
396                         } else if(collisionSystem == 3 || collisionSystem == 6) {
397                                 centMin = centMin*5;
398                                 centMax = centMax*5;
399                         } else if(collisionSystem == 4 || collisionSystem == 7) {
400                                 centMin = ((centMin*5)+45);
401                                 centMax = ((centMax*5)+45);
402                         }
403                         
404                         fBackList[iCut] = new TList();
405                         fBackList[iCut]->SetName(Form("%s_%s_%s Back histograms",cutstringEvent.Data(), cutstringPhoton.Data(),cutstringMeson.Data()));
406                         fBackList[iCut]->SetOwner(kTRUE);
407                         fCutFolder[iCut]->Add(fBackList[iCut]);
408
409                         sESDMotherBackInvMassPtZM[iCut] = new THnSparseF("Back_Back_InvMass_Pt_z_m","Back_Back_InvMass_Pt_z_m",nDim,nBins,xMin,xMax);
410                         fBackList[iCut]->Add(sESDMotherBackInvMassPtZM[iCut]);
411
412                         fMotherList[iCut] = new TList();
413                         fMotherList[iCut]->SetName(Form("%s_%s_%s Mother histograms",cutstringEvent.Data() ,cutstringPhoton.Data(),cutstringMeson.Data()));
414                         fMotherList[iCut]->SetOwner(kTRUE);
415                         fCutFolder[iCut]->Add(fMotherList[iCut]);
416
417                         sESDMotherInvMassPtZM[iCut] = new THnSparseF("Back_Mother_InvMass_Pt_z_m","Back_Mother_InvMass_Pt_z_m",nDim,nBins,xMin,xMax);
418                         fMotherList[iCut]->Add(sESDMotherInvMassPtZM[iCut]);
419
420                         if(((AliConversionMesonCuts*)fMesonCutArray->At(iCut))->BackgroundHandlerType() == 0){
421                                 fBGHandler[iCut] = new AliGammaConversionAODBGHandler(
422                                                                                                                                         collisionSystem,centMin,centMax,
423                                                                                                                                         ((AliConversionMesonCuts*)fMesonCutArray->At(iCut))->GetNumberOfBGEvents(),
424                                                                                                                                         ((AliConversionMesonCuts*)fMesonCutArray->At(iCut))->UseTrackMultiplicity(),
425                                                                                                                                         0,8,5);
426                                 fBGHandlerRP[iCut] = NULL;
427                         } else {
428                                 fBGHandlerRP[iCut] = new AliConversionAODBGHandlerRP(
429                                                                                                                                         ((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsHeavyIon(),
430                                                                                                                                         ((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->UseTrackMultiplicity(),
431                                                                                                                                         ((AliConversionMesonCuts*)fMesonCutArray->At(iCut))->GetNumberOfBGEvents());
432                                 fBGHandler[iCut] = NULL;
433                         }
434                 }
435         }
436 }
437 //________________________________________________________________________
438 void AliAnalysisTaskGammaConvV1::UserCreateOutputObjects(){
439
440         // Create histograms
441         if(fOutputContainer != NULL){
442                 delete fOutputContainer;
443                 fOutputContainer = NULL;
444         }
445         if(fOutputContainer == NULL){
446                 fOutputContainer = new TList();
447                 fOutputContainer->SetOwner(kTRUE);
448         }
449
450         // Array of current cut's gammas
451         fGammaCandidates = new TList();
452
453         fCutFolder = new TList*[fnCuts];
454         fESDList = new TList*[fnCuts];
455         fBackList = new TList*[fnCuts];
456         fMotherList = new TList*[fnCuts];
457         hNEvents = new TH1I*[fnCuts];
458         hNGoodESDTracks = new TH1I*[fnCuts];
459         hNGammaCandidates = new TH1I*[fnCuts];
460         hNGoodESDTracksVsNGammaCanditates = new TH2F*[fnCuts];
461         hNV0Tracks = new TH1I*[fnCuts];
462         hEtaShift = new TProfile*[fnCuts];
463         hESDConvGammaPt = new TH1F*[fnCuts];
464
465         if (fDoPhotonQA == 2){
466                 fPhotonDCAList = new TList*[fnCuts];
467                 tESDConvGammaPtDcazCat = new TTree*[fnCuts];
468         }
469         if (fDoPhotonQA > 0){
470                 hESDConvGammaR = new TH1F*[fnCuts];
471                 hESDConvGammaEta = new TH1F*[fnCuts];
472         } 
473         
474         if(fDoMesonAnalysis){
475                 hESDMotherInvMassPt = new TH2F*[fnCuts];
476                 hESDMotherBackInvMassPt = new TH2F*[fnCuts];
477                 hESDMotherInvMassEalpha = new TH2F*[fnCuts];
478                 if (fDoMesonQA == 2){
479                         fMesonDCAList = new TList*[fnCuts];
480                         tESDMesonsInvMassPtDcazMinDcazMaxFlag = new TTree*[fnCuts];
481                 }
482                 if (fDoMesonQA > 0){
483                         hESDMotherPi0PtY =  new TH2F*[fnCuts];
484                         hESDMotherEtaPtY =  new TH2F*[fnCuts];
485                         hESDMotherPi0PtAlpha =  new TH2F*[fnCuts];
486                         hESDMotherEtaPtAlpha =  new TH2F*[fnCuts];
487                         hESDMotherPi0PtOpenAngle =  new TH2F*[fnCuts];
488                         hESDMotherEtaPtOpenAngle =  new TH2F*[fnCuts];
489                 }   
490         }
491
492         for(Int_t iCut = 0; iCut<fnCuts;iCut++){
493
494                 TString cutstringEvent  = ((AliConvEventCuts*)fEventCutArray->At(iCut))->GetCutNumber();
495                 TString cutstringPhoton = ((AliConversionPhotonCuts*)fCutArray->At(iCut))->GetCutNumber();
496                 TString cutstringMeson = "NoMesonCut";
497                 if(fDoMesonAnalysis)cutstringMeson = ((AliConversionMesonCuts*)fMesonCutArray->At(iCut))->GetCutNumber();
498
499                 fCutFolder[iCut] = new TList();
500                 fCutFolder[iCut]->SetName(Form("Cut Number %s_%s_%s",cutstringEvent.Data() ,cutstringPhoton.Data(),cutstringMeson.Data()));
501                 fCutFolder[iCut]->SetOwner(kTRUE);
502                 fOutputContainer->Add(fCutFolder[iCut]);
503                 fESDList[iCut] = new TList();
504                 fESDList[iCut]->SetName(Form("%s_%s_%s ESD histograms",cutstringEvent.Data() ,cutstringPhoton.Data(),cutstringMeson.Data()));
505                 fESDList[iCut]->SetOwner(kTRUE);
506                 fCutFolder[iCut]->Add(fESDList[iCut]);
507
508                 hNEvents[iCut] = new TH1I("NEvents","NEvents",9,-0.5,8.5);
509                 hNEvents[iCut]->GetXaxis()->SetBinLabel(1,"Accepted");
510                 hNEvents[iCut]->GetXaxis()->SetBinLabel(2,"Centrality");
511                 hNEvents[iCut]->GetXaxis()->SetBinLabel(3,"Missing MC");
512                 if (((AliConvEventCuts*)fEventCutArray->At(iCut))->IsSpecialTrigger() > 1 ){ 
513                         TString TriggerNames = "Not Trigger: ";
514                         TriggerNames = TriggerNames+ ( (AliConvEventCuts*)fEventCutArray->At(iCut))->GetSpecialTriggerName();
515                         hNEvents[iCut]->GetXaxis()->SetBinLabel(4,TriggerNames.Data());
516                 } else {
517                         hNEvents[iCut]->GetXaxis()->SetBinLabel(4,"Trigger");
518                 }
519                 hNEvents[iCut]->GetXaxis()->SetBinLabel(5,"Vertex Z");
520                 hNEvents[iCut]->GetXaxis()->SetBinLabel(6,"Cont. Vertex");
521                 hNEvents[iCut]->GetXaxis()->SetBinLabel(7,"Pile-Up");
522                 hNEvents[iCut]->GetXaxis()->SetBinLabel(8,"no SDD");
523                 hNEvents[iCut]->GetXaxis()->SetBinLabel(9,"no V0AND");
524                 fESDList[iCut]->Add(hNEvents[iCut]);
525                 
526                 if(fIsHeavyIon == 1) hNGoodESDTracks[iCut] = new TH1I("GoodESDTracks","GoodESDTracks",4000,0,4000);
527                 else if(fIsHeavyIon == 2) hNGoodESDTracks[iCut] = new TH1I("GoodESDTracks","GoodESDTracks",400,0,400);
528                 else hNGoodESDTracks[iCut] = new TH1I("GoodESDTracks","GoodESDTracks",200,0,200);
529                 fESDList[iCut]->Add(hNGoodESDTracks[iCut]);
530                 if(fIsHeavyIon == 1) hNGammaCandidates[iCut] = new TH1I("GammaCandidates","GammaCandidates",100,0,100);
531                 else if(fIsHeavyIon == 2) hNGammaCandidates[iCut] = new TH1I("GammaCandidates","GammaCandidates",50,0,50);
532                 else hNGammaCandidates[iCut] = new TH1I("GammaCandidates","GammaCandidates",50,0,50);
533                 fESDList[iCut]->Add(hNGammaCandidates[iCut]);
534                 if(fIsHeavyIon == 1) hNGoodESDTracksVsNGammaCanditates[iCut] = new TH2F("GoodESDTracksVsGammaCandidates","GoodESDTracksVsGammaCandidates",4000,0,4000,100,0,100);
535                 else if(fIsHeavyIon == 2) hNGoodESDTracksVsNGammaCanditates[iCut] = new TH2F("GoodESDTracksVsGammaCandidates","GoodESDTracksVsGammaCandidates",400,0,400,50,0,50);
536                 else hNGoodESDTracksVsNGammaCanditates[iCut] = new TH2F("GoodESDTracksVsGammaCandidates","GoodESDTracksVsGammaCandidates",200,0,200,50,0,50);
537                 fESDList[iCut]->Add(hNGoodESDTracksVsNGammaCanditates[iCut]);
538
539                 
540                 if(fIsHeavyIon == 1) hNV0Tracks[iCut] = new TH1I("V0 Multiplicity","V0 Multiplicity",30000,0,30000);
541                 else if(fIsHeavyIon == 2) hNV0Tracks[iCut] = new TH1I("V0 Multiplicity","V0 Multiplicity",2500,0,2500);
542                 else hNV0Tracks[iCut] = new TH1I("V0 Multiplicity","V0 Multiplicity",1500,0,1500);
543                 fESDList[iCut]->Add(hNV0Tracks[iCut]);
544                 hEtaShift[iCut] = new TProfile("Eta Shift","Eta Shift",1, -0.5,0.5);
545                 fESDList[iCut]->Add(hEtaShift[iCut]);
546                 hESDConvGammaPt[iCut] = new TH1F("ESD_ConvGamma_Pt","ESD_ConvGamma_Pt",250,0,25);
547                 fESDList[iCut]->Add(hESDConvGammaPt[iCut]);
548
549                 if (fDoPhotonQA == 2){
550                         fPhotonDCAList[iCut] = new TList();
551                         fPhotonDCAList[iCut]->SetName(Form("%s_%s_%s Photon DCA tree",cutstringEvent.Data(),cutstringPhoton.Data(),cutstringMeson.Data()));
552                         fPhotonDCAList[iCut]->SetOwner(kTRUE);
553                         fCutFolder[iCut]->Add(fPhotonDCAList[iCut]);
554                                 
555                         tESDConvGammaPtDcazCat[iCut] = new TTree("ESD_ConvGamma_Pt_Dcaz_R_Eta","ESD_ConvGamma_Pt_Dcaz_R_Eta_Cat");   
556                         tESDConvGammaPtDcazCat[iCut]->Branch("Pt",&fPtGamma,"fPtGamma/F");
557                         tESDConvGammaPtDcazCat[iCut]->Branch("DcaZPhoton",&fDCAzPhoton,"fDCAzPhoton/F");
558         //          tESDConvGammaPtDcazCat[iCut]->Branch("R",&fRConvPhoton,"fRConvPhoton/F");
559         //          tESDConvGammaPtDcazCat[iCut]->Branch("Eta",&fEtaPhoton,"fEtaPhoton/F");
560                         
561                         tESDConvGammaPtDcazCat[iCut]->Branch("cat",&iCatPhoton,"iCatPhoton/b");
562                         if(fIsMC){
563                                 tESDConvGammaPtDcazCat[iCut]->Branch("photonMCInfo",&iPhotonMCInfo,"iPhotonMCInfo/b");
564                         }
565                         fPhotonDCAList[iCut]->Add(tESDConvGammaPtDcazCat[iCut]);
566                 }
567
568                 if (fDoPhotonQA > 0){
569                         hESDConvGammaR[iCut] = new TH1F("ESD_ConvGamma_R","ESD_ConvGamma_R",800,0,200);
570                         fESDList[iCut]->Add(hESDConvGammaR[iCut]);
571                         hESDConvGammaEta[iCut] = new TH1F("ESD_ConvGamma_Eta","ESD_ConvGamma_Eta",2000,-2,2);
572                         fESDList[iCut]->Add(hESDConvGammaEta[iCut]);
573                 }
574                 
575                 if(fDoMesonAnalysis){
576                         hESDMotherInvMassPt[iCut] = new TH2F("ESD_Mother_InvMass_Pt","ESD_Mother_InvMass_Pt",800,0,0.8,250,0,25);
577                         fESDList[iCut]->Add(hESDMotherInvMassPt[iCut]);
578                         hESDMotherBackInvMassPt[iCut] = new TH2F("ESD_Background_InvMass_Pt","ESD_Background_InvMass_Pt",800,0,0.8,250,0,25);
579                         fESDList[iCut]->Add(hESDMotherBackInvMassPt[iCut]);
580                         hESDMotherInvMassEalpha[iCut] = new TH2F("ESD_Mother_InvMass_vs_E_alpha","ESD_Mother_InvMass_vs_E_alpha",800,0,0.8,250,0,25);
581                         fESDList[iCut]->Add(hESDMotherInvMassEalpha[iCut]);
582                         if (fDoMesonQA == 2){    
583                                 fMesonDCAList[iCut] = new TList();
584                                 fMesonDCAList[iCut]->SetName(Form("%s_%s_%s Meson DCA tree",cutstringEvent.Data() ,cutstringPhoton.Data(),cutstringMeson.Data()));
585                                 fMesonDCAList[iCut]->SetOwner(kTRUE);
586                                 fCutFolder[iCut]->Add(fMesonDCAList[iCut]);
587                                 
588                                 tESDMesonsInvMassPtDcazMinDcazMaxFlag[iCut] = new TTree("ESD_Mesons_InvMass_Pt_DcazMin_DcazMax_Flag","ESD_Mesons_InvMass_Pt_DcazMin_DcazMax_Flag");   
589                                 tESDMesonsInvMassPtDcazMinDcazMaxFlag[iCut]->Branch("InvMass",&fInvMass,"fInvMass/F");
590                                 tESDMesonsInvMassPtDcazMinDcazMaxFlag[iCut]->Branch("Pt",&fPt,"fPt/F");
591                                 tESDMesonsInvMassPtDcazMinDcazMaxFlag[iCut]->Branch("DcaZMin",&fDCAzGammaMin,"fDCAzGammaMin/F");
592                                 tESDMesonsInvMassPtDcazMinDcazMaxFlag[iCut]->Branch("DcaZMax",&fDCAzGammaMax,"fDCAzGammaMax/F");
593                                 tESDMesonsInvMassPtDcazMinDcazMaxFlag[iCut]->Branch("kind",&iFlag,"iFlag/b");
594                                 if(fIsMC){
595                                 tESDMesonsInvMassPtDcazMinDcazMaxFlag[iCut]->Branch("mesonMCInfo",&iMesonMCInfo,"iMesonMCInfo/b");
596                                 }
597                                 fMesonDCAList[iCut]->Add(tESDMesonsInvMassPtDcazMinDcazMaxFlag[iCut]);
598                                 
599                         }
600                         if (fDoMesonQA > 0 ){
601                                 hESDMotherPi0PtY[iCut] = new TH2F("ESD_MotherPi0_Pt_Y","ESD_MotherPi0_Pt_Y",150,0.03,15.,150,-1.5,1.5);            
602                                 SetLogBinningXTH2(hESDMotherPi0PtY[iCut]);
603                                 fESDList[iCut]->Add(hESDMotherPi0PtY[iCut]);
604                                 hESDMotherEtaPtY[iCut] = new TH2F("ESD_MotherEta_Pt_Y","ESD_MotherEta_Pt_Y",150,0.03,15.,150,-1.5,1.5);
605                                 SetLogBinningXTH2(hESDMotherEtaPtY[iCut]);
606                                 fESDList[iCut]->Add(hESDMotherEtaPtY[iCut]);
607                                 hESDMotherPi0PtAlpha[iCut] = new TH2F("ESD_MotherPi0_Pt_Alpha","ESD_MotherPi0_Pt_Alpha",150,0.03,15.,100,0,1);            
608                                 SetLogBinningXTH2(hESDMotherPi0PtAlpha[iCut]);
609                                 fESDList[iCut]->Add(hESDMotherPi0PtAlpha[iCut]);
610                                 hESDMotherEtaPtAlpha[iCut] = new TH2F("ESD_MotherEta_Pt_Alpha","ESD_MotherEta_Pt_Alpha",150,0.03,15.,100,0,1);
611                                 SetLogBinningXTH2(hESDMotherEtaPtAlpha[iCut]);
612                                 fESDList[iCut]->Add(hESDMotherEtaPtAlpha[iCut]);
613                                 hESDMotherPi0PtOpenAngle[iCut] = new TH2F("ESD_MotherPi0_Pt_OpenAngle","ESD_MotherPi0_Pt_OpenAngle",150,0.03,15.,100,0,TMath::Pi());            
614                                 SetLogBinningXTH2(hESDMotherPi0PtOpenAngle[iCut]);
615                                 fESDList[iCut]->Add(hESDMotherPi0PtOpenAngle[iCut]);
616                                 hESDMotherEtaPtOpenAngle[iCut] = new TH2F("ESD_MotherEta_Pt_OpenAngle","ESD_MotherEta_Pt_OpenAngle",150,0.03,15.,100,0,TMath::Pi());            
617                                 SetLogBinningXTH2(hESDMotherEtaPtOpenAngle[iCut]);
618                                 fESDList[iCut]->Add(hESDMotherEtaPtOpenAngle[iCut]);
619                         }
620                         
621                                 
622                 }
623
624
625         }
626         if(fDoMesonAnalysis){
627                 InitBack(); // Init Background Handler
628         }
629
630         if(fIsMC){
631                 // MC Histogramms
632                 fMCList = new TList*[fnCuts];
633                 // True Histogramms
634                 fTrueList = new TList*[fnCuts];
635                 // Selected Header List
636                 fHeaderNameList = new TList*[fnCuts];
637                 hMCHeaders = new TH1I*[fnCuts];
638                 hMCAllGammaPt = new TH1F*[fnCuts];
639                 hMCDecayGammaPi0Pt = new TH1F*[fnCuts];
640                 hMCDecayGammaRhoPt = new TH1F*[fnCuts];
641                 hMCDecayGammaEtaPt = new TH1F*[fnCuts];
642                 hMCDecayGammaOmegaPt = new TH1F*[fnCuts];
643                 hMCDecayGammaEtapPt = new TH1F*[fnCuts];
644                 hMCDecayGammaPhiPt = new TH1F*[fnCuts];
645                 hMCDecayGammaSigmaPt = new TH1F*[fnCuts];
646                 hMCConvGammaPt = new TH1F*[fnCuts];
647                 hESDTrueConvGammaPt = new TH1F*[fnCuts];
648
649                 hESDCombinatorialPt = new TH2F*[fnCuts];
650                 hESDTruePrimaryConvGammaPt = new TH1F*[fnCuts];
651                 hESDTruePrimaryConvGammaESDPtMCPt = new TH2F*[fnCuts];
652                 hESDTrueSecondaryConvGammaPt = new TH1F*[fnCuts];
653                 hESDTrueSecondaryConvGammaFromXFromK0sPt = new TH1F*[fnCuts];
654                 hESDTrueSecondaryConvGammaFromXFromLambdaPt = new TH1F*[fnCuts];
655
656                 hESDTrueDalitzPsiPairDeltaPhi= new TH2F*[fnCuts];
657                 hESDTrueGammaPsiPairDeltaPhi= new TH2F*[fnCuts];
658
659                 if (fDoPhotonQA > 0){
660                         hMCConvGammaR = new TH1F*[fnCuts];
661                         hMCConvGammaEta = new TH1F*[fnCuts];
662                 hESDTrueConvGammaEta = new TH1F*[fnCuts];
663                 }
664
665                 if(fDoMesonAnalysis){
666                         hMCPi0Pt = new TH1F*[fnCuts];
667                         hMCPi0WOWeightPt = new TH1F*[fnCuts];
668                         hMCEtaPt = new TH1F*[fnCuts];
669                         hMCEtaWOWeightPt = new TH1F*[fnCuts];
670                         hMCPi0InAccPt = new TH1F*[fnCuts];
671                         hMCEtaInAccPt = new TH1F*[fnCuts];
672
673                         hESDTrueMotherInvMassPt = new TH2F*[fnCuts];
674                         hESDTruePrimaryMotherInvMassPt = new TH2F*[fnCuts];
675                         hESDTruePrimaryMotherW0WeightingInvMassPt = new TH2F*[fnCuts];
676                         pESDTruePrimaryMotherWeightsInvMassPt = new TProfile2D*[fnCuts];
677                         hESDTrueSecondaryMotherInvMassPt = new TH2F*[fnCuts];
678                         hESDTrueSecondaryMotherFromK0sInvMassPt = new TH2F*[fnCuts];
679                         hESDTrueSecondaryMotherFromEtaInvMassPt = new TH2F*[fnCuts];
680                         hESDTrueSecondaryMotherFromLambdaInvMassPt = new TH2F*[fnCuts];
681                         hESDTrueMotherDalitzInvMassPt = new TH2F*[fnCuts];
682                         if (fDoMesonQA > 0){
683                                 hMCPi0PtY = new TH2F*[fnCuts];
684                                 hMCEtaPtY = new TH2F*[fnCuts];
685                                 hMCPi0PtAlpha = new TH2F*[fnCuts];
686                                 hMCEtaPtAlpha = new TH2F*[fnCuts];
687                                 hMCK0sPt = new TH1F*[fnCuts];
688                                 hMCK0sWOWeightPt = new TH1F*[fnCuts];
689                                 hMCK0sPtY = new TH2F*[fnCuts];
690                                 hMCSecPi0PtvsSource= new TH2F*[fnCuts];
691                                 hMCSecPi0Source = new TH1F*[fnCuts];
692                                 hMCSecEtaPt = new TH1F*[fnCuts];
693                                 hMCSecEtaSource = new TH1F*[fnCuts];
694                                 hESDTruePrimaryPi0MCPtResolPt = new TH2F*[fnCuts];
695                                 hESDTruePrimaryEtaMCPtResolPt = new TH2F*[fnCuts];
696                                 hESDTrueK0sWithPi0DaughterMCPt = new TH1F*[fnCuts];
697                                 hESDTrueEtaWithPi0DaughterMCPt = new TH1F*[fnCuts];
698                                 hESDTrueLambdaWithPi0DaughterMCPt = new TH1F*[fnCuts];
699                                 hESDTrueBckGGInvMassPt = new TH2F*[fnCuts];
700                                 hESDTrueBckContInvMassPt = new TH2F*[fnCuts];
701                                 hESDTruePi0PtY = new TH2F*[fnCuts];
702                                 hESDTrueEtaPtY = new TH2F*[fnCuts];
703                                 hESDTruePi0PtAlpha = new TH2F*[fnCuts];
704                                 hESDTrueEtaPtAlpha = new TH2F*[fnCuts];
705                                 hESDTruePi0PtOpenAngle = new TH2F*[fnCuts];
706                                 hESDTrueEtaPtOpenAngle = new TH2F*[fnCuts];
707                         }
708                 }
709
710                 for(Int_t iCut = 0; iCut<fnCuts;iCut++){
711                         TString cutstringEvent  = ((AliConvEventCuts*)fEventCutArray->At(iCut))->GetCutNumber();
712                         TString cutstringPhoton = ((AliConversionPhotonCuts*)fCutArray->At(iCut))->GetCutNumber();
713                         TString cutstringMeson = "NoMesonCut";
714                         if(fDoMesonAnalysis)cutstringMeson = ((AliConversionMesonCuts*)fMesonCutArray->At(iCut))->GetCutNumber();
715
716                         fMCList[iCut] = new TList();
717                         fMCList[iCut]->SetName(Form("%s_%s_%s MC histograms",cutstringEvent.Data() ,cutstringPhoton.Data(),cutstringMeson.Data()));
718                         fMCList[iCut]->SetOwner(kTRUE);
719                         fCutFolder[iCut]->Add(fMCList[iCut]);
720                         hMCHeaders[iCut] = new TH1I("MC_Headers","MC_Headers",20,0,20);
721                         fMCList[iCut]->Add(hMCHeaders[iCut]);
722                         hMCAllGammaPt[iCut] = new TH1F("MC_AllGamma_Pt","MC_AllGamma_Pt",250,0,25);
723                         fMCList[iCut]->Add(hMCAllGammaPt[iCut]);
724                         hMCDecayGammaPi0Pt[iCut] = new TH1F("MC_DecayGammaPi0_Pt","MC_DecayGammaPi0_Pt",250,0,25);
725                         fMCList[iCut]->Add(hMCDecayGammaPi0Pt[iCut]);
726                         hMCDecayGammaRhoPt[iCut] = new TH1F("MC_DecayGammaRho_Pt","MC_DecayGammaRho_Pt",250,0,25);
727                         fMCList[iCut]->Add(hMCDecayGammaRhoPt[iCut]);
728                         hMCDecayGammaEtaPt[iCut] = new TH1F("MC_DecayGammaEta_Pt","MC_DecayGammaEta_Pt",250,0,25);
729                         fMCList[iCut]->Add(hMCDecayGammaEtaPt[iCut]);
730                         hMCDecayGammaOmegaPt[iCut] = new TH1F("MC_DecayGammaOmega_Pt","MC_DecayGammaOmmega_Pt",250,0,25);
731                         fMCList[iCut]->Add(hMCDecayGammaOmegaPt[iCut]);
732                         hMCDecayGammaEtapPt[iCut] = new TH1F("MC_DecayGammaEtap_Pt","MC_DecayGammaEtap_Pt",250,0,25);
733                         fMCList[iCut]->Add(hMCDecayGammaEtapPt[iCut]);
734                         hMCDecayGammaPhiPt[iCut] = new TH1F("MC_DecayGammaPhi_Pt","MC_DecayGammaPhi_Pt",250,0,25);
735                         fMCList[iCut]->Add(hMCDecayGammaPhiPt[iCut]);
736                         hMCDecayGammaSigmaPt[iCut] = new TH1F("MC_DecayGammaSigma_Pt","MC_DecayGammaSigma_Pt",250,0,25);
737                         fMCList[iCut]->Add(hMCDecayGammaSigmaPt[iCut]);
738                         hMCConvGammaPt[iCut] = new TH1F("MC_ConvGamma_Pt","MC_ConvGamma_Pt",250,0,25);
739                         fMCList[iCut]->Add(hMCConvGammaPt[iCut]);
740                 
741                         if (fDoPhotonQA > 0){
742                                 hMCConvGammaR[iCut] = new TH1F("MC_ConvGamma_R","MC_ConvGamma_R",800,0,200);
743                                 fMCList[iCut]->Add(hMCConvGammaR[iCut]);
744                                 hMCConvGammaEta[iCut] = new TH1F("MC_ConvGamma_Eta","MC_ConvGamma_Eta",2000,-2,2);
745                                 fMCList[iCut]->Add(hMCConvGammaEta[iCut]);
746                         }
747
748                         if(fDoMesonAnalysis){
749                                 hMCPi0Pt[iCut] = new TH1F("MC_Pi0_Pt","MC_Pi0_Pt",250,0,25);
750                                 hMCPi0Pt[iCut]->Sumw2();
751                                 fMCList[iCut]->Add(hMCPi0Pt[iCut]);
752                                 hMCPi0WOWeightPt[iCut] = new TH1F("MC_Pi0_WOWeights_Pt","MC_Pi0_WOWeights_Pt",250,0,25);
753                                 hMCPi0WOWeightPt[iCut]->Sumw2();
754                                 fMCList[iCut]->Add(hMCPi0WOWeightPt[iCut]);
755                                 
756                                 hMCEtaPt[iCut] = new TH1F("MC_Eta_Pt","MC_Eta_Pt",250,0,25);
757                                 hMCEtaPt[iCut]->Sumw2();
758                                 fMCList[iCut]->Add(hMCEtaPt[iCut]);
759                                 hMCEtaWOWeightPt[iCut] = new TH1F("MC_Eta_WOWeights_Pt","MC_Eta_WOWeights_Pt",250,0,25);
760                                 hMCEtaWOWeightPt[iCut]->Sumw2();
761                                 fMCList[iCut]->Add(hMCEtaWOWeightPt[iCut]);
762                                 hMCPi0InAccPt[iCut] = new TH1F("MC_Pi0InAcc_Pt","MC_Pi0InAcc_Pt",250,0,25);
763                                 hMCPi0InAccPt[iCut]->Sumw2();
764                                 fMCList[iCut]->Add(hMCPi0InAccPt[iCut]);
765                                 hMCEtaInAccPt[iCut] = new TH1F("MC_EtaInAcc_Pt","MC_EtaInAcc_Pt",250,0,25);
766                                 hMCEtaInAccPt[iCut]->Sumw2();
767                                 fMCList[iCut]->Add(hMCEtaInAccPt[iCut]);
768                                 if (fDoMesonQA > 0){
769                                         hMCPi0PtY[iCut] = new TH2F("MC_Pi0_Pt_Y","MC_Pi0_Pt_Y",150,0.03,15.,150,-1.5,1.5);
770                                         hMCPi0PtY[iCut]->Sumw2();
771                                         SetLogBinningXTH2(hMCPi0PtY[iCut]);
772                                         fMCList[iCut]->Add(hMCPi0PtY[iCut]);
773                                         hMCEtaPtY[iCut] = new TH2F("MC_Eta_Pt_Y","MC_Eta_Pt_Y",150,0.03,15.,150,-1.5,1.5);
774                                         hMCEtaPtY[iCut]->Sumw2();
775                                         SetLogBinningXTH2(hMCEtaPtY[iCut]);
776                                         fMCList[iCut]->Add(hMCEtaPtY[iCut]);
777                                         hMCPi0PtAlpha[iCut] = new TH2F("MC_Pi0_Pt_Alpha","MC_Pi0_Pt_Alpha",150,0.03,15.,100,0,1);
778                                         SetLogBinningXTH2(hMCPi0PtAlpha[iCut]);
779                                         fMCList[iCut]->Add(hMCPi0PtAlpha[iCut]);
780                                         hMCEtaPtAlpha[iCut] = new TH2F("MC_Eta_Pt_Alpha","MC_Eta_Pt_Alpha",150,0.03,15.,100,0,1);
781                                         SetLogBinningXTH2(hMCEtaPtAlpha[iCut]);
782                                         fMCList[iCut]->Add(hMCEtaPtAlpha[iCut]);
783
784                                         hMCK0sPt[iCut] = new TH1F("MC_K0s_Pt","MC_K0s_Pt",150,0,15);
785                                         hMCK0sPt[iCut]->Sumw2();
786                                         fMCList[iCut]->Add(hMCK0sPt[iCut]);
787                                         hMCK0sWOWeightPt[iCut] = new TH1F("MC_K0s_WOWeights_Pt","MC_K0s_WOWeights_Pt",150,0,15);
788                                         hMCK0sWOWeightPt[iCut]->Sumw2();
789                                         fMCList[iCut]->Add(hMCK0sWOWeightPt[iCut]);
790                                         hMCK0sPtY[iCut] = new TH2F("MC_K0s_Pt_Y","MC_K0s_Pt_Y",150,0.03,15.,150,-1.5,1.5);
791                                         hMCK0sPtY[iCut]->Sumw2();
792                                         SetLogBinningXTH2(hMCK0sPtY[iCut]);
793                                         fMCList[iCut]->Add(hMCK0sPtY[iCut]);
794                                         
795                                         hMCSecPi0Source[iCut] = new TH1F("MC_SecPi0_Source","MC_SecPi0_Source",5000,0.,5000);
796                                         fMCList[iCut]->Add(hMCSecPi0Source[iCut]);
797                                         hMCSecEtaSource[iCut] = new TH1F("MC_SecEta_Source","MC_SecEta_Source",5000,0,5000);
798                                         fMCList[iCut]->Add(hMCSecEtaSource[iCut]);
799                                         hMCSecPi0PtvsSource[iCut] = new TH2F("MC_SecPi0_Pt_Source","MC_SecPi0_Pt_Source",250,0.0,25.,16,-0.5,15.5);
800                                         hMCSecPi0PtvsSource[iCut]->Sumw2();
801                                         fMCList[iCut]->Add(hMCSecPi0PtvsSource[iCut]);
802                                         hMCSecEtaPt[iCut] = new TH1F("MC_SecEta_Pt","MC_SecEta_Pt",250,0,25);
803                                         hMCSecEtaPt[iCut]->Sumw2();
804                                         fMCList[iCut]->Add(hMCSecEtaPt[iCut]);
805                                 }
806
807                         }
808                         fTrueList[iCut] = new TList();
809                         fTrueList[iCut]->SetName(Form("%s_%s_%s True histograms",cutstringEvent.Data() ,cutstringPhoton.Data(),cutstringMeson.Data()));
810                         fTrueList[iCut]->SetOwner(kTRUE);
811                         fCutFolder[iCut]->Add(fTrueList[iCut]);
812
813                         hESDTrueConvGammaPt[iCut] = new TH1F("ESD_TrueConvGamma_Pt","ESD_TrueConvGamma_Pt",250,0,25);
814                         fTrueList[iCut]->Add(hESDTrueConvGammaPt[iCut]);
815
816                         hESDCombinatorialPt[iCut] = new TH2F("ESD_TrueCombinatorial_Pt","ESD_TrueCombinatorial_Pt",250,0,25,16,-0.5,15.5);
817                         hESDCombinatorialPt[iCut]->GetYaxis()->SetBinLabel( 1,"Elec+Elec");
818                         hESDCombinatorialPt[iCut]->GetYaxis()->SetBinLabel( 2,"Elec+Pion");
819                         hESDCombinatorialPt[iCut]->GetYaxis()->SetBinLabel( 3,"Elec+Kaon");
820                         hESDCombinatorialPt[iCut]->GetYaxis()->SetBinLabel( 4,"Elec+Proton");
821                         hESDCombinatorialPt[iCut]->GetYaxis()->SetBinLabel( 5,"Elec+Muon");
822                         hESDCombinatorialPt[iCut]->GetYaxis()->SetBinLabel( 6,"Pion+Pion");
823                         hESDCombinatorialPt[iCut]->GetYaxis()->SetBinLabel( 7,"Pion+Kaon");
824                         hESDCombinatorialPt[iCut]->GetYaxis()->SetBinLabel( 8,"Pion+Proton");
825                         hESDCombinatorialPt[iCut]->GetYaxis()->SetBinLabel( 9,"Pion+Muon");
826                         hESDCombinatorialPt[iCut]->GetYaxis()->SetBinLabel(10,"Kaon+Kaon");
827                         hESDCombinatorialPt[iCut]->GetYaxis()->SetBinLabel(11,"Kaon+Proton");
828                         hESDCombinatorialPt[iCut]->GetYaxis()->SetBinLabel(12,"Kaon+Muon");
829                         hESDCombinatorialPt[iCut]->GetYaxis()->SetBinLabel(13,"Proton+Proton");
830                         hESDCombinatorialPt[iCut]->GetYaxis()->SetBinLabel(14,"Proton+Muon");
831                         hESDCombinatorialPt[iCut]->GetYaxis()->SetBinLabel(15,"Muon+Muon");
832                         hESDCombinatorialPt[iCut]->GetYaxis()->SetBinLabel(16,"Rest");
833                         fTrueList[iCut]->Add(hESDCombinatorialPt[iCut]);
834                         hESDTruePrimaryConvGammaPt[iCut] = new TH1F("ESD_TruePrimaryConvGamma_Pt","ESD_TruePrimaryConvGamma_Pt",250,0,25);
835                         fTrueList[iCut]->Add(hESDTruePrimaryConvGammaPt[iCut]);
836                         hESDTrueSecondaryConvGammaPt[iCut] = new TH1F("ESD_TrueSecondaryConvGamma_Pt","ESD_TrueSecondaryConvGamma_Pt",250,0,25);
837                         fTrueList[iCut]->Add(hESDTrueSecondaryConvGammaPt[iCut]);
838
839                         hESDTrueSecondaryConvGammaFromXFromK0sPt[iCut]
840                                 = new TH1F("ESD_TrueSecondaryConvGammaFromXFromK0s_Pt", "ESD_TrueSecondaryConvGammaFromXFromK0s_Pt",250,0,25);
841                         fTrueList[iCut]->Add(hESDTrueSecondaryConvGammaFromXFromK0sPt[iCut]);
842                         hESDTrueSecondaryConvGammaFromXFromLambdaPt[iCut]
843                                 = new TH1F("ESD_TrueSecondaryConvGammaFromXFromLambda_Pt", "ESD_TrueSecondaryConvGammaFromXFromLambda_Pt",250,0,25);
844                         fTrueList[iCut]->Add(hESDTrueSecondaryConvGammaFromXFromLambdaPt[iCut]);
845
846                         hESDTrueDalitzPsiPairDeltaPhi[iCut]
847                                 = new TH2F("ESD_TrueDalitzPsiPairDeltaPhi_Pt", "ESD_TrueDalitzPsiPairDeltaPhi_Pt",100,-0.5,2,100,-0.5,0.5);
848                         fTrueList[iCut]->Add(hESDTrueDalitzPsiPairDeltaPhi[iCut]);
849                         
850                         hESDTrueGammaPsiPairDeltaPhi[iCut]
851                                 = new TH2F("ESD_TrueGammaPsiPairDeltaPhi_Pt", "ESD_TrueGammaPsiPairDeltaPhi_Pt",100,-0.5,2,100,-0.5,0.5);
852                         fTrueList[iCut]->Add(hESDTrueGammaPsiPairDeltaPhi[iCut]);
853
854                         hESDTruePrimaryConvGammaESDPtMCPt[iCut] = new TH2F("ESD_TruePrimaryConvGammaESD_PtMCPt", "ESD_TruePrimaryConvGammaESD_PtMCPt",250,0,25,250,0,25);
855                         fTrueList[iCut]->Add(hESDTruePrimaryConvGammaESDPtMCPt[iCut]);
856                         
857                         if(fDoMesonAnalysis){
858                                 hESDTrueMotherInvMassPt[iCut] = new TH2F("ESD_TrueMother_InvMass_Pt","ESD_TrueMother_InvMass_Pt",800,0,0.8,250,0,25);
859                                 fTrueList[iCut]->Add(hESDTrueMotherInvMassPt[iCut]);
860                                 hESDTruePrimaryMotherInvMassPt[iCut]
861                                 = new TH2F("ESD_TruePrimaryMother_InvMass_Pt", "ESD_TruePrimaryMother_InvMass_Pt", 800,0,0.8,250,0,25);
862                                 hESDTruePrimaryMotherInvMassPt[iCut]->Sumw2();
863                                 fTrueList[iCut]->Add(hESDTruePrimaryMotherInvMassPt[iCut]);
864                                 hESDTruePrimaryMotherW0WeightingInvMassPt[iCut]
865                                 = new TH2F("ESD_TruePrimaryMotherW0Weights_InvMass_Pt", "ESD_TruePrimaryMotherW0Weights_InvMass_Pt", 800,0,0.8,250,0,25);
866                                 hESDTruePrimaryMotherW0WeightingInvMassPt[iCut]->Sumw2();
867                                 fTrueList[iCut]->Add(hESDTruePrimaryMotherW0WeightingInvMassPt[iCut]);
868                                 pESDTruePrimaryMotherWeightsInvMassPt[iCut]
869                                 = new TProfile2D("ESD_TruePrimaryMotherWeights_InvMass_Pt", "ESD_TruePrimaryMotherWeights_InvMass_Pt", 800,0,0.8,250,0,25);
870                                 pESDTruePrimaryMotherWeightsInvMassPt[iCut]->Sumw2();
871                                 fTrueList[iCut]->Add(pESDTruePrimaryMotherWeightsInvMassPt[iCut]);
872                                 hESDTrueSecondaryMotherInvMassPt[iCut]
873                                 = new TH2F("ESD_TrueSecondaryMother_InvMass_Pt", "ESD_TrueSecondaryMother_InvMass_Pt", 800,0,0.8,250,0,25);
874                                 hESDTrueSecondaryMotherInvMassPt[iCut]->Sumw2();
875                                 fTrueList[iCut]->Add(hESDTrueSecondaryMotherInvMassPt[iCut]);
876                                 hESDTrueSecondaryMotherFromK0sInvMassPt[iCut]
877                                 = new TH2F("ESD_TrueSecondaryMotherFromK0s_InvMass_Pt","ESD_TrueSecondaryMotherFromK0s_InvMass_Pt",800,0,0.8,250,0,25);
878                                 hESDTrueSecondaryMotherFromK0sInvMassPt[iCut]->Sumw2();
879                                 fTrueList[iCut]->Add(hESDTrueSecondaryMotherFromK0sInvMassPt[iCut]);
880                                 hESDTrueSecondaryMotherFromEtaInvMassPt[iCut]
881                                 = new TH2F("ESD_TrueSecondaryMotherFromEta_InvMass_Pt","ESD_TrueSecondaryMotherFromEta_InvMass_Pt",800,0,0.8,250,0,25);
882                                 fTrueList[iCut]->Add(hESDTrueSecondaryMotherFromEtaInvMassPt[iCut]);
883                                 hESDTrueSecondaryMotherFromLambdaInvMassPt[iCut]
884                                 = new TH2F("ESD_TrueSecondaryMotherFromLambda_InvMass_Pt","ESD_TrueSecondaryMotherFromLambda_InvMass_Pt",800,0,0.8,250,0,25);
885                                 fTrueList[iCut]->Add(hESDTrueSecondaryMotherFromLambdaInvMassPt[iCut]);
886                                 hESDTrueMotherDalitzInvMassPt[iCut] = new TH2F("ESD_TrueDalitz_InvMass_Pt","ESD_TrueDalitz_InvMass_Pt",800,0,0.8,250,0,25);
887                                 fTrueList[iCut]->Add(hESDTrueMotherDalitzInvMassPt[iCut]);         
888                                 if (fDoMesonQA > 0){
889                                         hESDTruePrimaryPi0MCPtResolPt[iCut] = new TH2F("ESD_TruePrimaryPi0_MCPt_ResolPt","ESD_TruePrimaryPi0_ResolPt_MCPt",500,0.03,25,1000,-1.,1.);
890                                         hESDTruePrimaryPi0MCPtResolPt[iCut]->Sumw2();
891                                         SetLogBinningXTH2(hESDTruePrimaryPi0MCPtResolPt[iCut]);
892                                         fTrueList[iCut]->Add(hESDTruePrimaryPi0MCPtResolPt[iCut]);
893                                         hESDTruePrimaryEtaMCPtResolPt[iCut]  = new TH2F("ESD_TruePrimaryEta_MCPt_ResolPt","ESD_TruePrimaryEta_ResolPt_MCPt",500,0.03,25,1000,-1.,1.);
894                                         hESDTruePrimaryEtaMCPtResolPt[iCut]->Sumw2();
895                                         SetLogBinningXTH2(hESDTruePrimaryEtaMCPtResolPt[iCut]);
896                                         fTrueList[iCut]->Add(hESDTruePrimaryEtaMCPtResolPt[iCut]);
897                                         hESDTrueBckGGInvMassPt[iCut] = new TH2F("ESD_TrueBckGG_InvMass_Pt","ESD_TrueBckGG_InvMass_Pt",800,0,0.8,250,0,25);
898                                         fTrueList[iCut]->Add(hESDTrueBckGGInvMassPt[iCut]);
899                                         hESDTrueBckContInvMassPt[iCut] = new TH2F("ESD_TrueBckCont_InvMass_Pt","ESD_TrueBckCont_InvMass_Pt",800,0,0.8,250,0,25);
900                                         fTrueList[iCut]->Add(hESDTrueBckContInvMassPt[iCut]);
901                                         hESDTrueK0sWithPi0DaughterMCPt[iCut] = new TH1F("ESD_TrueK0sWithPi0Daughter_MCPt","ESD_TrueK0sWithPi0Daughter_MCPt",250,0,25);
902                                         fTrueList[iCut]->Add(hESDTrueK0sWithPi0DaughterMCPt[iCut]);
903                                         hESDTrueEtaWithPi0DaughterMCPt[iCut] = new TH1F("ESD_TrueEtaWithPi0Daughter_MCPt","ESD_TrueEtaWithPi0Daughter_MCPt",250,0,25);
904                                         fTrueList[iCut]->Add(hESDTrueEtaWithPi0DaughterMCPt[iCut]);
905                                         hESDTrueLambdaWithPi0DaughterMCPt[iCut] = new TH1F("ESD_TrueLambdaWithPi0Daughter_MCPt","ESD_TrueLambdaWithPi0Daughter_MCPt",250,0,25);
906                                         fTrueList[iCut]->Add(hESDTrueLambdaWithPi0DaughterMCPt[iCut]);
907
908                                         hESDTruePi0PtY[iCut] = new TH2F("ESD_TruePi0_Pt_Y","ESD_TruePi0_Pt_Y",150,0.03,15.,150,-1.5,1.5);
909                                         SetLogBinningXTH2(hESDTruePi0PtY[iCut]);
910                                         fTrueList[iCut]->Add(hESDTruePi0PtY[iCut]);
911                                         hESDTrueEtaPtY[iCut] = new TH2F("ESD_TrueEta_Pt_Y","ESD_TrueEta_Pt_Y",150,0.03,15.,150,-1.5,1.5);
912                                         SetLogBinningXTH2(hESDTrueEtaPtY[iCut]);
913                                         fTrueList[iCut]->Add(hESDTrueEtaPtY[iCut]);
914                                         hESDTruePi0PtAlpha[iCut] = new TH2F("ESD_TruePi0_Pt_Alpha","ESD_TruePi0_Pt_Alpha",150,0.03,15.,100,0,1);
915                                         SetLogBinningXTH2(hESDTruePi0PtAlpha[iCut]);
916                                         fTrueList[iCut]->Add(hESDTruePi0PtAlpha[iCut]);
917                                         hESDTrueEtaPtAlpha[iCut] = new TH2F("ESD_TrueEta_Pt_Alpha","ESD_TrueEta_Pt_Alpha",150,0.03,15.,100,0,1);
918                                         SetLogBinningXTH2(hESDTrueEtaPtAlpha[iCut]);
919                                         fTrueList[iCut]->Add(hESDTrueEtaPtAlpha[iCut]);
920                                         
921                                         hESDTruePi0PtOpenAngle[iCut] = new TH2F("ESD_TruePi0_Pt_OpenAngle","ESD_TruePi0_Pt_OpenAngle",150,0.03,15.,200,0,2*TMath::Pi());            
922                                         SetLogBinningXTH2(hESDTruePi0PtOpenAngle[iCut]);
923                                         fTrueList[iCut]->Add(hESDTruePi0PtOpenAngle[iCut]);
924                                         hESDTrueEtaPtOpenAngle[iCut] = new TH2F("ESD_TrueEta_Pt_OpenAngle","ESD_TrueEta_Pt_OpenAngle",150,0.03,15.,200,0,2*TMath::Pi());            
925                                         SetLogBinningXTH2(hESDTrueEtaPtOpenAngle[iCut]);
926                                         fTrueList[iCut]->Add(hESDTrueEtaPtOpenAngle[iCut]);
927
928                                         hESDTrueConvGammaEta[iCut] = new TH1F("ESD_TrueConvGamma_Eta","ESD_TrueConvGamma_Eta",2000,-2,2);
929                                         fTrueList[iCut]->Add(hESDTrueConvGammaEta[iCut]);
930
931                                 }
932                         }
933                 }
934         }
935
936         fV0Reader=(AliV0ReaderV1*)AliAnalysisManager::GetAnalysisManager()->GetTask("V0ReaderV1");
937         if(!fV0Reader){printf("Error: No V0 Reader");return;} // GetV0Reader
938
939         if(fV0Reader)
940                 if((AliConvEventCuts*)fV0Reader->GetEventCuts())
941                         if(((AliConvEventCuts*)fV0Reader->GetEventCuts())->GetCutHistograms())
942                                 fOutputContainer->Add(((AliConvEventCuts*)fV0Reader->GetEventCuts())->GetCutHistograms());
943
944         if(fV0Reader)
945                 if((AliConversionPhotonCuts*)fV0Reader->GetConversionCuts())
946                         if(((AliConversionPhotonCuts*)fV0Reader->GetConversionCuts())->GetCutHistograms())
947                                 fOutputContainer->Add(((AliConversionPhotonCuts*)fV0Reader->GetConversionCuts())->GetCutHistograms());
948
949         for(Int_t iCut = 0; iCut<fnCuts;iCut++){
950                 if(!((AliConvEventCuts*)fEventCutArray->At(iCut))) continue;
951                 if(((AliConvEventCuts*)fEventCutArray->At(iCut))->GetCutHistograms()){
952                         fCutFolder[iCut]->Add(((AliConvEventCuts*)fEventCutArray->At(iCut))->GetCutHistograms());
953                 }
954                 if(!((AliConversionPhotonCuts*)fCutArray->At(iCut))) continue;
955                 if(((AliConversionPhotonCuts*)fCutArray->At(iCut))->GetCutHistograms()){
956                         fCutFolder[iCut]->Add(((AliConversionPhotonCuts*)fCutArray->At(iCut))->GetCutHistograms());
957                 }
958                 if(fDoMesonAnalysis){
959                         if(!((AliConversionMesonCuts*)fMesonCutArray->At(iCut))) continue;
960                         if(((AliConversionMesonCuts*)fMesonCutArray->At(iCut))->GetCutHistograms()){
961                                 fCutFolder[iCut]->Add(((AliConversionMesonCuts*)fMesonCutArray->At(iCut))->GetCutHistograms());
962                         }
963                 }
964         }
965         PostData(1, fOutputContainer);
966 }
967 //_____________________________________________________________________________
968 Bool_t AliAnalysisTaskGammaConvV1::Notify()
969 {
970         for(Int_t iCut = 0; iCut<fnCuts;iCut++){
971                 if(!((AliConvEventCuts*)fEventCutArray->At(iCut))->GetDoEtaShift()){
972                         hEtaShift[iCut]->Fill(0.,(((AliConvEventCuts*)fEventCutArray->At(iCut))->GetEtaShift()));
973                         continue; // No Eta Shift requested, continue
974                 }
975                 if(((AliConvEventCuts*)fEventCutArray->At(iCut))->GetEtaShift() == 0.0){ // Eta Shift requested but not set, get shift automatically
976                         ((AliConvEventCuts*)fEventCutArray->At(iCut))->GetCorrectEtaShiftFromPeriod(fV0Reader->GetPeriodName());
977                         hEtaShift[iCut]->Fill(0.,(((AliConvEventCuts*)fEventCutArray->At(iCut))->GetEtaShift()));
978                         ((AliConvEventCuts*)fEventCutArray->At(iCut))->DoEtaShift(kFALSE); // Eta Shift Set, make sure that it is called only once   
979                         continue;
980                 }
981                 else{
982                         printf(" Gamma Conversion Task %s :: Eta Shift Manually Set to %f \n\n",
983                                         (((AliConvEventCuts*)fEventCutArray->At(iCut))->GetCutNumber()).Data(),((AliConvEventCuts*)fEventCutArray->At(iCut))->GetEtaShift());
984                         hEtaShift[iCut]->Fill(0.,(((AliConvEventCuts*)fEventCutArray->At(iCut))->GetEtaShift()));
985                         ((AliConvEventCuts*)fEventCutArray->At(iCut))->DoEtaShift(kFALSE); // Eta Shift Set, make sure that it is called only once   
986                 }
987         }
988         return kTRUE;
989 }
990 //_____________________________________________________________________________
991 void AliAnalysisTaskGammaConvV1::UserExec(Option_t *)
992 {
993         //
994         // Called for each event
995         //
996         Int_t eventQuality = ((AliConvEventCuts*)fV0Reader->GetEventCuts())->GetEventQuality();
997         if(eventQuality == 2 || eventQuality == 3){// Event Not Accepted due to MC event missing or wrong trigger for V0ReaderV1
998                 for(Int_t iCut = 0; iCut<fnCuts; iCut++){
999                         hNEvents[iCut]->Fill(eventQuality);
1000                 }
1001                 return;
1002         }
1003
1004         if(fIsMC) fMCEvent = MCEvent();
1005         if(fMCEvent == NULL) fIsMC = kFALSE;
1006         
1007         fInputEvent = InputEvent();
1008
1009         if(fIsMC && fInputEvent->IsA()==AliESDEvent::Class()){
1010                 fMCStack = fMCEvent->Stack();
1011                 if(fMCStack == NULL) fIsMC = kFALSE;
1012         }
1013
1014         fReaderGammas = fV0Reader->GetReconstructedGammas(); // Gammas from default Cut
1015         
1016         // ------------------- BeginEvent ----------------------------
1017
1018         AliEventplane *EventPlane = fInputEvent->GetEventplane();
1019         if(fIsHeavyIon ==1)fEventPlaneAngle = EventPlane->GetEventplane("V0",fInputEvent,2);
1020         else fEventPlaneAngle=0.0;
1021         
1022         if(fIsMC && fInputEvent->IsA()==AliAODEvent::Class() && !(fV0Reader->AreAODsRelabeled())){
1023                 RelabelAODPhotonCandidates(kTRUE);    // In case of AODMC relabeling MC
1024                 fV0Reader->RelabelAODs(kTRUE);
1025         }
1026         for(Int_t iCut = 0; iCut<fnCuts; iCut++){
1027                 fiCut = iCut;
1028                 Int_t eventNotAccepted =
1029                         ((AliConvEventCuts*)fEventCutArray->At(iCut))
1030                         ->IsEventAcceptedByCut(fV0Reader->GetEventCuts(),fInputEvent,fMCEvent,fIsHeavyIon);
1031                 if(eventNotAccepted){
1032                         // cout << "event rejected due to wrong trigger: " <<eventNotAccepted << endl;
1033                         hNEvents[iCut]->Fill(eventNotAccepted); // Check Centrality, PileUp, SDD and V0AND --> Not Accepted => eventQuality = 1
1034                         continue;
1035                 }
1036
1037                 if(eventQuality != 0){// Event Not Accepted
1038                         // cout << "event rejected due to: " <<eventQuality << endl;
1039                         hNEvents[iCut]->Fill(eventQuality);
1040                         continue;
1041                 }
1042
1043                 hNEvents[iCut]->Fill(eventQuality); // Should be 0 here
1044                 hNGoodESDTracks[iCut]->Fill(fV0Reader->GetNumberOfPrimaryTracks());
1045                 if(((AliConvEventCuts*)fEventCutArray->At(iCut))->IsHeavyIon() == 2)    hNV0Tracks[iCut]->Fill(fInputEvent->GetVZEROData()->GetMTotV0A());
1046                 else hNV0Tracks[iCut]->Fill(fInputEvent->GetVZEROData()->GetMTotV0A()+fInputEvent->GetVZEROData()->GetMTotV0C());
1047
1048                 if(fIsMC){
1049                         // Process MC Particle
1050                         if(((AliConvEventCuts*)fEventCutArray->At(iCut))->GetSignalRejection() != 0){
1051                                 if(fInputEvent->IsA()==AliESDEvent::Class()){
1052                                 ((AliConvEventCuts*)fEventCutArray->At(iCut))->GetNotRejectedParticles(((AliConvEventCuts*)fEventCutArray->At(iCut))->GetSignalRejection(),
1053                                                                                                                                                                         ((AliConvEventCuts*)fEventCutArray->At(iCut))->GetAcceptedHeader(),
1054                                                                                                                                                                         fMCEvent);
1055                                 }
1056                                 else if(fInputEvent->IsA()==AliAODEvent::Class()){
1057                                 ((AliConvEventCuts*)fEventCutArray->At(iCut))->GetNotRejectedParticles(((AliConvEventCuts*)fEventCutArray->At(iCut))->GetSignalRejection(),
1058                                                                                                                                                                         ((AliConvEventCuts*)fEventCutArray->At(iCut))->GetAcceptedHeader(),
1059                                                                                                                                                                         fInputEvent);
1060                                 }
1061
1062                                 if(((AliConvEventCuts*)fEventCutArray->At(iCut))->GetAcceptedHeader()){
1063                                 for(Int_t i = 0;i<(((AliConvEventCuts*)fEventCutArray->At(iCut))->GetAcceptedHeader())->GetEntries();i++){
1064                                         TString nameBin= hMCHeaders[iCut]->GetXaxis()->GetBinLabel(i+1);
1065                                         if (nameBin.CompareTo("")== 0){
1066                                                 TString nameHeader = ((TObjString*)((TList*)((AliConvEventCuts*)fEventCutArray->At(iCut))
1067                                                                                                                         ->GetAcceptedHeader())->At(i))->GetString();
1068 //                                              cout << nameHeader << endl;
1069                                                 hMCHeaders[iCut]->GetXaxis()->SetBinLabel(i+1,nameHeader.Data());
1070                                         }
1071                                 }
1072                                 }
1073                         }
1074                 }
1075                 if(fIsMC){
1076                         if(fInputEvent->IsA()==AliESDEvent::Class())
1077                                 ProcessMCParticles();
1078                         if(fInputEvent->IsA()==AliAODEvent::Class())
1079                                 ProcessAODMCParticles();
1080                 }
1081
1082                 ProcessPhotonCandidates(); // Process this cuts gammas
1083
1084                 hNGammaCandidates[iCut]->Fill(fGammaCandidates->GetEntries());
1085                 hNGoodESDTracksVsNGammaCanditates[iCut]->Fill(fV0Reader->GetNumberOfPrimaryTracks(),fGammaCandidates->GetEntries());
1086                 if(fDoMesonAnalysis){ // Meson Analysis
1087                         if(((AliConversionMesonCuts*)fMesonCutArray->At(iCut))->UseMCPSmearing() && fIsMC){
1088                                 fUnsmearedPx = new Double_t[fGammaCandidates->GetEntries()]; // Store unsmeared Momenta
1089                                 fUnsmearedPy = new Double_t[fGammaCandidates->GetEntries()];
1090                                 fUnsmearedPz = new Double_t[fGammaCandidates->GetEntries()];
1091                                 fUnsmearedE =  new Double_t[fGammaCandidates->GetEntries()];
1092
1093                                 for(Int_t gamma=0;gamma<fGammaCandidates->GetEntries();gamma++){ // Smear the AODPhotons in MC
1094                                 fUnsmearedPx[gamma] = ((AliAODConversionPhoton*)fGammaCandidates->At(gamma))->Px();
1095                                 fUnsmearedPy[gamma] = ((AliAODConversionPhoton*)fGammaCandidates->At(gamma))->Py();
1096                                 fUnsmearedPz[gamma] = ((AliAODConversionPhoton*)fGammaCandidates->At(gamma))->Pz();
1097                                 fUnsmearedE[gamma] =  ((AliAODConversionPhoton*)fGammaCandidates->At(gamma))->E();
1098                                 ((AliConversionMesonCuts*)fMesonCutArray->At(iCut))->SmearParticle(dynamic_cast<AliAODConversionPhoton*>(fGammaCandidates->At(gamma)));
1099                                 }
1100                         }
1101
1102                         CalculatePi0Candidates(); // Combine Gammas
1103                         if(((AliConversionMesonCuts*)fMesonCutArray->At(iCut))->DoBGCalculation()){
1104                                 if(((AliConversionMesonCuts*)fMesonCutArray->At(iCut))->BackgroundHandlerType() == 0){
1105                                 CalculateBackground(); // Combinatorial Background
1106                                 UpdateEventByEventData(); // Store Event for mixed Events
1107                                 }
1108                                 else{
1109                                 CalculateBackgroundRP(); // Combinatorial Background
1110                                 fBGHandlerRP[iCut]->AddEvent(fGammaCandidates,fInputEvent); // Store Event for mixed Events
1111                                 }
1112                         }
1113                         if(((AliConversionMesonCuts*)fMesonCutArray->At(iCut))->UseMCPSmearing() && fIsMC){
1114                                 for(Int_t gamma=0;gamma<fGammaCandidates->GetEntries();gamma++){ // Smear the AODPhotons in MC
1115                                 ((AliAODConversionPhoton*)fGammaCandidates->At(gamma))->SetPx(fUnsmearedPx[gamma]); // Reset Unsmeared Momenta
1116                                 ((AliAODConversionPhoton*)fGammaCandidates->At(gamma))->SetPy(fUnsmearedPy[gamma]);
1117                                 ((AliAODConversionPhoton*)fGammaCandidates->At(gamma))->SetPz(fUnsmearedPz[gamma]);
1118                                 ((AliAODConversionPhoton*)fGammaCandidates->At(gamma))->SetE(fUnsmearedE[gamma]);
1119                                 }
1120                                 delete[] fUnsmearedPx; fUnsmearedPx = 0x0;
1121                                 delete[] fUnsmearedPy; fUnsmearedPy = 0x0;
1122                                 delete[] fUnsmearedPz; fUnsmearedPz = 0x0;
1123                                 delete[] fUnsmearedE;  fUnsmearedE  = 0x0;
1124                         }
1125                 }
1126                 fGammaCandidates->Clear(); // delete this cuts good gammas
1127         }
1128
1129         if(fIsMC && fInputEvent->IsA()==AliAODEvent::Class() && !(fV0Reader->AreAODsRelabeled())){
1130                 RelabelAODPhotonCandidates(kFALSE); // Back to ESDMC Label
1131                 fV0Reader->RelabelAODs(kFALSE);
1132         }
1133
1134         PostData(1, fOutputContainer);
1135 }
1136 //________________________________________________________________________
1137 void AliAnalysisTaskGammaConvV1::ProcessPhotonCandidates()
1138 {
1139         Int_t nV0 = 0;
1140         TList *GammaCandidatesStepOne = new TList();
1141         TList *GammaCandidatesStepTwo = new TList();
1142         // Loop over Photon Candidates allocated by ReaderV1
1143         for(Int_t i = 0; i < fReaderGammas->GetEntriesFast(); i++){
1144                 AliAODConversionPhoton* PhotonCandidate = (AliAODConversionPhoton*) fReaderGammas->At(i);
1145                 if(!PhotonCandidate) continue;
1146                 fIsFromMBHeader = kTRUE;
1147                 if(fIsMC && ((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetSignalRejection() != 0){
1148                         Int_t isPosFromMBHeader
1149                                 = ((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsParticleFromBGEvent(PhotonCandidate->GetMCLabelPositive(), fMCStack, fInputEvent);
1150                         if(isPosFromMBHeader == 0 && ((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetSignalRejection() != 3) continue;
1151                         Int_t isNegFromMBHeader
1152                                 = ((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsParticleFromBGEvent(PhotonCandidate->GetMCLabelNegative(), fMCStack, fInputEvent);
1153                         if(isNegFromMBHeader == 0 && ((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetSignalRejection() != 3) continue;
1154
1155                         if( (isNegFromMBHeader+isPosFromMBHeader) != 4) fIsFromMBHeader = kFALSE;
1156                 }
1157
1158                 if(!((AliConversionPhotonCuts*)fCutArray->At(fiCut))->PhotonIsSelected(PhotonCandidate,fInputEvent)) continue;
1159                 if(!((AliConversionPhotonCuts*)fCutArray->At(fiCut))->InPlaneOutOfPlaneCut(PhotonCandidate->GetPhotonPhi(),fEventPlaneAngle)) continue;
1160                 if(!((AliConversionPhotonCuts*)fCutArray->At(fiCut))->UseElecSharingCut() &&
1161                         !((AliConversionPhotonCuts*)fCutArray->At(fiCut))->UseToCloseV0sCut()){
1162                         fGammaCandidates->Add(PhotonCandidate); // if no second loop is required add to events good gammas
1163
1164                         if(fIsFromMBHeader){
1165                                 hESDConvGammaPt[fiCut]->Fill(PhotonCandidate->Pt());
1166                                 if (fDoPhotonQA > 0){
1167                                 hESDConvGammaR[fiCut]->Fill(PhotonCandidate->GetConversionRadius());
1168                                 hESDConvGammaEta[fiCut]->Fill(PhotonCandidate->Eta());
1169                                 }   
1170                         }
1171                         if(fIsMC){
1172                                 if(fInputEvent->IsA()==AliESDEvent::Class())
1173                                 ProcessTruePhotonCandidates(PhotonCandidate);
1174                                 if(fInputEvent->IsA()==AliAODEvent::Class())
1175                                 ProcessTruePhotonCandidatesAOD(PhotonCandidate);
1176                         }
1177                         if (fIsFromMBHeader && fDoPhotonQA == 2){
1178                                 if (fIsHeavyIon == 1 && PhotonCandidate->Pt() > 0.399 && PhotonCandidate->Pt() < 12.){
1179                                 fPtGamma = PhotonCandidate->Pt();
1180                                 fDCAzPhoton = PhotonCandidate->GetDCAzToPrimVtx();
1181                                 fRConvPhoton = PhotonCandidate->GetConversionRadius();
1182                                 fEtaPhoton = PhotonCandidate->GetPhotonEta();
1183                                 iCatPhoton = PhotonCandidate->GetPhotonQuality();
1184                                 tESDConvGammaPtDcazCat[fiCut]->Fill();
1185                                 } else if ( PhotonCandidate->Pt() > 0.299 && PhotonCandidate->Pt() < 16.){
1186                                 fPtGamma = PhotonCandidate->Pt();
1187                                 fDCAzPhoton = PhotonCandidate->GetDCAzToPrimVtx();
1188                                 fRConvPhoton = PhotonCandidate->GetConversionRadius();
1189                                 fEtaPhoton = PhotonCandidate->GetPhotonEta();
1190                                 iCatPhoton = PhotonCandidate->GetPhotonQuality();
1191                                 tESDConvGammaPtDcazCat[fiCut]->Fill();
1192                                 }   
1193                         }   
1194                 } else if(((AliConversionPhotonCuts*)fCutArray->At(fiCut))->UseElecSharingCut()){ // if Shared Electron cut is enabled, Fill array, add to step one
1195                         ((AliConversionPhotonCuts*)fCutArray->At(fiCut))->FillElectonLabelArray(PhotonCandidate,nV0);
1196                         nV0++;
1197                         GammaCandidatesStepOne->Add(PhotonCandidate);
1198                 } else if(!((AliConversionPhotonCuts*)fCutArray->At(fiCut))->UseElecSharingCut() &&
1199                                 ((AliConversionPhotonCuts*)fCutArray->At(fiCut))->UseToCloseV0sCut()){ // shared electron is disabled, step one not needed -> step two
1200                         GammaCandidatesStepTwo->Add(PhotonCandidate);
1201                 }
1202         }
1203         if(((AliConversionPhotonCuts*)fCutArray->At(fiCut))->UseElecSharingCut()){
1204                 for(Int_t i = 0;i<GammaCandidatesStepOne->GetEntries();i++){
1205                         AliAODConversionPhoton *PhotonCandidate= (AliAODConversionPhoton*) GammaCandidatesStepOne->At(i);
1206                         if(!PhotonCandidate) continue;
1207                         fIsFromMBHeader = kTRUE;
1208                         if(fMCStack && ((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetSignalRejection() != 0){
1209                                 Int_t isPosFromMBHeader
1210                                 = ((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsParticleFromBGEvent(PhotonCandidate->GetMCLabelPositive(), fMCStack, fInputEvent);
1211                                 Int_t isNegFromMBHeader
1212                                 = ((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsParticleFromBGEvent(PhotonCandidate->GetMCLabelNegative(), fMCStack, fInputEvent);
1213                                 if( (isNegFromMBHeader+isPosFromMBHeader) != 4) fIsFromMBHeader = kFALSE;
1214                         }
1215                         if(!((AliConversionPhotonCuts*)fCutArray->At(fiCut))->RejectSharedElectronV0s(PhotonCandidate,i,GammaCandidatesStepOne->GetEntries())) continue;
1216                         if(!((AliConversionPhotonCuts*)fCutArray->At(fiCut))->UseToCloseV0sCut()){ // To Colse v0s cut diabled, step two not needed
1217                                 fGammaCandidates->Add(PhotonCandidate);
1218                                 if(fIsFromMBHeader){
1219                                 hESDConvGammaPt[fiCut]->Fill(PhotonCandidate->Pt());
1220                                 if (fDoPhotonQA > 0){
1221                                         hESDConvGammaR[fiCut]->Fill(PhotonCandidate->GetConversionRadius());
1222                                         hESDConvGammaEta[fiCut]->Fill(PhotonCandidate->Eta());
1223                                 }
1224                                 }
1225                         }
1226                         if(fIsMC){
1227                                 if(fInputEvent->IsA()==AliESDEvent::Class())
1228                                 ProcessTruePhotonCandidates(PhotonCandidate);
1229                                 if(fInputEvent->IsA()==AliAODEvent::Class())
1230                                 ProcessTruePhotonCandidatesAOD(PhotonCandidate);
1231                         } else GammaCandidatesStepTwo->Add(PhotonCandidate); // Close v0s cut enabled -> add to list two
1232                         
1233                         if (fIsFromMBHeader && fDoPhotonQA == 2){
1234                                 if (fIsHeavyIon ==1 && PhotonCandidate->Pt() > 0.399 && PhotonCandidate->Pt() < 12.){
1235                                 fPtGamma = PhotonCandidate->Pt();
1236                                 fDCAzPhoton = PhotonCandidate->GetDCAzToPrimVtx();
1237                                 fRConvPhoton = PhotonCandidate->GetConversionRadius();
1238                                 fEtaPhoton = PhotonCandidate->GetPhotonEta();
1239                                 iCatPhoton = PhotonCandidate->GetPhotonQuality();
1240                                 tESDConvGammaPtDcazCat[fiCut]->Fill();
1241                                 } else if ( PhotonCandidate->Pt() > 0.299 && PhotonCandidate->Pt() < 16.){
1242                                 fPtGamma = PhotonCandidate->Pt();
1243                                 fDCAzPhoton = PhotonCandidate->GetDCAzToPrimVtx();
1244                                 fRConvPhoton = PhotonCandidate->GetConversionRadius();
1245                                 fEtaPhoton = PhotonCandidate->GetPhotonEta();
1246                                 iCatPhoton = PhotonCandidate->GetPhotonQuality();
1247                                 tESDConvGammaPtDcazCat[fiCut]->Fill();
1248                                 }
1249                         }
1250                 }
1251         }
1252         if(((AliConversionPhotonCuts*)fCutArray->At(fiCut))->UseToCloseV0sCut()){
1253                 for(Int_t i = 0;i<GammaCandidatesStepTwo->GetEntries();i++){
1254                         AliAODConversionPhoton* PhotonCandidate = (AliAODConversionPhoton*) GammaCandidatesStepTwo->At(i);
1255                         if(!PhotonCandidate) continue;
1256                         fIsFromMBHeader = kTRUE;
1257                         if(fMCStack && ((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetSignalRejection() != 0){
1258                                 Int_t isPosFromMBHeader
1259                                 = ((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsParticleFromBGEvent(PhotonCandidate->GetMCLabelPositive(), fMCStack, fInputEvent);
1260                                 Int_t isNegFromMBHeader
1261                                 = ((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsParticleFromBGEvent(PhotonCandidate->GetMCLabelNegative(), fMCStack, fInputEvent);
1262                                 if( (isNegFromMBHeader+isPosFromMBHeader) != 4) fIsFromMBHeader = kFALSE;
1263                         }
1264                         if(!((AliConversionPhotonCuts*)fCutArray->At(fiCut))->RejectToCloseV0s(PhotonCandidate,GammaCandidatesStepTwo,i)) continue;
1265                         fGammaCandidates->Add(PhotonCandidate); // Add gamma to current cut TList
1266                         if(fIsFromMBHeader){
1267                                 hESDConvGammaPt[fiCut]->Fill(PhotonCandidate->Pt());
1268                                 if (fDoPhotonQA > 0){
1269                                 hESDConvGammaR[fiCut]->Fill(PhotonCandidate->GetConversionRadius());
1270                                 hESDConvGammaEta[fiCut]->Fill(PhotonCandidate->Eta());
1271                                 }
1272                         }
1273                         if(fIsMC){
1274                                 if(fInputEvent->IsA()==AliESDEvent::Class())
1275                                 ProcessTruePhotonCandidates(PhotonCandidate);
1276                                 if(fInputEvent->IsA()==AliAODEvent::Class())
1277                                 ProcessTruePhotonCandidatesAOD(PhotonCandidate);
1278                         }
1279                         if (fIsFromMBHeader){
1280                         if (fIsHeavyIon == 1 && PhotonCandidate->Pt() > 0.399 && PhotonCandidate->Pt() < 12.){
1281                                 fPtGamma = PhotonCandidate->Pt();
1282                                 fDCAzPhoton = PhotonCandidate->GetDCAzToPrimVtx();
1283                                 fRConvPhoton = PhotonCandidate->GetConversionRadius();
1284                                 fEtaPhoton = PhotonCandidate->GetPhotonEta();
1285                                 iCatPhoton = PhotonCandidate->GetPhotonQuality();
1286                                 tESDConvGammaPtDcazCat[fiCut]->Fill();
1287                                 } else if ( PhotonCandidate->Pt() > 0.299 && PhotonCandidate->Pt() < 16.){
1288                                 fPtGamma = PhotonCandidate->Pt();
1289                                 fDCAzPhoton = PhotonCandidate->GetDCAzToPrimVtx();
1290                                 fRConvPhoton = PhotonCandidate->GetConversionRadius();
1291                                 fEtaPhoton = PhotonCandidate->GetPhotonEta();
1292                                 iCatPhoton = PhotonCandidate->GetPhotonQuality();
1293                                 tESDConvGammaPtDcazCat[fiCut]->Fill();
1294                                 }
1295                         }
1296                 }
1297         }
1298
1299         delete GammaCandidatesStepOne;
1300         GammaCandidatesStepOne = 0x0;
1301         delete GammaCandidatesStepTwo;
1302         GammaCandidatesStepTwo = 0x0;
1303
1304 }
1305 //________________________________________________________________________
1306 void AliAnalysisTaskGammaConvV1::ProcessTruePhotonCandidatesAOD(AliAODConversionPhoton *TruePhotonCandidate)
1307 {
1308         Double_t magField = fInputEvent->GetMagneticField();
1309         if( magField  < 0.0 ){
1310                 magField =  1.0;
1311         }
1312         else {
1313                 magField =  -1.0;
1314         }
1315
1316         TClonesArray *AODMCTrackArray = dynamic_cast<TClonesArray*>(fInputEvent->FindListObject(AliAODMCParticle::StdBranchName()));
1317         if (AODMCTrackArray != NULL && TruePhotonCandidate != NULL){
1318
1319                 AliAODMCParticle *posDaughter = (AliAODMCParticle*) AODMCTrackArray->At(TruePhotonCandidate->GetMCLabelPositive());
1320                 AliAODMCParticle *negDaughter = (AliAODMCParticle*) AODMCTrackArray->At(TruePhotonCandidate->GetMCLabelNegative());
1321                 iPhotonMCInfo = 0;
1322                 
1323                 if(posDaughter == NULL || negDaughter == NULL) return; // One particle does not exist
1324                 Int_t pdgCode[2] = {abs(posDaughter->GetPdgCode()),abs(negDaughter->GetPdgCode())};
1325
1326                 if(posDaughter->GetMother() != negDaughter->GetMother()){
1327                         FillPhotonCombinatorialBackgroundHist(TruePhotonCandidate, pdgCode);
1328                         iPhotonMCInfo = 1;
1329                         return;
1330                 }
1331                 else if(posDaughter->GetMother() == -1){
1332                         FillPhotonCombinatorialBackgroundHist(TruePhotonCandidate, pdgCode);
1333                         iPhotonMCInfo = 1;
1334                         return;
1335                 }
1336
1337                 if(pdgCode[0]!=11 || pdgCode[1]!=11){
1338                         iPhotonMCInfo = 1;
1339                         return; //One Particle is not a electron
1340                 }
1341                 if(posDaughter->GetPdgCode()==negDaughter->GetPdgCode()){
1342                         iPhotonMCInfo = 1; 
1343                         return; // Same Charge
1344                 }
1345                 
1346                 AliAODMCParticle *Photon = (AliAODMCParticle*) AODMCTrackArray->At(posDaughter->GetMother());
1347                 AliVTrack * electronCandidate = ((AliConversionPhotonCuts*)fCutArray->At(fiCut))->GetTrack(fInputEvent,TruePhotonCandidate->GetTrackLabelNegative() );
1348                 AliVTrack * positronCandidate = ((AliConversionPhotonCuts*)fCutArray->At(fiCut))->GetTrack(fInputEvent,TruePhotonCandidate->GetTrackLabelPositive() );
1349                 Double_t deltaPhi = magField * TVector2::Phi_mpi_pi( electronCandidate->Phi()-positronCandidate->Phi());
1350
1351                 if(Photon->GetPdgCode() != 22){
1352                         hESDTrueDalitzPsiPairDeltaPhi[fiCut]->Fill(deltaPhi,TruePhotonCandidate->GetPsiPair());
1353                         iPhotonMCInfo = 1;
1354                         return; // Mother is no Photon
1355                 }
1356                 
1357                 if(((posDaughter->GetMCProcessCode())) != 5 || ((negDaughter->GetMCProcessCode())) != 5){
1358                         iPhotonMCInfo = 1;
1359                         return;// check if the daughters come from a conversion 
1360                 }   
1361                 // STILL A BUG IN ALIROOT >>8 HAS TPO BE REMOVED AFTER FIX
1362                 
1363                 
1364                 
1365                 // True Photon
1366                 if(fIsFromMBHeader){
1367                         hESDTrueConvGammaPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1368                         if (fDoPhotonQA > 0) hESDTrueConvGammaEta[fiCut]->Fill(TruePhotonCandidate->Eta());
1369                 }
1370                 hESDTrueGammaPsiPairDeltaPhi[fiCut]->Fill(deltaPhi,TruePhotonCandidate->GetPsiPair());  
1371                 if(Photon->IsPrimary()){ 
1372                         // Count just primary MC Gammas as true --> For Ratio esdtruegamma / mcconvgamma
1373                         if(fIsFromMBHeader){
1374                                 iPhotonMCInfo = 6;
1375                                 hESDTruePrimaryConvGammaPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1376                                 hESDTruePrimaryConvGammaESDPtMCPt[fiCut]->Fill(TruePhotonCandidate->Pt(),Photon->Pt()); // Allways Filled
1377                         }
1378                         // (Not Filled for i6, Extra Signal Gamma (parambox) are secondary)
1379                 }
1380                 else{
1381                         if(fIsFromMBHeader){
1382                                 hESDTrueSecondaryConvGammaPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1383                                 iPhotonMCInfo = 2;
1384                                 if(((AliAODMCParticle*)AODMCTrackArray->At(Photon->GetMother()))->GetMother() > -1 &&
1385                                         ((AliAODMCParticle*)AODMCTrackArray->At(((AliAODMCParticle*)AODMCTrackArray->At(Photon->GetMother()))->GetMother()))->GetPdgCode() == 3122){
1386                                         iPhotonMCInfo = 5;
1387                                         hESDTrueSecondaryConvGammaFromXFromLambdaPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1388                                 }
1389                                 if(((AliAODMCParticle*)AODMCTrackArray->At(Photon->GetMother()))->GetMother() > -1 &&
1390                                         ((AliAODMCParticle*)AODMCTrackArray->At(((AliAODMCParticle*)AODMCTrackArray->At(Photon->GetMother()))->GetMother()))->GetPdgCode() == 310){
1391                                         iPhotonMCInfo = 4;
1392                                         hESDTrueSecondaryConvGammaFromXFromK0sPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1393                                 }
1394                                 if(((AliAODMCParticle*)AODMCTrackArray->At(Photon->GetMother()))->GetMother() > -1 &&
1395                                         ((AliAODMCParticle*)AODMCTrackArray->At(((AliAODMCParticle*)AODMCTrackArray->At(Photon->GetMother()))->GetMother()))->GetPdgCode() == 221){
1396                                         iPhotonMCInfo = 3;
1397                                 }
1398                         }
1399                 }
1400         }
1401         return;
1402 }
1403 //________________________________________________________________________
1404 void AliAnalysisTaskGammaConvV1::ProcessTruePhotonCandidates(AliAODConversionPhoton *TruePhotonCandidate)
1405 {
1406         
1407         Double_t magField = fInputEvent->GetMagneticField();
1408         if( magField  < 0.0 ){
1409                 magField =  1.0;
1410         }
1411         else {
1412                 magField =  -1.0;
1413         }
1414
1415         // Process True Photons
1416         TParticle *posDaughter = TruePhotonCandidate->GetPositiveMCDaughter(fMCStack);
1417         TParticle *negDaughter = TruePhotonCandidate->GetNegativeMCDaughter(fMCStack);
1418
1419         iPhotonMCInfo = 0;
1420         
1421         if(posDaughter == NULL || negDaughter == NULL) return; // One particle does not exist
1422         Int_t pdgCode[2] = {abs(posDaughter->GetPdgCode()),abs(negDaughter->GetPdgCode())};
1423         iPhotonMCInfo = 1;
1424         if(posDaughter->GetMother(0) != negDaughter->GetMother(0)){
1425                 FillPhotonCombinatorialBackgroundHist(TruePhotonCandidate, pdgCode);
1426                 return;
1427         } else if(posDaughter->GetMother(0) == -1){
1428                 FillPhotonCombinatorialBackgroundHist(TruePhotonCandidate, pdgCode);
1429                 return;
1430         }
1431
1432         if(pdgCode[0]!=11 || pdgCode[1]!=11) return; //One Particle is not a electron
1433         
1434         if(posDaughter->GetPdgCode()==negDaughter->GetPdgCode()) return; // Same Charge
1435
1436         TParticle *Photon = TruePhotonCandidate->GetMCParticle(fMCStack);
1437         AliVTrack * electronCandidate = ((AliConversionPhotonCuts*)fCutArray->At(fiCut))->GetTrack(fInputEvent,TruePhotonCandidate->GetTrackLabelNegative() );
1438         AliVTrack * positronCandidate = ((AliConversionPhotonCuts*)fCutArray->At(fiCut))->GetTrack(fInputEvent,TruePhotonCandidate->GetTrackLabelPositive() );
1439         Double_t deltaPhi = magField * TVector2::Phi_mpi_pi( electronCandidate->Phi()-positronCandidate->Phi());
1440         
1441         if(Photon->GetPdgCode() != 22){
1442                 hESDTrueDalitzPsiPairDeltaPhi[fiCut]->Fill(deltaPhi,TruePhotonCandidate->GetPsiPair());
1443                 return; // Mother is no Photon
1444         }
1445         
1446         if(posDaughter->GetUniqueID() != 5 || negDaughter->GetUniqueID() !=5) return;// check if the daughters come from a conversion
1447
1448         
1449         
1450         // True Photon
1451         if(fIsFromMBHeader){
1452                 hESDTrueConvGammaPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1453                 if (fDoPhotonQA > 0) hESDTrueConvGammaEta[fiCut]->Fill(TruePhotonCandidate->Eta());
1454         }
1455         hESDTrueGammaPsiPairDeltaPhi[fiCut]->Fill(deltaPhi,TruePhotonCandidate->GetPsiPair());  
1456         if(posDaughter->GetMother(0) <= fMCStack->GetNprimary()){
1457                 // Count just primary MC Gammas as true --> For Ratio esdtruegamma / mcconvgamma
1458                 if(fIsFromMBHeader){
1459                         iPhotonMCInfo = 6;
1460                         hESDTruePrimaryConvGammaPt[fiCut]->Fill(TruePhotonCandidate->Pt());  
1461                         hESDTruePrimaryConvGammaESDPtMCPt[fiCut]->Fill(TruePhotonCandidate->Pt(),Photon->Pt()); // Allways Filled
1462
1463                 }
1464                 // (Not Filled for i6, Extra Signal Gamma (parambox) are secondary)
1465         } else {
1466                 if(fIsFromMBHeader){
1467                         iPhotonMCInfo = 2;
1468                         hESDTrueSecondaryConvGammaPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1469                         if(fMCStack->Particle(Photon->GetMother(0))->GetMother(0) > -1 &&
1470                                 fMCStack->Particle(fMCStack->Particle(Photon->GetMother(0))->GetMother(0))->GetPdgCode() == 3122){
1471                                 hESDTrueSecondaryConvGammaFromXFromLambdaPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1472                                 iPhotonMCInfo = 5;
1473                         }
1474                         if(fMCStack->Particle(Photon->GetMother(0))->GetMother(0) > -1 &&
1475                                 fMCStack->Particle(fMCStack->Particle(Photon->GetMother(0))->GetMother(0))->GetPdgCode() == 310){
1476                                 hESDTrueSecondaryConvGammaFromXFromK0sPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1477                                 iPhotonMCInfo = 4;
1478                         }
1479                         if(fMCStack->Particle(Photon->GetMother(0))->GetMother(0) > -1 &&
1480                                 fMCStack->Particle(fMCStack->Particle(Photon->GetMother(0))->GetMother(0))->GetPdgCode() == 221){
1481                                 iPhotonMCInfo = 3;
1482                         }
1483                 }
1484         }
1485 }
1486 //________________________________________________________________________
1487 void AliAnalysisTaskGammaConvV1::ProcessAODMCParticles()
1488 {
1489
1490         TClonesArray *AODMCTrackArray = dynamic_cast<TClonesArray*>(fInputEvent->FindListObject(AliAODMCParticle::StdBranchName()));
1491
1492         if (AODMCTrackArray){
1493                 // Loop over all primary MC particle
1494                 for(Int_t i = 0; i < AODMCTrackArray->GetEntriesFast(); i++) {
1495
1496                         AliAODMCParticle* particle = static_cast<AliAODMCParticle*>(AODMCTrackArray->At(i));
1497                         if (!particle) continue;
1498                         if (!particle->IsPrimary()) continue;
1499
1500                         Int_t isMCFromMBHeader = -1;
1501                         if(((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetSignalRejection() != 0){
1502                                 isMCFromMBHeader
1503                                         = ((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsParticleFromBGEvent(i, fMCStack, fInputEvent);
1504                                 if(isMCFromMBHeader == 0 && ((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetSignalRejection() != 3) continue;
1505                         }
1506                         
1507                         if(!((AliConversionPhotonCuts*)fCutArray->At(fiCut))->InPlaneOutOfPlaneCut(particle->Phi(),fEventPlaneAngle,kFALSE)) continue;
1508                         if(((AliConversionPhotonCuts*)fCutArray->At(fiCut))->PhotonIsSelectedAODMC(particle,AODMCTrackArray,kFALSE)){
1509                                 hMCAllGammaPt[fiCut]->Fill(particle->Pt()); // All MC Gamma
1510                                 if(particle->GetMother() >-1){ // Meson Decay Gamma
1511                                         switch((static_cast<AliAODMCParticle*>(AODMCTrackArray->At(particle->GetMother())))->GetPdgCode()){
1512                                         case 111: // Pi0
1513                                         hMCDecayGammaPi0Pt[fiCut]->Fill(particle->Pt());
1514                                         break;
1515                                         case 113: // Rho0
1516                                         hMCDecayGammaRhoPt[fiCut]->Fill(particle->Pt());
1517                                         break;
1518                                         case 221: // Eta
1519                                         hMCDecayGammaEtaPt[fiCut]->Fill(particle->Pt());
1520                                         break;
1521                                         case 223: // Omega
1522                                         hMCDecayGammaOmegaPt[fiCut]->Fill(particle->Pt());
1523                                         break;
1524                                         case 331: // Eta'
1525                                         hMCDecayGammaEtapPt[fiCut]->Fill(particle->Pt());
1526                                         break;
1527                                         case 333: // Phi
1528                                         hMCDecayGammaPhiPt[fiCut]->Fill(particle->Pt());
1529                                         break;
1530                                         case 3212: // Sigma
1531                                         hMCDecayGammaSigmaPt[fiCut]->Fill(particle->Pt());
1532                                         break;
1533                                         }
1534                                 }
1535                         }
1536                         if(((AliConversionPhotonCuts*)fCutArray->At(fiCut))->PhotonIsSelectedAODMC(particle,AODMCTrackArray,kTRUE)){
1537                                 Double_t rConv = 0;
1538                                 for(Int_t daughterIndex=particle->GetDaughter(0);daughterIndex<=particle->GetDaughter(1);daughterIndex++){
1539                                         AliAODMCParticle *tmpDaughter = static_cast<AliAODMCParticle*>(AODMCTrackArray->At(daughterIndex));
1540                                         if(!tmpDaughter) continue;
1541                                         if(abs(tmpDaughter->GetPdgCode()) == 11){
1542                                                 rConv = sqrt( (tmpDaughter->Xv()*tmpDaughter->Xv()) + (tmpDaughter->Yv()*tmpDaughter->Yv()) );
1543                                         }
1544                                 }
1545                                 hMCConvGammaPt[fiCut]->Fill(particle->Pt());
1546                                 if (fDoPhotonQA > 0){
1547                                         hMCConvGammaR[fiCut]->Fill(rConv);
1548                                         hMCConvGammaEta[fiCut]->Fill(particle->Eta());
1549                                 }
1550                         }
1551                         // Converted MC Gamma
1552                         if(fDoMesonAnalysis){
1553                                 if(particle->GetPdgCode() == 310 && fDoMesonQA > 0){
1554                                         Double_t mesonY = 10.;
1555                                         if(particle->E() - particle->Pz() == 0 || particle->E() + particle->Pz() == 0){
1556                                                 mesonY=10.-((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetEtaShift();
1557                                         } else{
1558                                                 mesonY = 0.5*(TMath::Log((particle->E()+particle->Pz()) / (particle->E()-particle->Pz())))
1559                                                 -((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetEtaShift();
1560                                         }
1561                                         Float_t weightedK0s= 1;
1562                                         if(((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsParticleFromBGEvent(i, fMCStack, fInputEvent)){
1563                                                 if (particle->Pt()>0.005){
1564                                                         weightedK0s= ((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetWeightForMeson(fV0Reader->GetPeriodName(),i, 0x0, fInputEvent);
1565                                                         //cout << "MC input \t"<<i << "\t" <<  particle->Pt()<<"\t"<<weighted << endl;
1566                                                 }
1567                                         }
1568                                         hMCK0sPt[fiCut]->Fill(particle->Pt(),weightedK0s);
1569                                         hMCK0sWOWeightPt[fiCut]->Fill(particle->Pt());
1570                                         hMCK0sPtY[fiCut]->Fill(particle->Pt(),mesonY,weightedK0s);
1571                                 }         
1572                                 if(((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->MesonIsSelectedAODMC(particle,AODMCTrackArray,((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetEtaShift())){
1573                                         AliAODMCParticle* daughter0 = static_cast<AliAODMCParticle*>(AODMCTrackArray->At(particle->GetDaughter(0)));
1574                                         AliAODMCParticle* daughter1 = static_cast<AliAODMCParticle*>(AODMCTrackArray->At(particle->GetDaughter(1)));
1575                                         Float_t weighted= 1;
1576                                         if(((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsParticleFromBGEvent(i, fMCStack, fInputEvent)){
1577                                                 if (particle->Pt()>0.005){
1578                                                         weighted= ((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetWeightForMeson(fV0Reader->GetPeriodName(),i, 0x0, fInputEvent);
1579                         //                   if(particle->GetPdgCode() == 221){
1580                         //                      cout << "MC input \t"<<i << "\t" <<  particle->Pt()<<"\t"<<weighted << endl;
1581                         //                   }
1582                                                 }
1583                                         }
1584                                         Double_t mesonY = 10.;
1585                                         if(particle->E() - particle->Pz() == 0 || particle->E() + particle->Pz() == 0){
1586                                                 mesonY=10.-((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetEtaShift();
1587                                         } else{
1588                                                 mesonY = 0.5*(TMath::Log((particle->E()+particle->Pz()) / (particle->E()-particle->Pz())))
1589                                                 -((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetEtaShift();
1590                                         }
1591
1592                                         Double_t alpha = -1;
1593                                         if (particle->GetPdgCode() == 111 || particle->GetPdgCode() == 221){
1594                                                 alpha = TMath::Abs((daughter0->E() - daughter1->E()))/(daughter0->E() + daughter1->E());
1595                                         }
1596
1597                                         if(particle->GetPdgCode() == 111){
1598                                                 hMCPi0Pt[fiCut]->Fill(particle->Pt(),weighted); // All MC Pi0
1599                                                 hMCPi0WOWeightPt[fiCut]->Fill(particle->Pt());
1600                                                 if (fDoMesonQA > 0){
1601                                                         hMCPi0PtY[fiCut]->Fill(particle->Pt(),mesonY,weighted); // All MC Pi0
1602                                                         hMCPi0PtAlpha[fiCut]->Fill(particle->Pt(),alpha); // All MC Pi0
1603                                                 }       
1604                                         } else if(particle->GetPdgCode() == 221){
1605                                                 hMCEtaPt[fiCut]->Fill(particle->Pt(),weighted); // All MC Eta
1606                                                 hMCEtaWOWeightPt[fiCut]->Fill(particle->Pt());
1607                                                 if (fDoMesonQA > 0){
1608                                                         hMCEtaPtY[fiCut]->Fill(particle->Pt(),mesonY,weighted); // All MC Pi0
1609                                                         hMCEtaPtAlpha[fiCut]->Fill(particle->Pt(),alpha); // All MC Pi0
1610                                                 }       
1611                                         }
1612                                         
1613                                         // Check the acceptance for both gammas
1614                                         if(((AliConversionPhotonCuts*)fCutArray->At(fiCut))->PhotonIsSelectedAODMC(daughter0,AODMCTrackArray,kFALSE) &&
1615                                         ((AliConversionPhotonCuts*)fCutArray->At(fiCut))->PhotonIsSelectedAODMC(daughter1,AODMCTrackArray,kFALSE)  &&
1616                                         ((AliConversionPhotonCuts*)fCutArray->At(fiCut))->InPlaneOutOfPlaneCut(daughter0->Phi(),fEventPlaneAngle,kFALSE) &&
1617                                         ((AliConversionPhotonCuts*)fCutArray->At(fiCut))->InPlaneOutOfPlaneCut(daughter1->Phi(),fEventPlaneAngle,kFALSE)){
1618
1619                                                 if(particle->GetPdgCode() == 111){
1620                                                         hMCPi0InAccPt[fiCut]->Fill(particle->Pt(),weighted); // MC Pi0 with gamma in acc
1621                                                 } else if(particle->GetPdgCode() == 221){
1622                                                         hMCEtaInAccPt[fiCut]->Fill(particle->Pt(),weighted); // MC Eta with gamma in acc
1623                                                 }
1624                                         }
1625                                 }
1626                         }
1627                 }
1628         }
1629         return;
1630 }
1631 //________________________________________________________________________
1632 void AliAnalysisTaskGammaConvV1::ProcessMCParticles()
1633 {
1634         // Loop over all primary MC particle
1635         for(Int_t i = 0; i < fMCStack->GetNprimary(); i++) {
1636                 TParticle* particle = (TParticle *)fMCStack->Particle(i);
1637                 if (!particle) continue;
1638
1639                 Int_t isMCFromMBHeader = -1;
1640                 if(((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetSignalRejection() != 0){
1641                         isMCFromMBHeader
1642                                 = ((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsParticleFromBGEvent(i, fMCStack, fInputEvent);
1643                         if(isMCFromMBHeader == 0 && ((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetSignalRejection() != 3) continue;
1644                 }
1645
1646                 if(!((AliConversionPhotonCuts*)fCutArray->At(fiCut))->InPlaneOutOfPlaneCut(particle->Phi(),fEventPlaneAngle,kFALSE)) continue;
1647                 if(((AliConversionPhotonCuts*)fCutArray->At(fiCut))->PhotonIsSelectedMC(particle,fMCStack,kFALSE)){
1648                         hMCAllGammaPt[fiCut]->Fill(particle->Pt()); // All MC Gamma
1649                         if(particle->GetMother(0) >-1){ // Meson Decay Gamma
1650                                 switch(fMCStack->Particle(particle->GetMother(0))->GetPdgCode()){
1651                                 case 111: // Pi0
1652                                 hMCDecayGammaPi0Pt[fiCut]->Fill(particle->Pt());
1653                                 break;
1654                                 case 113: // Rho0
1655                                 hMCDecayGammaRhoPt[fiCut]->Fill(particle->Pt());
1656                                 break;
1657                                 case 221: // Eta
1658                                 hMCDecayGammaEtaPt[fiCut]->Fill(particle->Pt());
1659                                 break;
1660                                 case 223: // Omega
1661                                 hMCDecayGammaOmegaPt[fiCut]->Fill(particle->Pt());
1662                                 break;
1663                                 case 331: // Eta'
1664                                 hMCDecayGammaEtapPt[fiCut]->Fill(particle->Pt());
1665                                 break;
1666                                 case 333: // Phi
1667                                 hMCDecayGammaPhiPt[fiCut]->Fill(particle->Pt());
1668                                 break;
1669                                 case 3212: // Sigma
1670                                 hMCDecayGammaSigmaPt[fiCut]->Fill(particle->Pt());
1671                                 break;
1672                                 }
1673                         }
1674                 }
1675                 if(((AliConversionPhotonCuts*)fCutArray->At(fiCut))->PhotonIsSelectedMC(particle,fMCStack,kTRUE)){
1676                         hMCConvGammaPt[fiCut]->Fill(particle->Pt());
1677                         if (fDoPhotonQA > 0){
1678                                 hMCConvGammaR[fiCut]->Fill(((TParticle*)fMCStack->Particle(particle->GetFirstDaughter()))->R());
1679                                 hMCConvGammaEta[fiCut]->Fill(particle->Eta());
1680                         }
1681                 } // Converted MC Gamma
1682                 if(fDoMesonAnalysis){
1683                         if(particle->GetPdgCode() == 310 && fDoMesonQA > 0){
1684                                 Double_t mesonY = 10.;
1685                                 if(particle->Energy() - particle->Pz() == 0 || particle->Energy() + particle->Pz() == 0){
1686                                         mesonY=10.-((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetEtaShift();
1687                                 } else{
1688                                         mesonY = 0.5*(TMath::Log((particle->Energy()+particle->Pz()) / (particle->Energy()-particle->Pz())))
1689                                         -((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetEtaShift();
1690                                 }
1691                                 Float_t weightedK0s= 1;
1692                                 if(((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsParticleFromBGEvent(i, fMCStack, fInputEvent)){
1693                                         if (particle->Pt()>0.005){
1694                                                 weightedK0s= ((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetWeightForMeson(fV0Reader->GetPeriodName(),i, fMCStack, fInputEvent);
1695                                                 //cout << "MC input \t"<<i << "\t" <<  particle->Pt()<<"\t"<<weighted << endl;
1696                                         }
1697                                 }
1698                                 if (fMCStack->IsPhysicalPrimary(i)){
1699                                         hMCK0sPt[fiCut]->Fill(particle->Pt(),weightedK0s);
1700                                         hMCK0sWOWeightPt[fiCut]->Fill(particle->Pt());
1701                                         hMCK0sPtY[fiCut]->Fill(particle->Pt(),mesonY,weightedK0s);
1702                                 }       
1703                         }
1704                         if(((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))
1705                                 ->MesonIsSelectedMC(particle,fMCStack,((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetEtaShift())){
1706                                 TParticle* daughter0 = (TParticle*)fMCStack->Particle(particle->GetFirstDaughter());
1707                                 TParticle* daughter1 = (TParticle*)fMCStack->Particle(particle->GetLastDaughter());
1708
1709                                 Float_t weighted= 1;
1710                                 if(((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsParticleFromBGEvent(i, fMCStack, fInputEvent)){
1711                                         if (particle->Pt()>0.005){
1712                                                 weighted= ((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetWeightForMeson(fV0Reader->GetPeriodName(),i, fMCStack, fInputEvent);
1713                 //                   if(particle->GetPdgCode() == 221){
1714                 //                      cout << "MC input \t"<<i << "\t" <<  particle->Pt()<<"\t"<<weighted << endl;
1715                 //                   }
1716                                         }
1717                                 }
1718                                 Double_t mesonY = 10.;
1719                                 if(particle->Energy() - particle->Pz() == 0 || particle->Energy() + particle->Pz() == 0){
1720                                         mesonY=10.-((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetEtaShift();
1721                                 } else {
1722                                         mesonY = 0.5*(TMath::Log((particle->Energy()+particle->Pz()) / (particle->Energy()-particle->Pz())))
1723                                         -((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetEtaShift();
1724                                 }
1725                                 Double_t alpha = -1;
1726                                 if (particle->GetPdgCode() == 111 || particle->GetPdgCode() == 221){
1727                                         alpha = TMath::Abs((daughter0->Energy() - daughter1->Energy()))/(daughter0->Energy() + daughter1->Energy());
1728                                 }
1729
1730                                 if(particle->GetPdgCode() == 111){
1731                                         hMCPi0Pt[fiCut]->Fill(particle->Pt(),weighted); // All MC Pi0
1732                                         hMCPi0WOWeightPt[fiCut]->Fill(particle->Pt());
1733                                         if (fDoMesonQA > 0){
1734                                                 hMCPi0PtY[fiCut]->Fill(particle->Pt(),mesonY,weighted); // All MC Pi0
1735                                                 hMCPi0PtAlpha[fiCut]->Fill(particle->Pt(),alpha); // All MC Pi0
1736                                         }       
1737                                 } else if(particle->GetPdgCode() == 221){
1738                                         hMCEtaPt[fiCut]->Fill(particle->Pt(),weighted); // All MC Eta
1739                                         hMCEtaWOWeightPt[fiCut]->Fill(particle->Pt());
1740                                         if (fDoMesonQA > 0){
1741                                                 hMCEtaPtY[fiCut]->Fill(particle->Pt(),mesonY,weighted); // All MC Pi0
1742                                                 hMCEtaPtAlpha[fiCut]->Fill(particle->Pt(),alpha); // All MC Pi0
1743                                         }       
1744                                 } 
1745
1746                                 // Check the acceptance for both gammas
1747                                 if(((AliConversionPhotonCuts*)fCutArray->At(fiCut))->PhotonIsSelectedMC(daughter0,fMCStack,kFALSE) &&
1748                                 ((AliConversionPhotonCuts*)fCutArray->At(fiCut))->PhotonIsSelectedMC(daughter1,fMCStack,kFALSE)  &&
1749                                 ((AliConversionPhotonCuts*)fCutArray->At(fiCut))->InPlaneOutOfPlaneCut(daughter0->Phi(),fEventPlaneAngle,kFALSE) &&
1750                                 ((AliConversionPhotonCuts*)fCutArray->At(fiCut))->InPlaneOutOfPlaneCut(daughter1->Phi(),fEventPlaneAngle,kFALSE)){
1751
1752                                         if(particle->GetPdgCode() == 111){
1753                                                 hMCPi0InAccPt[fiCut]->Fill(particle->Pt(),weighted); // MC Pi0 with gamma in acc
1754                                         } else if(particle->GetPdgCode() == 221){
1755                                                 hMCEtaInAccPt[fiCut]->Fill(particle->Pt(),weighted); // MC Eta with gamma in acc
1756                                         }
1757                                 }
1758                         }
1759                 }
1760         }
1761
1762         if (fDoMesonQA){
1763                 for(Int_t i = fMCStack->GetNprimary(); i < fMCStack->GetNtrack(); i++) {
1764                         TParticle* particle = (TParticle *)fMCStack->Particle(i);                       
1765                         if (!particle) continue;
1766
1767                         Int_t isMCFromMBHeader = -1;
1768                         if(((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetSignalRejection() != 0){
1769                                 isMCFromMBHeader
1770                                         = ((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsParticleFromBGEvent(i, fMCStack, fInputEvent);
1771                                 if(isMCFromMBHeader == 0 && ((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetSignalRejection() != 3) continue;
1772                         }
1773
1774                         if(fDoMesonAnalysis){
1775                                 if(((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->MesonIsSelectedMC(particle,fMCStack,((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetEtaShift())){
1776                                         Float_t weighted= 1;
1777                                         if(((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsParticleFromBGEvent(i, fMCStack, fInputEvent)){
1778                                                 if (particle->Pt()>0.005){
1779                                                         weighted= ((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetWeightForMeson(fV0Reader->GetPeriodName(),i, fMCStack, fInputEvent);
1780                         //                   if(particle->GetPdgCode() == 221){
1781                         //                      cout << "MC input \t"<<i << "\t" <<  particle->Pt()<<"\t"<<weighted << endl;
1782                         //                   }
1783                                                 }
1784                                         }
1785                                         
1786                                         if(particle->GetPdgCode() == 111){      
1787                                                 Int_t pdgCode = ((TParticle*)fMCStack->Particle( particle->GetFirstMother() ))->GetPdgCode();
1788                                                 Int_t source = GetSourceClassification(111,pdgCode);
1789                                                 hMCSecPi0PtvsSource[fiCut]->Fill(particle->Pt(),source,weighted); // All MC Pi0
1790                                                 hMCSecPi0Source[fiCut]->Fill(pdgCode);                                          
1791                                         } else if(particle->GetPdgCode() == 221){
1792                                                 Int_t pdgCode = ((TParticle*)fMCStack->Particle( particle->GetFirstMother() ))->GetPdgCode();
1793                                                 hMCSecEtaPt[fiCut]->Fill(particle->Pt(),weighted); // All MC Pi0
1794                                                 hMCSecEtaSource[fiCut]->Fill(pdgCode);
1795                                         } 
1796                                 }
1797                         }
1798                 }
1799         }
1800 }
1801
1802 //________________________________________________________________________
1803 void AliAnalysisTaskGammaConvV1::CalculatePi0Candidates(){
1804
1805         // Conversion Gammas
1806         if(fGammaCandidates->GetEntries()>1){
1807                 for(Int_t firstGammaIndex=0;firstGammaIndex<fGammaCandidates->GetEntries()-1;firstGammaIndex++){
1808                         AliAODConversionPhoton *gamma0=dynamic_cast<AliAODConversionPhoton*>(fGammaCandidates->At(firstGammaIndex));
1809                         if (gamma0==NULL) continue;
1810                         for(Int_t secondGammaIndex=firstGammaIndex+1;secondGammaIndex<fGammaCandidates->GetEntries();secondGammaIndex++){
1811                                 AliAODConversionPhoton *gamma1=dynamic_cast<AliAODConversionPhoton*>(fGammaCandidates->At(secondGammaIndex));
1812                                 //Check for same Electron ID
1813                                 if (gamma1==NULL) continue;
1814                                 if(gamma0->GetTrackLabelPositive() == gamma1->GetTrackLabelPositive() ||
1815                                 gamma0->GetTrackLabelNegative() == gamma1->GetTrackLabelNegative() ||
1816                                 gamma0->GetTrackLabelNegative() == gamma1->GetTrackLabelPositive() ||
1817                                 gamma0->GetTrackLabelPositive() == gamma1->GetTrackLabelNegative() ) continue;
1818
1819                                 AliAODConversionMother *pi0cand = new AliAODConversionMother(gamma0,gamma1);
1820                                 pi0cand->SetLabels(firstGammaIndex,secondGammaIndex);
1821
1822                                 pi0cand->CalculateDistanceOfClossetApproachToPrimVtx(fInputEvent->GetPrimaryVertex());
1823                                 if((((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->MesonIsSelected(pi0cand,kTRUE,((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetEtaShift()))){
1824                                         hESDMotherInvMassPt[fiCut]->Fill(pi0cand->M(),pi0cand->Pt());
1825                                         if(pi0cand->GetAlpha()<0.1) hESDMotherInvMassEalpha[fiCut]->Fill(pi0cand->M(),pi0cand->E());
1826                                 
1827                                         if (fDoMesonQA > 0){
1828                                                 if ( pi0cand->M() > 0.05 && pi0cand->M() < 0.17){
1829                                                         hESDMotherPi0PtY[fiCut]->Fill(pi0cand->Pt(),pi0cand->Rapidity()-((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetEtaShift());  
1830                                                         hESDMotherPi0PtAlpha[fiCut]->Fill(pi0cand->Pt(),pi0cand->GetAlpha());  
1831                                                         hESDMotherPi0PtOpenAngle[fiCut]->Fill(pi0cand->Pt(),pi0cand->GetOpeningAngle()); 
1832                                                         
1833                                                 } 
1834                                                 if ( pi0cand->M() > 0.45 && pi0cand->M() < 0.65){
1835                                                         hESDMotherEtaPtY[fiCut]->Fill(pi0cand->Pt(),pi0cand->Rapidity()-((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetEtaShift());  
1836                                                         hESDMotherEtaPtAlpha[fiCut]->Fill(pi0cand->Pt(),pi0cand->GetAlpha());  
1837                                                         hESDMotherEtaPtOpenAngle[fiCut]->Fill(pi0cand->Pt(),pi0cand->GetOpeningAngle());       
1838                                                 }
1839                                         }   
1840                                         if(((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->DoBGCalculation()){
1841                                                 Int_t zbin = 0;
1842                                                 Int_t mbin = 0;
1843
1844                                                 if(((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->BackgroundHandlerType() == 0){
1845                                                         zbin = fBGHandler[fiCut]->GetZBinIndex(fInputEvent->GetPrimaryVertex()->GetZ());
1846                                                         if(((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->UseTrackMultiplicity()){
1847                                                                 mbin = fBGHandler[fiCut]->GetMultiplicityBinIndex(fV0Reader->GetNumberOfPrimaryTracks());
1848                                                         } else {
1849                                                                 mbin = fBGHandler[fiCut]->GetMultiplicityBinIndex(fGammaCandidates->GetEntries());
1850                                                         }
1851                                                 }
1852                                                 else{
1853                                                         zbin = fBGHandlerRP[fiCut]->GetZBinIndex(fInputEvent->GetPrimaryVertex()->GetZ());
1854                                                         if(((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->UseTrackMultiplicity()){
1855                                                                 mbin = fBGHandlerRP[fiCut]->GetMultiplicityBinIndex(fV0Reader->GetNumberOfPrimaryTracks());
1856                                                         } else {
1857                                                                 mbin = fBGHandlerRP[fiCut]->GetMultiplicityBinIndex(fGammaCandidates->GetEntries());
1858                                                         }
1859                                                 }
1860                                                 Double_t sparesFill[4] = {pi0cand->M(),pi0cand->Pt(),(Double_t)zbin,(Double_t)mbin};
1861                                                 sESDMotherInvMassPtZM[fiCut]->Fill(sparesFill,1);
1862                                         }
1863                                         
1864
1865                                         if(fIsMC){
1866                                                 if(fInputEvent->IsA()==AliESDEvent::Class())
1867                                                         ProcessTrueMesonCandidates(pi0cand,gamma0,gamma1);
1868                                                 if(fInputEvent->IsA()==AliAODEvent::Class())
1869                                                         ProcessTrueMesonCandidatesAOD(pi0cand,gamma0,gamma1);
1870                                         }
1871                                         if (fDoMesonQA == 2){
1872                                                 fInvMass = pi0cand->M();
1873                                                 fPt  = pi0cand->Pt();
1874                                                 if (abs(gamma0->GetDCAzToPrimVtx()) < abs(gamma1->GetDCAzToPrimVtx())){
1875                                                         fDCAzGammaMin = gamma0->GetDCAzToPrimVtx();
1876                                                         fDCAzGammaMax = gamma1->GetDCAzToPrimVtx();
1877                                                 } else {
1878                                                         fDCAzGammaMin = gamma1->GetDCAzToPrimVtx();
1879                                                         fDCAzGammaMax = gamma0->GetDCAzToPrimVtx();
1880                                                 }
1881                                                 iFlag = pi0cand->GetMesonQuality();
1882                 //                   cout << "gamma 0: " << gamma0->GetV0Index()<< "\t" << gamma0->GetPx() << "\t" << gamma0->GetPy() << "\t" <<  gamma0->GetPz() << "\t" << endl; 
1883                 //                   cout << "gamma 1: " << gamma1->GetV0Index()<< "\t"<< gamma1->GetPx() << "\t" << gamma1->GetPy() << "\t" <<  gamma1->GetPz() << "\t" << endl; 
1884                 //                    cout << "pi0: "<<fInvMass << "\t" << fPt <<"\t" << fDCAzGammaMin << "\t" << fDCAzGammaMax << "\t" << (Int_t)iFlag << "\t" << (Int_t)iMesonMCInfo <<endl;
1885                                                 if (fIsHeavyIon == 1 && fPt > 0.399 && fPt < 20. ) {
1886                                                         if (fInvMass > 0.08 && fInvMass < 0.2) tESDMesonsInvMassPtDcazMinDcazMaxFlag[fiCut]->Fill();
1887                                                         if ((fInvMass > 0.45 && fInvMass < 0.6) &&  (fPt > 0.999 && fPt < 20.) )tESDMesonsInvMassPtDcazMinDcazMaxFlag[fiCut]->Fill();
1888                                                 } else if (fPt > 0.299 && fPt < 20. )  {
1889                                                         if ( (fInvMass > 0.08 && fInvMass < 0.2) || (fInvMass > 0.45 && fInvMass < 0.6)) tESDMesonsInvMassPtDcazMinDcazMaxFlag[fiCut]->Fill();
1890                                                 }   
1891                                         }
1892                                 }
1893                                 delete pi0cand;
1894                                 pi0cand=0x0;
1895                         }
1896                 }
1897         }
1898 }
1899
1900 //______________________________________________________________________
1901 void AliAnalysisTaskGammaConvV1::ProcessTrueMesonCandidates(AliAODConversionMother *Pi0Candidate, AliAODConversionPhoton *TrueGammaCandidate0, AliAODConversionPhoton *TrueGammaCandidate1)
1902 {
1903         // Process True Mesons
1904         AliStack *MCStack = fMCEvent->Stack();
1905         iMesonMCInfo = 0;
1906         if(TrueGammaCandidate0->GetV0Index()<fInputEvent->GetNumberOfV0s()){
1907                 Bool_t isTruePi0 = kFALSE;
1908                 Bool_t isTrueEta = kFALSE;
1909                 Bool_t isTruePi0Dalitz = kFALSE;
1910                 Bool_t isTrueEtaDalitz = kFALSE;
1911                 Bool_t gamma0DalitzCand = kFALSE;
1912                 Bool_t gamma1DalitzCand = kFALSE;
1913                 Int_t gamma0MCLabel = TrueGammaCandidate0->GetMCParticleLabel(MCStack);
1914                 Int_t gamma0MotherLabel = -1;
1915                 if(gamma0MCLabel != -1){ // Gamma is Combinatorial; MC Particles don't belong to the same Mother
1916                         // Daughters Gamma 0
1917                         TParticle * negativeMC = (TParticle*)TrueGammaCandidate0->GetNegativeMCDaughter(MCStack);
1918                         TParticle * positiveMC = (TParticle*)TrueGammaCandidate0->GetPositiveMCDaughter(MCStack);
1919                         TParticle * gammaMC0 = (TParticle*)MCStack->Particle(gamma0MCLabel);
1920                         if(abs(negativeMC->GetPdgCode())==11 && abs(positiveMC->GetPdgCode())==11){  // Electrons ...
1921                                 if(negativeMC->GetUniqueID() == 5 && positiveMC->GetUniqueID() ==5){ // ... From Conversion ...
1922                                         if(gammaMC0->GetPdgCode() == 22){ // ... with Gamma Mother
1923                                                 gamma0MotherLabel=gammaMC0->GetFirstMother();
1924                                         }
1925                                 }
1926                                 if(gammaMC0->GetPdgCode() ==111){ // Dalitz candidate
1927                                         gamma0DalitzCand = kTRUE;
1928                                         gamma0MotherLabel=-111;
1929                                 }
1930                                 if(gammaMC0->GetPdgCode() ==221){ // Dalitz candidate
1931                                         gamma0DalitzCand = kTRUE;
1932                                         gamma0MotherLabel=-221;
1933                                 }
1934                         }
1935                 }
1936                 if(TrueGammaCandidate1->GetV0Index()<fInputEvent->GetNumberOfV0s()){
1937                         Int_t gamma1MCLabel = TrueGammaCandidate1->GetMCParticleLabel(MCStack);
1938                         Int_t gamma1MotherLabel = -1;
1939                         if(gamma1MCLabel != -1){ // Gamma is Combinatorial; MC Particles don't belong to the same Mother
1940                                 // Daughters Gamma 1
1941                                 TParticle * negativeMC = (TParticle*)TrueGammaCandidate1->GetNegativeMCDaughter(MCStack);
1942                                 TParticle * positiveMC = (TParticle*)TrueGammaCandidate1->GetPositiveMCDaughter(MCStack);
1943                                 TParticle * gammaMC1 = (TParticle*)MCStack->Particle(gamma1MCLabel);
1944                                 if(abs(negativeMC->GetPdgCode())==11 && abs(positiveMC->GetPdgCode())==11){  // Electrons ...
1945                                         if(negativeMC->GetUniqueID() == 5 && positiveMC->GetUniqueID() ==5){ // ... From Conversion ...
1946                                                 if(gammaMC1->GetPdgCode() == 22){ // ... with Gamma Mother
1947                                                         gamma1MotherLabel=gammaMC1->GetFirstMother();
1948                                                 }
1949                                         }
1950                                         if(gammaMC1->GetPdgCode() ==111 ){ // Dalitz candidate
1951                                                 gamma1DalitzCand = kTRUE;
1952                                                 gamma1MotherLabel=-111;
1953                                         }
1954                                         if(gammaMC1->GetPdgCode() ==221){ // Dalitz candidate
1955                                                 gamma1DalitzCand = kTRUE;
1956                                                 gamma1MotherLabel=-221;
1957                                         }
1958                                 }
1959                         }
1960                         if(gamma0MotherLabel>=0 && gamma0MotherLabel==gamma1MotherLabel){
1961                                 if(((TParticle*)MCStack->Particle(gamma1MotherLabel))->GetPdgCode() == 111){
1962                                         isTruePi0=kTRUE;
1963                                 }
1964                                 if(((TParticle*)MCStack->Particle(gamma1MotherLabel))->GetPdgCode() == 221){
1965                                         isTrueEta=kTRUE;
1966                                 }
1967                         }
1968                         
1969                         //Identify Dalitz candidate
1970                         if (gamma1DalitzCand || gamma0DalitzCand){
1971                                 if (gamma0DalitzCand && gamma0MCLabel >=0 && gamma0MCLabel==gamma1MotherLabel){
1972                                         if (gamma0MotherLabel == -111) isTruePi0Dalitz = kTRUE;
1973                                         if (gamma0MotherLabel == -221) isTrueEtaDalitz = kTRUE;   
1974                                 }   
1975                                 if (gamma1DalitzCand && gamma1MCLabel >=0 && gamma1MCLabel==gamma0MotherLabel){
1976                                         if (gamma1MotherLabel == -111) isTruePi0Dalitz = kTRUE;
1977                                         if (gamma1MotherLabel == -221) isTrueEtaDalitz = kTRUE;   
1978                                 }
1979                         }
1980                         
1981                         
1982                         if(isTruePi0 || isTrueEta){// True Pion or Eta
1983                                 hESDTrueMotherInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt()); 
1984                                 if (fDoMesonQA > 0){
1985                                         if (isTruePi0){
1986                                                 if ( Pi0Candidate->M() > 0.05 && Pi0Candidate->M() < 0.17){
1987                                                         hESDTruePi0PtY[fiCut]->Fill(Pi0Candidate->Pt(),Pi0Candidate->Rapidity()-((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetEtaShift()); 
1988                                                         hESDTruePi0PtAlpha[fiCut]->Fill(Pi0Candidate->Pt(),Pi0Candidate->GetAlpha()); 
1989                                                         hESDTruePi0PtOpenAngle[fiCut]->Fill(Pi0Candidate->Pt(),Pi0Candidate->GetOpeningAngle()); 
1990                                                 }
1991                                         } else if (isTrueEta){   
1992                                                 if ( Pi0Candidate->M() > 0.45 && Pi0Candidate->M() < 0.65){
1993                                                         hESDTrueEtaPtY[fiCut]->Fill(Pi0Candidate->Pt(),Pi0Candidate->Rapidity()-((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetEtaShift()); 
1994                                                         hESDTrueEtaPtAlpha[fiCut]->Fill(Pi0Candidate->Pt(),Pi0Candidate->GetAlpha()); 
1995                                                         hESDTrueEtaPtOpenAngle[fiCut]->Fill(Pi0Candidate->Pt(),Pi0Candidate->GetOpeningAngle()); 
1996                                                 }
1997                                         }
1998                                 }   
1999                                 if(gamma0MotherLabel >= MCStack->GetNprimary()){ // Secondary Meson
2000                                         Int_t secMotherLabel = ((TParticle*)MCStack->Particle(gamma1MotherLabel))->GetMother(0);
2001                                         Float_t weightedSec= 1;
2002                                         if(((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsParticleFromBGEvent(secMotherLabel, fMCStack, fInputEvent) && MCStack->Particle(secMotherLabel)->GetPdgCode()==310){
2003                                                 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
2004                                                 //cout << "MC input \t"<<i << "\t" <<  particle->Pt()<<"\t"<<weighted << endl;
2005                                         }
2006                                         hESDTrueSecondaryMotherInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt(),weightedSec);
2007                                         iMesonMCInfo = 2;
2008                                         if (secMotherLabel >-1){
2009                                                 if(MCStack->Particle(secMotherLabel)->GetPdgCode()==310){
2010                                                         iMesonMCInfo = 4;
2011                                                         hESDTrueSecondaryMotherFromK0sInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt(),weightedSec);
2012                                                         if (fDoMesonQA > 0)hESDTrueK0sWithPi0DaughterMCPt[fiCut]
2013                                                                                         ->Fill(MCStack->Particle(secMotherLabel)->Pt());
2014                                                 }
2015                                                 if(MCStack->Particle(secMotherLabel)->GetPdgCode()==221){
2016                                                         iMesonMCInfo = 3;
2017                                                         hESDTrueSecondaryMotherFromEtaInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt(),weightedSec);
2018                                                         if (fDoMesonQA > 0)hESDTrueEtaWithPi0DaughterMCPt[fiCut]
2019                                                                                         ->Fill(MCStack->Particle(secMotherLabel)->Pt());
2020                                                 }
2021                                                 if(MCStack->Particle(secMotherLabel)->GetPdgCode()==3122){
2022                                                         iMesonMCInfo = 7;
2023                                                         hESDTrueSecondaryMotherFromLambdaInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt(),weightedSec);
2024                                                         if (fDoMesonQA > 0)hESDTrueLambdaWithPi0DaughterMCPt[fiCut]
2025                                                                                         ->Fill(MCStack->Particle(secMotherLabel)->Pt());
2026                                                 }
2027                                         }
2028                                 } else { // Only primary pi0 for efficiency calculation
2029                                         iMesonMCInfo = 6;
2030                                         Float_t weighted= 1;
2031                                         if(((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsParticleFromBGEvent(gamma1MotherLabel, fMCStack, fInputEvent)){
2032                                                 if (((TParticle*)MCStack->Particle(gamma1MotherLabel))->Pt()>0.005){
2033                                                         weighted= ((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetWeightForMeson(fV0Reader->GetPeriodName(),gamma1MotherLabel, fMCStack, fInputEvent);
2034                 //                      cout << "rec \t " <<gamma1MotherLabel << "\t" <<  weighted << endl;
2035                                                 }
2036                                         }
2037                                         hESDTruePrimaryMotherInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt(),weighted);
2038                                         hESDTruePrimaryMotherW0WeightingInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt());
2039                                         pESDTruePrimaryMotherWeightsInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt(),weighted);
2040                                         
2041                                         
2042                                         if (fDoMesonQA > 0){
2043                                                 if(isTruePi0){ // Only primary pi0 for resolution
2044                                                         hESDTruePrimaryPi0MCPtResolPt[fiCut]->Fill(((TParticle*)MCStack->Particle(gamma1MotherLabel))->Pt(),(Pi0Candidate->Pt()-((TParticle*)MCStack->Particle(gamma1MotherLabel))->Pt())/((TParticle*)MCStack->Particle(gamma1MotherLabel))->Pt(),weighted);
2045                                                 }
2046                                                 if (isTrueEta){ // Only primary eta for resolution
2047                                                         hESDTruePrimaryEtaMCPtResolPt[fiCut]->Fill(((TParticle*)MCStack->Particle(gamma1MotherLabel))->Pt(),(Pi0Candidate->Pt()-((TParticle*)MCStack->Particle(gamma1MotherLabel))->Pt())/((TParticle*)MCStack->Particle(gamma1MotherLabel))->Pt(),weighted);
2048                                                 }
2049                                         }
2050                                 }
2051                         } else if(!isTruePi0 && !isTrueEta){ // Background
2052                                 if (fDoMesonQA > 0){
2053                                         if(gamma0MotherLabel>-1 && gamma1MotherLabel>-1){ // Both Tracks are Photons and have a mother but not Pi0 or Eta
2054                                                 hESDTrueBckGGInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt());
2055                                                 iMesonMCInfo = 1;
2056                                         } else { // No photon or without mother
2057                                                 hESDTrueBckContInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt());
2058                                         }
2059                                 }
2060                                 if( isTruePi0Dalitz || isTrueEtaDalitz ){
2061                                 // Dalitz
2062                                         iMesonMCInfo = 5;
2063                                         hESDTrueMotherDalitzInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt());
2064                                 } else if (gamma0DalitzCand || gamma1DalitzCand){
2065                                         if (fDoMesonQA > 0)hESDTrueBckContInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt());
2066                                 }   
2067                         }
2068                 }
2069         }
2070 }
2071 //______________________________________________________________________
2072 void AliAnalysisTaskGammaConvV1::ProcessTrueMesonCandidatesAOD(AliAODConversionMother *Pi0Candidate, AliAODConversionPhoton *TrueGammaCandidate0, AliAODConversionPhoton *TrueGammaCandidate1)
2073 {
2074
2075         // Process True Mesons
2076         TClonesArray *AODMCTrackArray = dynamic_cast<TClonesArray*>(fInputEvent->FindListObject(AliAODMCParticle::StdBranchName()));
2077         Bool_t isTruePi0 = kFALSE;
2078         Bool_t isTrueEta = kFALSE;
2079         Bool_t isTruePi0Dalitz = kFALSE;
2080         Bool_t isTrueEtaDalitz = kFALSE;
2081         Bool_t gamma0DalitzCand = kFALSE;
2082         Bool_t gamma1DalitzCand = kFALSE;
2083    
2084         if (AODMCTrackArray!=NULL && TrueGammaCandidate0 != NULL){
2085                 AliAODMCParticle *positiveMC = static_cast<AliAODMCParticle*>(AODMCTrackArray->At(TrueGammaCandidate0->GetMCLabelPositive()));
2086                 AliAODMCParticle *negativeMC = static_cast<AliAODMCParticle*>(AODMCTrackArray->At(TrueGammaCandidate0->GetMCLabelNegative()));
2087
2088                 iMesonMCInfo = 0;
2089                 Int_t gamma0MCLabel = -1;
2090                 Int_t gamma0MotherLabel = -1;
2091                 if(!positiveMC||!negativeMC)
2092                         return;
2093                 
2094                 if(positiveMC->GetMother()>-1&&(negativeMC->GetMother() == positiveMC->GetMother())){
2095                         gamma0MCLabel = positiveMC->GetMother();
2096                 }
2097
2098                 if(gamma0MCLabel != -1){ // Gamma is Combinatorial; MC Particles don't belong to the same Mother
2099                         // Daughters Gamma 0
2100                         AliAODMCParticle * gammaMC0 = static_cast<AliAODMCParticle*>(AODMCTrackArray->At(gamma0MCLabel));
2101                         if(abs(negativeMC->GetPdgCode())==11 && abs(positiveMC->GetPdgCode())==11){  // Electrons ...
2102                                 if(((positiveMC->GetMCProcessCode())) == 5 && ((negativeMC->GetMCProcessCode())) == 5){ // ... From Conversion ...     
2103                                         if(gammaMC0->GetPdgCode() == 22){ // ... with Gamma Mother
2104                                         gamma0MotherLabel=gammaMC0->GetMother();
2105                                         }
2106                                 }
2107                                 if(gammaMC0->GetPdgCode() ==111){ // Dalitz candidate
2108                                         gamma0DalitzCand = kTRUE;
2109                                         gamma0MotherLabel=-111;
2110                                 }
2111                                 if(gammaMC0->GetPdgCode() ==221){ // Dalitz candidate
2112                                         gamma0DalitzCand = kTRUE;
2113                                         gamma0MotherLabel=-221;
2114                                 }
2115                         }
2116                 }
2117                 positiveMC = static_cast<AliAODMCParticle*>(AODMCTrackArray->At(TrueGammaCandidate1->GetMCLabelPositive()));
2118                 negativeMC = static_cast<AliAODMCParticle*>(AODMCTrackArray->At(TrueGammaCandidate1->GetMCLabelNegative()));
2119                 
2120                 Int_t gamma1MCLabel = -1;
2121                 Int_t gamma1MotherLabel = -1;
2122                 if(!positiveMC||!negativeMC)
2123                         return;
2124                 
2125                 if(positiveMC->GetMother()>-1&&(negativeMC->GetMother() == positiveMC->GetMother())){
2126                         gamma1MCLabel = positiveMC->GetMother();
2127                 }
2128                 if(gamma1MCLabel != -1){ // Gamma is Combinatorial; MC Particles don't belong to the same Mother
2129                         // Daughters Gamma 1
2130                         AliAODMCParticle * gammaMC1 = static_cast<AliAODMCParticle*>(AODMCTrackArray->At(gamma1MCLabel));
2131                         if(abs(negativeMC->GetPdgCode())==11 && abs(positiveMC->GetPdgCode())==11){  // Electrons ...
2132                                 if(((positiveMC->GetMCProcessCode())) == 5 && ((negativeMC->GetMCProcessCode())) == 5){ // ... From Conversion ...     
2133                                         if(gammaMC1->GetPdgCode() == 22){ // ... with Gamma Mother
2134                                         gamma1MotherLabel=gammaMC1->GetMother();
2135                                         }
2136                                 }
2137                                 if(gammaMC1->GetPdgCode() ==111 ){ // Dalitz candidate
2138                                                 gamma1DalitzCand = kTRUE;
2139                                                 gamma1MotherLabel=-111;
2140                                 }
2141                                 if(gammaMC1->GetPdgCode() ==221){ // Dalitz candidate
2142                                         gamma1DalitzCand = kTRUE;
2143                                         gamma1MotherLabel=-221;
2144                                 }
2145                         }
2146                 }
2147                 if(gamma0MotherLabel>=0 && gamma0MotherLabel==gamma1MotherLabel){
2148                         if(static_cast<AliAODMCParticle*>(AODMCTrackArray->At(gamma1MotherLabel))->GetPdgCode() == 111){
2149                                 isTruePi0=kTRUE;
2150                         }
2151                         if(static_cast<AliAODMCParticle*>(AODMCTrackArray->At(gamma1MotherLabel))->GetPdgCode() == 221){
2152                                 isTrueEta=kTRUE;
2153                         }
2154                 }
2155                 
2156                 //Identify Dalitz candidate
2157                 if (gamma1DalitzCand || gamma0DalitzCand){
2158                         if (gamma0DalitzCand && gamma0MCLabel >=0 && gamma0MCLabel==gamma1MotherLabel){
2159                                 if (gamma0MotherLabel == -111) isTruePi0Dalitz = kTRUE;
2160                                 if (gamma0MotherLabel == -221) isTrueEtaDalitz = kTRUE;   
2161                         }   
2162                         if (gamma1DalitzCand && gamma1MCLabel >=0 && gamma1MCLabel==gamma0MotherLabel){
2163                                 if (gamma1MotherLabel == -111) isTruePi0Dalitz = kTRUE;
2164                                 if (gamma1MotherLabel == -221) isTrueEtaDalitz = kTRUE;   
2165                         }
2166                 }
2167                                         
2168                 if(isTruePi0 || isTrueEta){// True Pion or Eta
2169                         hESDTrueMotherInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt());
2170                         if (fDoMesonQA > 0){
2171                                 if (isTruePi0){
2172                                         if ( Pi0Candidate->M() > 0.05 && Pi0Candidate->M() < 0.17){
2173                                         hESDTruePi0PtY[fiCut]->Fill(Pi0Candidate->Pt(),Pi0Candidate->Rapidity()-((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetEtaShift());
2174                                         hESDTruePi0PtAlpha[fiCut]->Fill(Pi0Candidate->Pt(),Pi0Candidate->GetAlpha()); 
2175                                         hESDTruePi0PtOpenAngle[fiCut]->Fill(Pi0Candidate->Pt(),Pi0Candidate->GetOpeningAngle()); 
2176                                         }
2177                                 } else if (isTrueEta){   
2178                                         if ( Pi0Candidate->M() > 0.45 && Pi0Candidate->M() < 0.65){
2179                                         hESDTrueEtaPtY[fiCut]->Fill(Pi0Candidate->Pt(),Pi0Candidate->Rapidity()-((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetEtaShift()); 
2180                                         hESDTrueEtaPtAlpha[fiCut]->Fill(Pi0Candidate->Pt(),Pi0Candidate->GetAlpha()); 
2181                                         hESDTrueEtaPtOpenAngle[fiCut]->Fill(Pi0Candidate->Pt(),Pi0Candidate->GetOpeningAngle()); 
2182                                         }
2183                                 }
2184                         }
2185                         if(!(static_cast<AliAODMCParticle*>(AODMCTrackArray->At(gamma0MotherLabel))->IsPrimary())){ // Secondary Meson
2186                                 Int_t secMotherLabel = static_cast<AliAODMCParticle*>(AODMCTrackArray->At(gamma1MotherLabel))->GetMother();
2187                                 Float_t weightedSec= 1;
2188                                 if(((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsParticleFromBGEvent(secMotherLabel, 0x0, fInputEvent) && static_cast<AliAODMCParticle*>(AODMCTrackArray->At(secMotherLabel))->GetPdgCode()==310){
2189                                         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
2190                                         //cout << "MC input \t"<<i << "\t" <<  particle->Pt()<<"\t"<<weighted << endl;
2191                                 }
2192                                 hESDTrueSecondaryMotherInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt(),weightedSec);
2193                                 iMesonMCInfo = 2;   
2194                                 if (secMotherLabel >-1){
2195                                         if(static_cast<AliAODMCParticle*>(AODMCTrackArray->At(secMotherLabel))->GetPdgCode()==310){
2196                                                 iMesonMCInfo = 4;
2197                                                 hESDTrueSecondaryMotherFromK0sInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt(),weightedSec);
2198                                                 if (fDoMesonQA > 0)hESDTrueK0sWithPi0DaughterMCPt[fiCut]
2199                                                                                         ->Fill(static_cast<AliAODMCParticle*>(AODMCTrackArray->At(secMotherLabel))->Pt());
2200                                         }
2201                                         if(static_cast<AliAODMCParticle*>(AODMCTrackArray->At(secMotherLabel))->GetPdgCode()==221){
2202                                                 iMesonMCInfo = 3;
2203                                                 hESDTrueSecondaryMotherFromEtaInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt(),weightedSec);
2204                                                 if (fDoMesonQA > 0)hESDTrueEtaWithPi0DaughterMCPt[fiCut]
2205                                                                                         ->Fill(static_cast<AliAODMCParticle*>(AODMCTrackArray->At(secMotherLabel))->Pt());
2206                                         }
2207                                         if(static_cast<AliAODMCParticle*>(AODMCTrackArray->At(secMotherLabel))->GetPdgCode()==3122){
2208                                                 iMesonMCInfo = 7;
2209                                                 hESDTrueSecondaryMotherFromLambdaInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt(),weightedSec);
2210                                                 if (fDoMesonQA > 0)hESDTrueLambdaWithPi0DaughterMCPt[fiCut]
2211                                                                                         ->Fill(static_cast<AliAODMCParticle*>(AODMCTrackArray->At(secMotherLabel))->Pt());
2212                                         }
2213                                 }
2214                         }else{ // Only primary pi0 for efficiency calculation
2215                                 Float_t weighted= 1;
2216                                 iMesonMCInfo = 6;
2217                                 if(((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsParticleFromBGEvent(gamma1MotherLabel, 0x0, fInputEvent)){
2218                                         if (static_cast<AliAODMCParticle*>(AODMCTrackArray->At(gamma1MotherLabel))->Pt()>0.005){
2219                                                 weighted= ((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetWeightForMeson(fV0Reader->GetPeriodName(),gamma1MotherLabel, 0x0, fInputEvent);
2220                                         //                      cout << "rec \t " <<gamma1MotherLabel << "\t" <<  weighted << endl;
2221                                         }
2222                                 }
2223                                 hESDTruePrimaryMotherInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt(),weighted);
2224                                 hESDTruePrimaryMotherW0WeightingInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt());
2225                                 pESDTruePrimaryMotherWeightsInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt(),weighted);
2226                                         
2227                                 if (fDoMesonQA > 0){
2228                                         if(isTruePi0){ // Only primary pi0 for resolution
2229                                                 hESDTruePrimaryPi0MCPtResolPt[fiCut]->Fill(static_cast<AliAODMCParticle*>(AODMCTrackArray->At(gamma1MotherLabel))->Pt(),
2230                                                                                                                                 (Pi0Candidate->Pt()-static_cast<AliAODMCParticle*>(AODMCTrackArray->At(gamma1MotherLabel))->Pt())/static_cast<AliAODMCParticle*>(AODMCTrackArray->At(gamma1MotherLabel))->Pt(),weighted);
2231                                         
2232                                         }
2233                                         if (isTrueEta){ // Only primary eta for resolution
2234                                                 hESDTruePrimaryEtaMCPtResolPt[fiCut]->Fill(static_cast<AliAODMCParticle*>(AODMCTrackArray->At(gamma1MotherLabel))->Pt(),
2235                                                                                                                                 (Pi0Candidate->Pt()-static_cast<AliAODMCParticle*>(AODMCTrackArray->At(gamma1MotherLabel))->Pt())/static_cast<AliAODMCParticle*>(AODMCTrackArray->At(gamma1MotherLabel))->Pt(),weighted);
2236                                         }
2237                                 }
2238                         }
2239                 } else if(!isTruePi0 && !isTrueEta) { // Background
2240                         if (fDoMesonQA > 0){
2241                                 if(gamma0MotherLabel>-1 && gamma1MotherLabel>-1){ // Both Tracks are Photons and have a mother but not Pi0 or Eta
2242                                         hESDTrueBckGGInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt());
2243                                         iMesonMCInfo = 1;
2244                                 } else { // No photon or without mother
2245                                         hESDTrueBckContInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt());
2246                                 }
2247                         }
2248                         if( isTruePi0Dalitz || isTrueEtaDalitz ){
2249                                 // Dalitz
2250                                 iMesonMCInfo = 5;
2251                                 hESDTrueMotherDalitzInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt());
2252                         } else if (gamma0DalitzCand || gamma1DalitzCand){
2253                                 if (fDoMesonQA > 0)hESDTrueBckContInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt());
2254                         }
2255                 }
2256         }
2257         return;
2258 }
2259 //________________________________________________________________________
2260 void AliAnalysisTaskGammaConvV1::CalculateBackground(){
2261
2262         Int_t zbin= fBGHandler[fiCut]->GetZBinIndex(fInputEvent->GetPrimaryVertex()->GetZ());
2263         Int_t mbin = 0;
2264
2265         if(((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->UseTrackMultiplicity()){
2266                 mbin = fBGHandler[fiCut]->GetMultiplicityBinIndex(fV0Reader->GetNumberOfPrimaryTracks());
2267         } else {
2268                 mbin = fBGHandler[fiCut]->GetMultiplicityBinIndex(fGammaCandidates->GetEntries());
2269         }
2270
2271         if(((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->UseRotationMethod()){
2272
2273                 for(Int_t iCurrent=0;iCurrent<fGammaCandidates->GetEntries();iCurrent++){
2274                         AliAODConversionPhoton currentEventGoodV0 = *(AliAODConversionPhoton*)(fGammaCandidates->At(iCurrent));
2275                         for(Int_t iCurrent2=iCurrent+1;iCurrent2<fGammaCandidates->GetEntries();iCurrent2++){
2276                                 for(Int_t nRandom=0;nRandom<((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->GetNumberOfBGEvents();nRandom++){
2277                                 AliAODConversionPhoton currentEventGoodV02 = *(AliAODConversionPhoton*)(fGammaCandidates->At(iCurrent2));
2278
2279                                 if(((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->DoBGProbability()){
2280                                         AliAODConversionMother *backgroundCandidateProb = new AliAODConversionMother(&currentEventGoodV0,&currentEventGoodV02);
2281                                         Double_t massBGprob = backgroundCandidateProb->M();
2282                                         if(massBGprob>0.1 && massBGprob<0.14){
2283                                                 if(fRandom.Rndm()>fBGHandler[fiCut]->GetBGProb(zbin,mbin)){
2284                                                         delete backgroundCandidateProb;
2285                                                         continue;
2286                                                 }
2287                                         }
2288                                         delete backgroundCandidateProb;
2289                                         backgroundCandidateProb = 0x0;
2290                                 }
2291
2292                                 RotateParticle(&currentEventGoodV02);
2293                                 AliAODConversionMother *backgroundCandidate = new AliAODConversionMother(&currentEventGoodV0,&currentEventGoodV02);
2294                                 backgroundCandidate->CalculateDistanceOfClossetApproachToPrimVtx(fInputEvent->GetPrimaryVertex());
2295                                 if((((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))
2296                                         ->MesonIsSelected(backgroundCandidate,kFALSE,((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetEtaShift()))){
2297                                         hESDMotherBackInvMassPt[fiCut]->Fill(backgroundCandidate->M(),backgroundCandidate->Pt());
2298                                         Double_t sparesFill[4] = {backgroundCandidate->M(),backgroundCandidate->Pt(),(Double_t)zbin,(Double_t)mbin};
2299                                         sESDMotherBackInvMassPtZM[fiCut]->Fill(sparesFill,1);
2300                                 }
2301                                 delete backgroundCandidate;
2302                                 backgroundCandidate = 0x0;
2303                                 }
2304                         }
2305                 }
2306         }else{
2307                 AliGammaConversionAODBGHandler::GammaConversionVertex *bgEventVertex = NULL;
2308
2309                 if(((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->UseTrackMultiplicity()){
2310                         for(Int_t nEventsInBG=0;nEventsInBG<fBGHandler[fiCut]->GetNBGEvents();nEventsInBG++){
2311                                 AliGammaConversionAODVector *previousEventV0s = fBGHandler[fiCut]->GetBGGoodV0s(zbin,mbin,nEventsInBG);
2312                                 if(fMoveParticleAccordingToVertex == kTRUE){
2313                                 bgEventVertex = fBGHandler[fiCut]->GetBGEventVertex(zbin,mbin,nEventsInBG);
2314                                 }
2315
2316                                 for(Int_t iCurrent=0;iCurrent<fGammaCandidates->GetEntries();iCurrent++){
2317                                 AliAODConversionPhoton currentEventGoodV0 = *(AliAODConversionPhoton*)(fGammaCandidates->At(iCurrent));
2318                                 for(UInt_t iPrevious=0;iPrevious<previousEventV0s->size();iPrevious++){
2319                                         AliAODConversionPhoton previousGoodV0 = (AliAODConversionPhoton)(*(previousEventV0s->at(iPrevious)));
2320                                         if(fMoveParticleAccordingToVertex == kTRUE){
2321                                                 MoveParticleAccordingToVertex(&previousGoodV0,bgEventVertex);
2322                                         }