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