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