AliPhysics  master (3d17d9d)
AliAnalysisTaskHeavyNeutralMesonToGG.cxx
Go to the documentation of this file.
1 /**************************************************************************
2  * Copyright(c) 1998-2020, ALICE Experiment at CERN, All rights reserved. *
3  * *
4  * Author: Friederike Bock, Sandro Wenzel *
5  * Version 1.0 *
6  * *
7  * *
8  * Permission to use, copy, modify and distribute this software and its *
9  * documentation strictly for non-commercial purposes is hereby granted *
10  * without fee, provided that the above copyright notice appears in all *
11  * copies and that both the copyright notice and this permission notice *
12  * appear in the supporting documentation. The authors make no claims *
13  * about the suitability of this software for any purpose. It is *
14  * provided "as is" without express or implied warranty. *
15  **************************************************************************/
16 
18 //----------------------------------------------------------------
19 // Class used to do analysis on conversion photons + calo photons
20 //----------------------------------------------------------------
22 #include "TChain.h"
23 #include "TTree.h"
24 #include "TBranch.h"
25 #include "TFile.h"
26 #include "TH1F.h"
27 #include "TH2F.h"
28 #include "TH3F.h"
29 #include "THnSparse.h"
30 #include "TCanvas.h"
31 #include "TNtuple.h"
32 #include "AliAnalysisTask.h"
33 #include "AliAnalysisManager.h"
34 #include "AliESDEvent.h"
35 #include "AliESDInputHandler.h"
36 #include "AliMCEventHandler.h"
37 #include "AliMCEvent.h"
38 #include "AliMCParticle.h"
39 #include "AliCentrality.h"
40 #include "AliESDVZERO.h"
41 #include "AliESDpid.h"
43 #include "AliVParticle.h"
44 #include "AliESDtrack.h"
45 #include "AliESDtrackCuts.h"
46 #include "AliKFVertex.h"
47 #include "AliGenCocktailEventHeader.h"
48 #include "AliConversionAODBGHandlerRP.h"
49 #include "AliAODMCParticle.h"
50 #include "AliAODMCHeader.h"
51 #include "AliEventplane.h"
53 #include "AliAODEvent.h"
54 #include "AliAODInputHandler.h"
55 #include "AliESDEvent.h"
56 #include "AliESDInputHandler.h"
57 #include "AliInputEventHandler.h"
58 #include "AliCaloTrackMatcher.h"
59 #include <vector>
60 #include <map>
61 
63 
64 //________________________________________________________________________
66  fRandom(0),
67  fV0Reader(NULL),
68  fBGHandler(NULL),
69  fBGHandlerRP(NULL),
70  fBGClusHandler(NULL),
71  fBGClusHandlerRP(NULL),
72  fInputEvent(NULL),
73  fMCEvent(NULL),
74  fEventCuts(NULL),
75  fConversionCuts(NULL),
76  fCaloPhotonCuts(NULL),
77  fMesonCuts(NULL),
78  fCutFolder(NULL),
79  fESDList(NULL),
80  fBackList(NULL),
81  fMotherList(NULL),
82  fTrueList(NULL),
83  fMCList(NULL),
84  fOutputContainer(NULL),
85  fGammaCandidates(NULL),
86  fClusterCandidates(NULL),
87  fEventCutArray(NULL),
88  fCutArray(NULL),
89  fClusterCutArray(NULL),
90  fMesonCutArray(NULL),
91  fReaderGammas(NULL),
92  fV0ReaderName("V0ReaderV1"),
93  fCorrTaskSetting(""),
94  fFileNameBroken(NULL),
95  fSparseMotherInvMassPtZM(NULL),
96  fSparseMotherBackInvMassPtZM(NULL),
97  fTreeBrokenFiles(NULL),
98  fProfileTruePrimaryMesonWeightsInvMassPt(NULL),
99  fHistoMotherInvMassPt(NULL),
100  fHistoMotherMatchedInvMassPt(NULL),
101  fHistoMotherBackInvMassPt(NULL),
102  fHistoMotherMesonPtY(NULL),
103  fHistoMotherMesonPtAlpha(NULL),
104  fHistoMotherMesonPtOpenAngle(NULL),
105  fHistoMotherMesonConvPhotonEtaPhi(NULL),
106  fHistoTrueMesonInvMassPt(NULL),
107  fHistoTrueMesonMatchedInvMassPt(NULL),
108  fHistoTrueMesonCaloPhotonInvMassPt(NULL),
109  fHistoTrueMesonCaloConvertedPhotonInvMassPt(NULL),
110  fHistoTrueMesonCaloMixedPhotonConvPhotonInvMassPt(NULL),
111  fHistoTrueMesonCaloConvertedPhotonMatchedInvMassPt(NULL),
112  fHistoTrueMesonCaloElectronInvMassPt(NULL),
113  fHistoTrueMesonCaloMergedClusterInvMassPt(NULL),
114  fHistoTrueMesonCaloMergedClusterPartConvInvMassPt(NULL),
115  fHistoTruePrimaryMesonInvMassPt(NULL),
116  fHistoTruePrimaryMesonW0WeightingInvMassPt(NULL),
117  fHistoTruePrimaryMesonMCPtResolPt(NULL),
118  fHistoTrueMotherMesonConvPhotonEtaPhi(NULL),
119  fHistoTrueBckGGInvMassPt(NULL),
120  fHistoTrueBckFullMesonContainedInOneClusterInvMassPt(NULL),
121  fHistoTrueBckAsymEClustersInvMassPt(NULL),
122  fHistoTrueBckContInvMassPt(NULL),
123  fHistoTrueMesonPtY(NULL),
124  fHistoTrueMesonPtAlpha(NULL),
125  fHistoTrueMesonPtOpenAngle(NULL),
126  fHistoMCMesonPtY(NULL),
127  fHistoMCMesonPtAlpha(NULL),
128  fHistoMCMesonPtJetPt(NULL),
129  fHistoTrueNLabelsInClusPt(NULL),
130  fHistoDoubleCountTrueMesonInvMassPt(NULL),
131  fHistoDoubleCountTrueConvGammaRPt(NULL),
132  fHistoDoubleCountTrueClusterGammaPt(NULL),
133  fHistoSPDClusterTrackletBackground(NULL),
134  fProfileEtaShift(NULL),
135  fProfileJetJetXSection(NULL),
136  fHistoMCHeaders(NULL),
137  fHistoConvGammaPt(NULL),
138  fHistoClusGammaPt(NULL),
139  fHistoClusGammaE(NULL),
140  fHistoClusOverlapHeadersGammaPt(NULL),
141  fHistoClusAllHeadersGammaPt(NULL),
142  fHistoClusRejectedHeadersGammaPt(NULL),
143  fHistoMotherInvMassRejected(NULL),
144  fHistoMCMesonPt(NULL),
145  fHistoMCMesonWOWeightPt(NULL),
146  fHistoMCMesonWOEvtWeightPt(NULL),
147  fHistoMCMesonInAccPt(NULL),
148  fHistoMCMesonWOWeightInAccPt(NULL),
149  fHistoMCMesonWOEvtWeightInAccPt(NULL),
150  fHistoTrueConvGammaPt(NULL),
151  fHistoTruePrimaryConvGammaPt(NULL),
152  fHistoTrueClusGammaPt(NULL),
153  fHistoTrueClusConvGammaPt(NULL),
154  fHistoTrueClusConvGammaFullyPt(NULL),
155  fHistoTruePrimaryClusGammaPt(NULL),
156  fHistoTruePrimaryClusConvGammaPt(NULL),
157  fHistoMultipleCountTrueMeson(NULL),
158  fHistoMultipleCountTrueConvGamma(NULL),
159  fHistoMultipleCountTrueClusterGamma(NULL),
160  fHistoNEvents(NULL),
161  fHistoNEventsWOWeight(NULL),
162  fHistoNGoodESDTracks(NULL),
163  fHistoVertexZ(NULL),
164  fHistoVertexX(NULL),
165  fHistoVertexY(NULL),
166  fHistoNGammaConvCandidates(NULL),
167  fHistoNGammaCaloCandidates(NULL),
168  fHistoNV0Tracks(NULL),
169  fHistoJetJetNTrials(NULL),
170  fMapMultipleCountTrueMesons(),
171  fMapMultipleCountTrueConvGammas(),
172  fMapMultipleCountTrueClusterGammas(),
173  fVectorRecTrueMesons(0),
174  fVectorDoubleCountTrueMesons(0),
175  fVectorDoubleCountTrueConvGammas(0),
176  fVectorDoubleCountTrueClusterGammas(0),
177  fUnsmearedPx(NULL),
178  fUnsmearedPy(NULL),
179  fUnsmearedPz(NULL),
180  fUnsmearedE(NULL),
181  fMesonInvMassWindow(NULL),
182  fMCEventPos(NULL),
183  fMCEventNeg(NULL),
184  fESDArrayPos(NULL),
185  fESDArrayNeg(NULL),
186  fEventPlaneAngle(-100),
187  fMesonInvMassMin(0),
188  fMesonInvMassMax(0),
189  fMesonInvMassNBins(0),
190  fWeightJetJetMC(1),
191  fNGammaCandidates(0),
192  fnCuts(0),
193  fiCut(0),
194  fIsHeavyIon(0),
195  fMesonRecoMode(-1),
196  fMesonType(-1),
197  fMesonPDG(0),
198  fDoMesonQA(0),
199  fDoPhotonQA(0),
200  fDoClusterQA(0),
201  fIsMC(0),
202  fMoveParticleAccordingToVertex(kTRUE),
203  fDoLightOutput(kFALSE),
204  fIsFromDesiredHeader(kTRUE),
205  fIsOverlappingWithOtherHeader(kFALSE),
206  fDoTHnSparse(kTRUE),
207  fSetPlotHistsExtQA(kFALSE),
208  fDoConvGammaShowerShapeTree(kFALSE),
209  fEnableSortForClusMC(kFALSE),
210  fDoPrimaryTrackMatching(kFALSE),
211  fDoInvMassShowerShapeTree(kFALSE),
212  fAllowOverlapHeaders(kTRUE),
213  fEnableClusterCutsForTrigger(kFALSE)
214 {
215 
216 }
217 
218 //________________________________________________________________________
220  AliAnalysisTaskSE(name),
221  fRandom(0),
222  fV0Reader(NULL),
223  fBGHandler(NULL),
224  fBGHandlerRP(NULL),
225  fBGClusHandler(NULL),
226  fBGClusHandlerRP(NULL),
227  fInputEvent(NULL),
228  fMCEvent(NULL),
229  fEventCuts(NULL),
230  fConversionCuts(NULL),
231  fCaloPhotonCuts(NULL),
232  fMesonCuts(NULL),
233  fCutFolder(NULL),
234  fESDList(NULL),
235  fBackList(NULL),
236  fMotherList(NULL),
237  fTrueList(NULL),
238  fMCList(NULL),
239  fOutputContainer(NULL),
240  fGammaCandidates(NULL),
241  fClusterCandidates(NULL),
242  fEventCutArray(NULL),
243  fCutArray(NULL),
244  fClusterCutArray(NULL),
245  fMesonCutArray(NULL),
246  fReaderGammas(NULL),
247  fV0ReaderName("V0ReaderV1"),
248  fCorrTaskSetting(""),
249  fFileNameBroken(NULL),
250  fSparseMotherInvMassPtZM(NULL),
251  fSparseMotherBackInvMassPtZM(NULL),
252  fTreeBrokenFiles(NULL),
253  fProfileTruePrimaryMesonWeightsInvMassPt(NULL),
254  fHistoMotherInvMassPt(NULL),
255  fHistoMotherMatchedInvMassPt(NULL),
256  fHistoMotherBackInvMassPt(NULL),
257  fHistoMotherMesonPtY(NULL),
258  fHistoMotherMesonPtAlpha(NULL),
259  fHistoMotherMesonPtOpenAngle(NULL),
260  fHistoMotherMesonConvPhotonEtaPhi(NULL),
261  fHistoTrueMesonInvMassPt(NULL),
262  fHistoTrueMesonMatchedInvMassPt(NULL),
263  fHistoTrueMesonCaloPhotonInvMassPt(NULL),
264  fHistoTrueMesonCaloConvertedPhotonInvMassPt(NULL),
265  fHistoTrueMesonCaloMixedPhotonConvPhotonInvMassPt(NULL),
266  fHistoTrueMesonCaloConvertedPhotonMatchedInvMassPt(NULL),
267  fHistoTrueMesonCaloElectronInvMassPt(NULL),
268  fHistoTrueMesonCaloMergedClusterInvMassPt(NULL),
269  fHistoTrueMesonCaloMergedClusterPartConvInvMassPt(NULL),
270  fHistoTruePrimaryMesonInvMassPt(NULL),
271  fHistoTruePrimaryMesonW0WeightingInvMassPt(NULL),
272  fHistoTruePrimaryMesonMCPtResolPt(NULL),
273  fHistoTrueMotherMesonConvPhotonEtaPhi(NULL),
274  fHistoTrueBckGGInvMassPt(NULL),
275  fHistoTrueBckFullMesonContainedInOneClusterInvMassPt(NULL),
276  fHistoTrueBckAsymEClustersInvMassPt(NULL),
277  fHistoTrueBckContInvMassPt(NULL),
278  fHistoTrueMesonPtY(NULL),
279  fHistoTrueMesonPtAlpha(NULL),
280  fHistoTrueMesonPtOpenAngle(NULL),
281  fHistoMCMesonPtY(NULL),
282  fHistoMCMesonPtAlpha(NULL),
283  fHistoMCMesonPtJetPt(NULL),
284  fHistoTrueNLabelsInClusPt(NULL),
285  fHistoDoubleCountTrueMesonInvMassPt(NULL),
286  fHistoDoubleCountTrueConvGammaRPt(NULL),
287  fHistoDoubleCountTrueClusterGammaPt(NULL),
288  fHistoSPDClusterTrackletBackground(NULL),
289  fProfileEtaShift(NULL),
290  fProfileJetJetXSection(NULL),
291  fHistoMCHeaders(NULL),
292  fHistoConvGammaPt(NULL),
293  fHistoClusGammaPt(NULL),
294  fHistoClusGammaE(NULL),
295  fHistoClusOverlapHeadersGammaPt(NULL),
296  fHistoClusAllHeadersGammaPt(NULL),
297  fHistoClusRejectedHeadersGammaPt(NULL),
298  fHistoMotherInvMassRejected(NULL),
299  fHistoMCMesonPt(NULL),
300  fHistoMCMesonWOWeightPt(NULL),
301  fHistoMCMesonWOEvtWeightPt(NULL),
302  fHistoMCMesonInAccPt(NULL),
303  fHistoMCMesonWOWeightInAccPt(NULL),
304  fHistoMCMesonWOEvtWeightInAccPt(NULL),
305  fHistoTrueConvGammaPt(NULL),
306  fHistoTruePrimaryConvGammaPt(NULL),
307  fHistoTrueClusGammaPt(NULL),
308  fHistoTrueClusConvGammaPt(NULL),
309  fHistoTrueClusConvGammaFullyPt(NULL),
310  fHistoTruePrimaryClusGammaPt(NULL),
311  fHistoTruePrimaryClusConvGammaPt(NULL),
312  fHistoMultipleCountTrueMeson(NULL),
313  fHistoMultipleCountTrueConvGamma(NULL),
314  fHistoMultipleCountTrueClusterGamma(NULL),
315  fHistoNEvents(NULL),
316  fHistoNEventsWOWeight(NULL),
317  fHistoNGoodESDTracks(NULL),
318  fHistoVertexZ(NULL),
319  fHistoVertexX(NULL),
320  fHistoVertexY(NULL),
321  fHistoNGammaConvCandidates(NULL),
322  fHistoNGammaCaloCandidates(NULL),
323  fHistoNV0Tracks(NULL),
324  fHistoJetJetNTrials(NULL),
325  fMapMultipleCountTrueMesons(),
326  fMapMultipleCountTrueConvGammas(),
327  fMapMultipleCountTrueClusterGammas(),
328  fVectorRecTrueMesons(0),
329  fVectorDoubleCountTrueMesons(0),
330  fVectorDoubleCountTrueConvGammas(0),
331  fVectorDoubleCountTrueClusterGammas(0),
332  fUnsmearedPx(NULL),
333  fUnsmearedPy(NULL),
334  fUnsmearedPz(NULL),
335  fUnsmearedE(NULL),
336  fMesonInvMassWindow(NULL),
337  fMCEventPos(NULL),
338  fMCEventNeg(NULL),
339  fESDArrayPos(NULL),
340  fESDArrayNeg(NULL),
341  fEventPlaneAngle(-100),
342  fMesonInvMassMin(0),
343  fMesonInvMassMax(0),
344  fMesonInvMassNBins(0),
345  fWeightJetJetMC(1),
346  fNGammaCandidates(0),
347  fnCuts(0),
348  fiCut(0),
349  fIsHeavyIon(0),
350  fMesonRecoMode(-1),
351  fMesonType(-1),
352  fMesonPDG(0),
353  fDoMesonQA(0),
354  fDoPhotonQA(0),
355  fDoClusterQA(0),
356  fIsMC(0),
357  fMoveParticleAccordingToVertex(kTRUE),
358  fDoLightOutput(kFALSE),
359  fIsFromDesiredHeader(kTRUE),
360  fIsOverlappingWithOtherHeader(kFALSE),
361  fDoTHnSparse(kTRUE),
362  fSetPlotHistsExtQA(kFALSE),
363  fDoConvGammaShowerShapeTree(kFALSE),
364  fEnableSortForClusMC(kFALSE),
365  fDoPrimaryTrackMatching(kFALSE),
366  fDoInvMassShowerShapeTree(kFALSE),
367  fAllowOverlapHeaders(kTRUE),
368  fEnableClusterCutsForTrigger(kFALSE)
369 {
370  // Define output slots here
371  DefineOutput(1, TList::Class());
372 }
373 
375 {
376  if(fGammaCandidates){
377  delete fGammaCandidates;
378  fGammaCandidates = 0x0;
379  }
380  if(fClusterCandidates){
381  delete fClusterCandidates;
382  fClusterCandidates = 0x0;
383  }
384  if(fBGHandler){
385  delete[] fBGHandler;
386  fBGHandler = 0x0;
387  }
388  if(fBGHandlerRP){
389  delete[] fBGHandlerRP;
390  fBGHandlerRP = 0x0;
391  }
392  if(fBGClusHandler){
393  delete[] fBGClusHandler;
394  fBGClusHandler = 0x0;
395  }
396  if(fBGClusHandlerRP){
397  delete[] fBGClusHandlerRP;
398  fBGClusHandlerRP = 0x0;
399  }
400 }
401 //___________________________________________________________
403 
404  const Int_t nDim = 4;
405  Int_t nBins[nDim] = {800,300,7,4};
406  Double_t xMin[nDim] = {0,0, 0,0};
407  Double_t xMax[nDim] = {1.2,30,7,4};
408 
409  if(fDoTHnSparse){
410  fSparseMotherInvMassPtZM = new THnSparseF*[fnCuts];
411  fSparseMotherBackInvMassPtZM = new THnSparseF*[fnCuts];
412  }
413 
415  fBGHandlerRP = new AliConversionAODBGHandlerRP*[fnCuts];
416 
418  fBGClusHandlerRP = new AliConversionAODBGHandlerRP*[fnCuts];
419 
420  for(Int_t iCut = 0; iCut<fnCuts;iCut++){
421  if (((AliConversionMesonCuts*)fMesonCutArray->At(iCut))->DoBGCalculation()){
422  TString cutstringEvent = ((AliConvEventCuts*)fEventCutArray->At(iCut))->GetCutNumber();
423  TString cutstringConvGamma = "";
424  TString cutstringCaloGamma = "";
425  if (fMesonRecoMode < 2) cutstringConvGamma = ((AliConversionPhotonCuts*)fCutArray->At(iCut))->GetCutNumber();
426  if (fMesonRecoMode > 0 || fEnableClusterCutsForTrigger ) cutstringCaloGamma = ((AliCaloPhotonCuts*)fClusterCutArray->At(iCut))->GetCutNumber();
427  TString cutstringMeson = ((AliConversionMesonCuts*)fMesonCutArray->At(iCut))->GetCutNumber();
428 
429  TString fullCutString = "";
430  if (fMesonRecoMode == 0) fullCutString = Form("%i_%s_%s_%s",fMesonRecoMode, cutstringEvent.Data(), cutstringConvGamma.Data(), cutstringMeson.Data());
431  else if (fMesonRecoMode == 1) fullCutString = Form("%i_%s_%s_%s_%s",fMesonRecoMode, cutstringEvent.Data(),cutstringConvGamma.Data(),cutstringCaloGamma.Data(), cutstringMeson.Data());
432  else if (fMesonRecoMode == 2) fullCutString = Form("%i_%s_%s_%s",fMesonRecoMode, cutstringEvent.Data(), cutstringCaloGamma.Data(), cutstringMeson.Data());
433 
434  Int_t collisionSystem = atoi((TString)(((AliConvEventCuts*)fEventCutArray->At(iCut))->GetCutNumber())(0,1));
435  Int_t centMin = atoi((TString)(((AliConvEventCuts*)fEventCutArray->At(iCut))->GetCutNumber())(1,1));
436  Int_t centMax = atoi((TString)(((AliConvEventCuts*)fEventCutArray->At(iCut))->GetCutNumber())(2,1));
437 
438  if( collisionSystem == 1 || collisionSystem == 2 ||
439  collisionSystem == 5 || collisionSystem == 8 ||
440  collisionSystem == 9){
441  centMin = centMin*10;
442  centMax = centMax*10;
443  if(centMax ==0 && centMax!=centMin) centMax=100;
444  }else if(collisionSystem == 3 || collisionSystem == 6){
445  centMin = centMin*5;
446  centMax = centMax*5;
447  }else if(collisionSystem == 4 || collisionSystem == 7){
448  centMin = ((centMin*5)+45);
449  centMax = ((centMax*5)+45);
450  }
451 
452  if(fDoTHnSparse){
453  fBackList[iCut] = new TList();
454  fBackList[iCut]->SetName(Form("%s Back histograms",fullCutString.Data()));
455  fBackList[iCut]->SetOwner(kTRUE);
456  fCutFolder[iCut]->Add(fBackList[iCut]);
457 
458  fSparseMotherBackInvMassPtZM[iCut] = new THnSparseF("Back_Back_InvMass_Pt_z_m", "Back_Back_InvMass_Pt_z_m",nDim,nBins,xMin,xMax);
459  fBackList[iCut]->Add(fSparseMotherBackInvMassPtZM[iCut]);
460 
461  fMotherList[iCut] = new TList();
462  fMotherList[iCut]->SetName(Form("%s Mother histograms",fullCutString.Data()));
463  fMotherList[iCut]->SetOwner(kTRUE);
464  fCutFolder[iCut]->Add(fMotherList[iCut]);
465 
466  fSparseMotherInvMassPtZM[iCut] = new THnSparseF("Back_Mother_InvMass_Pt_z_m", "Back_Mother_InvMass_Pt_z_m",nDim,nBins,xMin,xMax);
467  fMotherList[iCut]->Add(fSparseMotherInvMassPtZM[iCut]);
468  }
469 
470  if(((AliConversionMesonCuts*)fMesonCutArray->At(iCut))->BackgroundHandlerType() == 0){
471  fBGHandler[iCut] = new AliGammaConversionAODBGHandler( collisionSystem,centMin,centMax,
472  ((AliConversionMesonCuts*)fMesonCutArray->At(iCut))->GetNumberOfBGEvents(),
473  ((AliConversionMesonCuts*)fMesonCutArray->At(iCut))->UseTrackMultiplicity(),
474  2,8,5);
475  fBGClusHandler[iCut] = new AliGammaConversionAODBGHandler( collisionSystem,centMin,centMax,
476  ((AliConversionMesonCuts*)fMesonCutArray->At(iCut))->GetNumberOfBGEvents(),
477  ((AliConversionMesonCuts*)fMesonCutArray->At(iCut))->UseTrackMultiplicity(),
478  2,8,5);
479  fBGHandlerRP[iCut] = NULL;
480  }else{
481  fBGHandlerRP[iCut] = new AliConversionAODBGHandlerRP(((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsHeavyIon(),
482  ((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->UseTrackMultiplicity(),
483  ((AliConversionMesonCuts*)fMesonCutArray->At(iCut))->GetNumberOfBGEvents());
484  fBGClusHandlerRP[iCut] = new AliConversionAODBGHandlerRP( ((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsHeavyIon(),
485  ((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->UseTrackMultiplicity(),
486  ((AliConversionMesonCuts*)fMesonCutArray->At(iCut))->GetNumberOfBGEvents());
487  fBGHandler[iCut] = NULL;
488  }
489  }
490  }
491 }
492 
493 //________________________________________________________________________
495 
496  fV0Reader = (AliV0ReaderV1*)AliAnalysisManager::GetAnalysisManager()->GetTask(fV0ReaderName.Data());
497  if(!fV0Reader){printf("Error: No V0 Reader");return;}// GetV0Reader
498  if (fIsMC == 2){
499  fDoClusterQA = 0;
500  fDoTHnSparse = kFALSE;
501  } else if (fIsMC == 3){
502  fDoTHnSparse = kFALSE;
503  }
504 
505  if (fMesonRecoMode == 2)
507 
508  // Create histograms
509  if(fOutputContainer != NULL){
510  delete fOutputContainer;
511  fOutputContainer = NULL;
512  }
513  if(fOutputContainer == NULL){
514  fOutputContainer = new TList();
515  fOutputContainer->SetOwner(kTRUE);
516  }
517 
518  // Array of current cut's gammas
519  fGammaCandidates = new TList();
520  fClusterCandidates = new TList();
521  fClusterCandidates->SetOwner(kTRUE);
522 
523  fCutFolder = new TList*[fnCuts];
524  fESDList = new TList*[fnCuts];
525 
526  if(fDoTHnSparse){
527  fBackList = new TList*[fnCuts];
528  fMotherList = new TList*[fnCuts];
529  }
530 
531  fHistoNEvents = new TH1F*[fnCuts];
532  if (fIsMC > 1){
533  fHistoNEventsWOWeight = new TH1F*[fnCuts];
534  }
535  if (fIsMC == 2){
536  fProfileJetJetXSection = new TProfile*[fnCuts];
537  fHistoJetJetNTrials = new TH1F*[fnCuts];
538  }
539  fHistoNGoodESDTracks = new TH1F*[fnCuts];
540  fHistoVertexZ = new TH1F*[fnCuts];
541  if(!fDoLightOutput){
542  fHistoVertexX = new TH1F*[fnCuts];
543  fHistoVertexY = new TH1F*[fnCuts];
544  }
545  if (fMesonRecoMode < 2) fHistoNGammaConvCandidates = new TH1F*[fnCuts];
546  if(fIsHeavyIon==2) fProfileEtaShift = new TProfile*[fnCuts];
547  if(!fDoLightOutput){
549  fHistoNV0Tracks = new TH1F*[fnCuts];
550  if (fMesonRecoMode < 2)
551  fHistoConvGammaPt = new TH1F*[fnCuts];
552  }
553 
555  fHistoNGammaCaloCandidates = new TH1F*[fnCuts];
556  if (!fDoLightOutput ){
557  fHistoClusGammaPt = new TH1F*[fnCuts];
558  fHistoClusGammaE = new TH1F*[fnCuts];
559  if (fIsMC > 0){
561  fHistoClusAllHeadersGammaPt = new TH1F*[fnCuts];
563  }
564  }
565  }
566 
568  fHistoMotherInvMassRejected = new TH1F*[fnCuts];
570  if(!fDoLightOutput && fMesonRecoMode == 1){
572  }
573  if (fDoMesonQA > 0){
578  }
579 
580  fMesonInvMassWindow = new Double_t[2];
581  if (fMesonType < 0 || fMesonType > 2){
582  if(!fV0Reader){printf("Error: No V0 Reader");return;}// GetV0Reader
583  } else if (fMesonType == 0){ // pi0 case 134.9770 ± 0.0005 MeV
584  fMesonPDG = 111;
585  fMesonInvMassMin = 0.;
586  fMesonInvMassMax = 0.400;
587  fMesonInvMassNBins = 400;
588  fMesonInvMassWindow[0] = 0.05;
589  fMesonInvMassWindow[1] = 0.17;
590  } else if (fMesonType == 1){ // eta case 547.862 ± 0.017 MeV
591  fMesonPDG = 221;
592  fMesonInvMassMin = 0.300;
593  fMesonInvMassMax = 0.800;
594  fMesonInvMassNBins = 500;
595  fMesonInvMassWindow[0] = 0.45;
596  fMesonInvMassWindow[1] = 0.65;
597  } else if (fMesonType == 2){ // eta' case 957.78 ± 0.06 MeV
598  fMesonPDG = 331;
599  fMesonInvMassMin = 0.700;
600  fMesonInvMassMax = 1.200;
601  fMesonInvMassNBins = 500;
602  fMesonInvMassWindow[0] = 0.85;
603  fMesonInvMassWindow[1] = 1.05;
604  }
605  // set common binning in pT for mesons and photons
606  Int_t nBinsPt = 200;
607  Float_t binWidthPt = 0.1;
608  Float_t minPt = 0;
609  Float_t maxPt = 20;
610  Int_t nBinsQAPt = 170;
611  Float_t maxQAPt = 20;
612  Int_t nBinsClusterPt = 500;
613  Float_t minClusterPt = 0;
614  Float_t maxClusterPt = 50;
615  Double_t *arrPtBinning = new Double_t[1200];
616  Double_t *arrQAPtBinning = new Double_t[1200];
617  Double_t *arrClusPtBinning = new Double_t[1200];
618  // Set special pt binning for pp 8TeV
619  if (((AliConvEventCuts*)fV0Reader->GetEventCuts())->GetEnergyEnum() == AliConvEventCuts::k8TeV ){
620  nBinsQAPt = 190;
621  maxQAPt = 40;
622  for(Int_t i=0; i<nBinsQAPt+1;i++){
623  if(i<60) arrQAPtBinning[i] = 0.05*i;
624  else if(i<130) arrQAPtBinning[i] = 3.+0.1*(i-60);
625  else if(i<170) arrQAPtBinning[i] = 10.+0.25*(i-130);
626  else if(i<190) arrQAPtBinning[i] = 20.+1.0*(i-170);
627  else arrQAPtBinning[i] = maxQAPt;
628  }
629  nBinsPt = 400;
630  minPt = 0;
631  maxPt = 40;
632  for(Int_t i=0; i<nBinsPt+1;i++){
633  arrPtBinning[i] = ((maxPt-minPt)/nBinsPt)*i;
634  }
635  nBinsClusterPt = 800;
636  minClusterPt = 0;
637  maxClusterPt = 80;
638  for(Int_t i=0; i<nBinsPt+1;i++){
639  arrClusPtBinning[i] = ((maxClusterPt-minClusterPt)/nBinsClusterPt)*i;
640  }
641  // Set special pt binning for pPb 5TeV
642  } else if ( ((AliConvEventCuts*)fV0Reader->GetEventCuts())->GetEnergyEnum() == AliConvEventCuts::kpPb5TeV ||
643  ((AliConvEventCuts*)fV0Reader->GetEventCuts())->GetEnergyEnum() == AliConvEventCuts::kpPb5TeVR2 ||
644  ((AliConvEventCuts*)fV0Reader->GetEventCuts())->GetEnergyEnum() == AliConvEventCuts::k5TeV ){
645  binWidthPt = 0.05;
646  nBinsPt = 205;
647  minPt = 0;
648  maxPt = 60;
649  for(Int_t i=0; i<nBinsPt+1;i++){
650  if (i < 1) arrPtBinning[i] = 0.3*i;
651  else if(i<55) arrPtBinning[i] = 0.3+0.05*(i-1);
652  else if(i<125) arrPtBinning[i] = 3.+0.1*(i-55);
653  else if(i<165) arrPtBinning[i] = 10.+0.25*(i-125);
654  else if(i<205) arrPtBinning[i] = 20.+1.0*(i-165);
655  else arrPtBinning[i] = maxPt;
656  }
657  nBinsQAPt = 210;
658  maxQAPt = 60;
659  for(Int_t i=0; i<nBinsQAPt+1;i++){
660  if(i<60) arrQAPtBinning[i] = 0.05*i;
661  else if(i<130) arrQAPtBinning[i] = 3.+0.1*(i-60);
662  else if(i<170) arrQAPtBinning[i] = 10.+0.25*(i-130);
663  else if(i<210) arrQAPtBinning[i] = 20.+1.0*(i-170);
664  else arrQAPtBinning[i] = maxQAPt;
665  }
666  nBinsClusterPt = 301;
667  minClusterPt = 0;
668  maxClusterPt = 100;
669  for(Int_t i=0; i<nBinsClusterPt+1;i++){
670  if (i < 1) arrClusPtBinning[i] = 0.3*i;
671  else if(i<55) arrClusPtBinning[i] = 0.3+0.05*(i-1);
672  else if(i<125) arrClusPtBinning[i] = 3.+0.1*(i-55);
673  else if(i<155) arrClusPtBinning[i] = 10.+0.2*(i-125);
674  else if(i<211) arrClusPtBinning[i] = 16.+0.25*(i-155);
675  else if(i<251) arrClusPtBinning[i] = 30.+0.5*(i-211);
676  else if(i<301) arrClusPtBinning[i] = 50.+1.0*(i-251);
677  else arrClusPtBinning[i] = maxClusterPt;
678  }
679  // Set special pt binning for pp 13TeV, pPb 8TeV
680  } else if ( ((AliConvEventCuts*)fV0Reader->GetEventCuts())->GetEnergyEnum() == AliConvEventCuts::k13TeV ||
681  ((AliConvEventCuts*)fV0Reader->GetEventCuts())->GetEnergyEnum() == AliConvEventCuts::k13TeVLowB ||
682  ((AliConvEventCuts*)fV0Reader->GetEventCuts())->GetEnergyEnum() == AliConvEventCuts::kpPb8TeV ){
683  nBinsPt = 285;
684  minPt = 0;
685  maxPt = 100;
686  binWidthPt = 0.05;
687  for(Int_t i=0; i<nBinsPt+1;i++){
688  if (i < 1) arrPtBinning[i] = 0.3*i;
689  else if(i<55) arrPtBinning[i] = 0.3+0.05*(i-1);
690  else if(i<125) arrPtBinning[i] = 3.+0.1*(i-55);
691  else if(i<185) arrPtBinning[i] = 10.+0.25*(i-125);
692  else if(i<235) arrPtBinning[i] = 25.+0.5*(i-185);
693  else if(i<285) arrPtBinning[i] = 50.+1.0*(i-235);
694  else arrPtBinning[i] = maxPt;
695  }
696  nBinsQAPt = 270;
697  maxQAPt = 100;
698  for(Int_t i=0; i<nBinsQAPt+1;i++){
699  if(i<60) arrQAPtBinning[i] = 0.05*i;
700  else if(i<130) arrQAPtBinning[i] = 3.+0.1*(i-60);
701  else if(i<170) arrQAPtBinning[i] = 10.+0.25*(i-130);
702  else if(i<210) arrQAPtBinning[i] = 20.+0.5*(i-170);
703  else if(i<270) arrQAPtBinning[i] = 40.+1.0*(i-210);
704  else arrQAPtBinning[i] = maxQAPt;
705  }
706  nBinsClusterPt = 301;
707  minClusterPt = 0;
708  maxClusterPt = 100;
709  for(Int_t i=0; i<nBinsClusterPt+1;i++){
710  if (i < 1) arrClusPtBinning[i] = 0.3*i;
711  else if(i<55) arrClusPtBinning[i] = 0.3+0.05*(i-1);
712  else if(i<125) arrClusPtBinning[i] = 3.+0.1*(i-55);
713  else if(i<155) arrClusPtBinning[i] = 10.+0.2*(i-125);
714  else if(i<211) arrClusPtBinning[i] = 16.+0.25*(i-155);
715  else if(i<251) arrClusPtBinning[i] = 30.+0.5*(i-211);
716  else if(i<301) arrClusPtBinning[i] = 50.+1.0*(i-251);
717  else arrClusPtBinning[i] = maxClusterPt;
718  }
719  } else if (((AliConvEventCuts*)fV0Reader->GetEventCuts())->GetEnergyEnum() == AliConvEventCuts::kXeXe5440GeV ){
720  nBinsPt = 90;
721  minPt = 0;
722  maxPt = 20;
723  for(Int_t i=0; i<nBinsPt+1;i++){
724  if (i < 1) arrPtBinning[i] = 0.3*i;
725  else if(i<58) arrPtBinning[i] = 0.3+0.1*(i-1);
726  else if(i<82) arrPtBinning[i] = 6.+0.25*(i-58);
727  else if(i<90) arrPtBinning[i] = 12.+1.0*(i-82);
728  else arrPtBinning[i] = maxPt;
729  }
730  nBinsQAPt = 92;
731  maxQAPt = 20;
732  for(Int_t i=0; i<nBinsQAPt+1;i++){
733  if(i<60) arrQAPtBinning[i] = 0.1*i;
734  else if(i<84) arrQAPtBinning[i] = 6.+0.25*(i-60);
735  else if(i<92) arrQAPtBinning[i] = 12.+1.0*(i-84);
736  else arrQAPtBinning[i] = maxQAPt;
737  }
738  nBinsClusterPt = 148;
739  minClusterPt = 0;
740  maxClusterPt = 40;
741  for(Int_t i=0; i<nBinsClusterPt+1;i++){
742  if (i < 1) arrClusPtBinning[i] = 0.3*i;
743  else if(i<98) arrClusPtBinning[i] = 0.3+0.1*(i-1);
744  else if(i<123) arrClusPtBinning[i] = 10.+0.2*(i-98);
745  else if(i<148) arrClusPtBinning[i] = 15.+1.0*(i-123);
746  else arrClusPtBinning[i] = maxClusterPt;
747  }
748  } else if (((AliConvEventCuts*)fV0Reader->GetEventCuts())->GetEnergyEnum() == AliConvEventCuts::kPbPb5TeV ){
749  nBinsPt = 90;
750  minPt = 0;
751  maxPt = 20;
752  for(Int_t i=0; i<nBinsPt+1;i++){
753  if (i < 1) arrPtBinning[i] = 0.3*i;
754  else if(i<58) arrPtBinning[i] = 0.3+0.1*(i-1);
755  else if(i<82) arrPtBinning[i] = 6.+0.25*(i-58);
756  else if(i<90) arrPtBinning[i] = 12.+1.0*(i-82);
757  else arrPtBinning[i] = maxPt;
758  }
759  nBinsQAPt = 92;
760  maxQAPt = 20;
761  for(Int_t i=0; i<nBinsQAPt+1;i++){
762  if(i<60) arrQAPtBinning[i] = 0.1*i;
763  else if(i<84) arrQAPtBinning[i] = 6.+0.25*(i-60);
764  else if(i<92) arrQAPtBinning[i] = 12.+1.0*(i-84);
765  else arrQAPtBinning[i] = maxQAPt;
766  }
767  nBinsClusterPt = 148;
768  minClusterPt = 0;
769  maxClusterPt = 40;
770  for(Int_t i=0; i<nBinsClusterPt+1;i++){
771  if (i < 1) arrClusPtBinning[i] = 0.3*i;
772  else if(i<98) arrClusPtBinning[i] = 0.3+0.1*(i-1);
773  else if(i<123) arrClusPtBinning[i] = 10.+0.2*(i-98);
774  else if(i<148) arrClusPtBinning[i] = 15.+1.0*(i-123);
775  else arrClusPtBinning[i] = maxClusterPt;
776  }
777  // default binning
778  } else {
779  for(Int_t i=0; i<nBinsPt+1;i++){
780  arrPtBinning[i] = ((maxPt-minPt)/nBinsPt)*i;
781  }
782  for(Int_t i=0; i<nBinsClusterPt+1;i++){
783  arrClusPtBinning[i] = ((maxClusterPt-minClusterPt)/nBinsClusterPt)*i;
784  }
785  for(Int_t i=0; i<nBinsQAPt+1;i++){
786  if(i<60) arrQAPtBinning[i] = 0.05*i;
787  else if(i<130) arrQAPtBinning[i] = 3.+0.1*(i-60);
788  else if(i<170) arrQAPtBinning[i] = 10.+0.25*(i-130);
789  else arrQAPtBinning[i] = maxQAPt;
790  }
791  }
792 
793  for(Int_t iCut = 0; iCut<fnCuts;iCut++){
794  TString cutstringEvent = ((AliConvEventCuts*)fEventCutArray->At(iCut))->GetCutNumber();
795  TString cutstringConvGamma = "";
796  TString cutstringCaloGamma = "";
797  if (fMesonRecoMode < 2) cutstringConvGamma = ((AliConversionPhotonCuts*)fCutArray->At(iCut))->GetCutNumber();
798  if (fMesonRecoMode > 0 || fEnableClusterCutsForTrigger) cutstringCaloGamma = ((AliCaloPhotonCuts*)fClusterCutArray->At(iCut))->GetCutNumber();
799  TString cutstringMeson = ((AliConversionMesonCuts*)fMesonCutArray->At(iCut))->GetCutNumber();
800 
801  TString fullCutString = "";
802  if (fMesonRecoMode == 0) fullCutString = Form("%i_%s_%s_%s",fMesonRecoMode, cutstringEvent.Data(), cutstringConvGamma.Data(), cutstringMeson.Data());
803  else if (fMesonRecoMode == 1) fullCutString = Form("%i_%s_%s_%s_%s",fMesonRecoMode, cutstringEvent.Data(),cutstringConvGamma.Data(),cutstringCaloGamma.Data(), cutstringMeson.Data());
804  else if (fMesonRecoMode == 2) fullCutString = Form("%i_%s_%s_%s",fMesonRecoMode, cutstringEvent.Data(), cutstringCaloGamma.Data(), cutstringMeson.Data());
805 
806  fCutFolder[iCut] = new TList();
807  fCutFolder[iCut]->SetName(Form("Cut Number %s",fullCutString.Data()));
808  fCutFolder[iCut]->SetOwner(kTRUE);
809  fOutputContainer->Add(fCutFolder[iCut]);
810  fESDList[iCut] = new TList();
811  fESDList[iCut]->SetName(Form("%s ESD histograms",fullCutString.Data()));
812  fESDList[iCut]->SetOwner(kTRUE);
813  fCutFolder[iCut]->Add(fESDList[iCut]);
814 
815  fHistoNEvents[iCut] = new TH1F("NEvents", "NEvents", 14, -0.5, 13.5);
816  fHistoNEvents[iCut]->GetXaxis()->SetBinLabel(1,"Accepted");
817  fHistoNEvents[iCut]->GetXaxis()->SetBinLabel(2,"Centrality");
818  fHistoNEvents[iCut]->GetXaxis()->SetBinLabel(3,"Miss. MC or inc. ev.");
819  if (((AliConvEventCuts*)fEventCutArray->At(iCut))->IsSpecialTrigger() > 1 ){
820  TString TriggerNames = "Not Trigger: ";
821  TriggerNames = TriggerNames+ ( (AliConvEventCuts*)fEventCutArray->At(iCut))->GetSpecialTriggerName();
822  fHistoNEvents[iCut]->GetXaxis()->SetBinLabel(4,TriggerNames.Data());
823  }else {
824  fHistoNEvents[iCut]->GetXaxis()->SetBinLabel(4,"Trigger");
825  }
826  fHistoNEvents[iCut]->GetXaxis()->SetBinLabel(5,"Vertex Z");
827  fHistoNEvents[iCut]->GetXaxis()->SetBinLabel(6,"Cont. Vertex");
828  fHistoNEvents[iCut]->GetXaxis()->SetBinLabel(7,"Pile-Up");
829  fHistoNEvents[iCut]->GetXaxis()->SetBinLabel(8,"no SDD");
830  fHistoNEvents[iCut]->GetXaxis()->SetBinLabel(9,"no V0AND");
831  fHistoNEvents[iCut]->GetXaxis()->SetBinLabel(10,"EMCAL problem");
832  fHistoNEvents[iCut]->GetXaxis()->SetBinLabel(11,"rejectedForJetJetMC");
833  fHistoNEvents[iCut]->GetXaxis()->SetBinLabel(12,"SPD hits vs tracklet");
834  fHistoNEvents[iCut]->GetXaxis()->SetBinLabel(13,"Out-of-Bunch pileup Past-Future");
835  fHistoNEvents[iCut]->GetXaxis()->SetBinLabel(14,"Pileup V0M-TPCout Tracks");
836  fESDList[iCut]->Add(fHistoNEvents[iCut]);
837 
838  if (fIsMC > 1){
839  fHistoNEventsWOWeight[iCut] = new TH1F("NEventsWOWeight", "NEventsWOWeight", 14, -0.5, 13.5);
840  fHistoNEventsWOWeight[iCut]->GetXaxis()->SetBinLabel(1,"Accepted");
841  fHistoNEventsWOWeight[iCut]->GetXaxis()->SetBinLabel(2,"Centrality");
842  fHistoNEventsWOWeight[iCut]->GetXaxis()->SetBinLabel(3,"Miss. MC or inc. ev.");
843  if (((AliConvEventCuts*)fEventCutArray->At(iCut))->IsSpecialTrigger() > 1 ){
844  TString TriggerNames = "Not Trigger: ";
845  TriggerNames = TriggerNames+ ( (AliConvEventCuts*)fEventCutArray->At(iCut))->GetSpecialTriggerName();
846  fHistoNEventsWOWeight[iCut]->GetXaxis()->SetBinLabel(4,TriggerNames.Data());
847  }else {
848  fHistoNEventsWOWeight[iCut]->GetXaxis()->SetBinLabel(4,"Trigger");
849  }
850  fHistoNEventsWOWeight[iCut]->GetXaxis()->SetBinLabel(5,"Vertex Z");
851  fHistoNEventsWOWeight[iCut]->GetXaxis()->SetBinLabel(6,"Cont. Vertex");
852  fHistoNEventsWOWeight[iCut]->GetXaxis()->SetBinLabel(7,"Pile-Up");
853  fHistoNEventsWOWeight[iCut]->GetXaxis()->SetBinLabel(8,"no SDD");
854  fHistoNEventsWOWeight[iCut]->GetXaxis()->SetBinLabel(9,"no V0AND");
855  fHistoNEventsWOWeight[iCut]->GetXaxis()->SetBinLabel(10,"EMCAL problem");
856  fHistoNEventsWOWeight[iCut]->GetXaxis()->SetBinLabel(11,"rejectedForJetJetMC");
857  fHistoNEventsWOWeight[iCut]->GetXaxis()->SetBinLabel(12,"SPD hits vs tracklet");
858  fHistoNEventsWOWeight[iCut]->GetXaxis()->SetBinLabel(13,"Out-of-Bunch pileup Past-Future");
859  fHistoNEventsWOWeight[iCut]->GetXaxis()->SetBinLabel(14,"Pileup V0M-TPCout Tracks");
860  fESDList[iCut]->Add(fHistoNEventsWOWeight[iCut]);
861  }
862 
863  if (fIsMC == 2){
864  fProfileJetJetXSection[iCut] = new TProfile("XSection", "XSection", 1, -0.5, 0.5);
865  fESDList[iCut]->Add(fProfileJetJetXSection[iCut]);
866  fHistoJetJetNTrials[iCut] = new TH1F("NTrials", "#sum{NTrials}", 1, 0, 1);
867  fHistoJetJetNTrials[iCut]->GetXaxis()->SetBinLabel(1,"#sum{NTrials}");
868  fESDList[iCut]->Add(fHistoJetJetNTrials[iCut]);
869  }
870 
871  if(fIsHeavyIon == 1)
872  fHistoNGoodESDTracks[iCut] = new TH1F("GoodESDTracks", "GoodESDTracks; # TPC tracks", 4000, 0, 4000);
873  else if(fIsHeavyIon == 2)
874  fHistoNGoodESDTracks[iCut] = new TH1F("GoodESDTracks", "GoodESDTracks; # TPC tracks", 400, 0, 400);
875  else
876  fHistoNGoodESDTracks[iCut] = new TH1F("GoodESDTracks", "GoodESDTracks; # TPC tracks", 200, 0, 200);
877  fESDList[iCut]->Add(fHistoNGoodESDTracks[iCut]);
878 
879  fHistoVertexZ[iCut] = new TH1F("VertexZ", "VertexZ", 200, -10, 10);
880  fESDList[iCut]->Add(fHistoVertexZ[iCut]);
881  if(!fDoLightOutput){
882  fHistoVertexX[iCut] = new TH1F("VertexX", "VertexX", 100, -5, 5);
883  fESDList[iCut]->Add(fHistoVertexX[iCut]);
884  fHistoVertexY[iCut] = new TH1F("VertexY", "VertexY", 100, -5, 5);
885  fESDList[iCut]->Add(fHistoVertexY[iCut]);
886  }
887 
888  if (fMesonRecoMode < 2){
889  if(fIsHeavyIon == 1)
890  fHistoNGammaConvCandidates[iCut] = new TH1F("GammaConvCandidates", "GammaConvCandidates; # accepted #gamma_{conv}", 100, 0, 100);
891  else if(fIsHeavyIon == 2)
892  fHistoNGammaConvCandidates[iCut] = new TH1F("GammaConvCandidates", "GammaConvCandidates; # accepted #gamma_{conv}", 50, 0, 50);
893  else
894  fHistoNGammaConvCandidates[iCut] = new TH1F("GammaConvCandidates", "GammaConvCandidates; # accepted #gamma_{conv}", 50, 0, 50);
895  fESDList[iCut]->Add(fHistoNGammaConvCandidates[iCut]);
896  }
898  if(fIsHeavyIon == 1)
899  fHistoNGammaCaloCandidates[iCut] = new TH1F("GammaCaloCandidates", "GammaCaloCandidates; # accepted #gamma_{conv}", 100, 0, 100);
900  else if(fIsHeavyIon == 2)
901  fHistoNGammaCaloCandidates[iCut] = new TH1F("GammaCaloCandidates", "GammaCaloCandidates; # accepted #gamma_{conv}", 50, 0, 50);
902  else
903  fHistoNGammaCaloCandidates[iCut] = new TH1F("GammaCaloCandidates", "GammaCaloCandidates; # accepted #gamma_{conv}", 50, 0, 50);
904  fESDList[iCut]->Add(fHistoNGammaCaloCandidates[iCut]);
905  }
906 
907  if(!fDoLightOutput){
908  fHistoSPDClusterTrackletBackground[iCut] = new TH2F("SPD tracklets vs SPD clusters", "SPD tracklets vs SPD clusters", 100, 0, 200, 250, 0, 1000);
910 
911  if(fIsHeavyIon == 1)
912  fHistoNV0Tracks[iCut] = new TH1F("V0 Multiplicity", "V0 Multiplicity; VZERO amp [arb. units]", 30000, 0, 30000);
913  else if(fIsHeavyIon == 2)
914  fHistoNV0Tracks[iCut] = new TH1F("V0 Multiplicity", "V0 Multiplicity; VZERO amp [arb. units]", 2500, 0, 2500);
915  else
916  fHistoNV0Tracks[iCut] = new TH1F("V0 Multiplicity", "V0 Multiplicity; VZERO amp [arb. units]", 1500, 0, 1500);
917  fESDList[iCut]->Add(fHistoNV0Tracks[iCut]);
918 
919  if (fMesonRecoMode < 2){
920  fHistoConvGammaPt[iCut] = new TH1F("ESD_ConvGamma_Pt", "ESD_ConvGamma_Pt; p_{T,conv}(GeV/c)", (Int_t)((maxPt-minPt)/binWidthPt), minPt, maxPt);
921  fESDList[iCut]->Add(fHistoConvGammaPt[iCut]);
922  }
924  fHistoClusGammaPt[iCut] = new TH1F("ClusGamma_Pt", "ClusGamma_Pt; p_{T,clus} (GeV/c)", nBinsClusterPt, arrClusPtBinning);
925  fESDList[iCut]->Add(fHistoClusGammaPt[iCut]);
926  fHistoClusGammaE[iCut] = new TH1F("ClusGamma_E", "ClusGamma_E; E_{clus} (GeV)", nBinsClusterPt, arrClusPtBinning);
927  fESDList[iCut]->Add(fHistoClusGammaE[iCut]);
928  if (fIsMC > 0){
929  fHistoClusOverlapHeadersGammaPt[iCut] = new TH1F("ClusGammaOverlapHeaders_Pt", "ClusGammaOverlapHeaders_Pt; p_{T,clus} (GeV/c), selected header w/ overlap",
930  nBinsClusterPt, arrClusPtBinning);
931  fESDList[iCut]->Add(fHistoClusOverlapHeadersGammaPt[iCut]);
932  fHistoClusAllHeadersGammaPt[iCut] = new TH1F("ClusGammaAllHeaders_Pt", "ClusGammaAllHeaders_Pt; p_{T,clus} (GeV/c), all headers",
933  nBinsClusterPt, arrClusPtBinning);
934  fESDList[iCut]->Add(fHistoClusAllHeadersGammaPt[iCut]);
935  fHistoClusRejectedHeadersGammaPt[iCut] = new TH1F("ClusGammaRejectedHeaders_Pt", "ClusGammaRejectedHeaders_Pt; p_{T,clus} (GeV/c), rejected headers",
936  nBinsClusterPt, arrClusPtBinning);
937  fESDList[iCut]->Add(fHistoClusRejectedHeadersGammaPt[iCut]);
938  }
939  }
940  }
941 
942  if(fIsHeavyIon == 2){
943  fProfileEtaShift[iCut] = new TProfile("Eta Shift", "Eta Shift", 1, -0.5, 0.5);
944  fESDList[iCut]->Add(fProfileEtaShift[iCut]);
945  }
946 
947  if (fIsMC > 1){
948  fHistoNEvents[iCut]->Sumw2();
949  fHistoNGoodESDTracks[iCut]->Sumw2();
950  fHistoVertexZ[iCut]->Sumw2();
951  if (fMesonRecoMode < 2) fHistoNGammaConvCandidates[iCut]->Sumw2();
953  if(!fDoLightOutput){
954  fHistoVertexX[iCut]->Sumw2();
955  fHistoVertexY[iCut]->Sumw2();
956  fHistoSPDClusterTrackletBackground[iCut]->Sumw2();
957  fHistoNV0Tracks[iCut]->Sumw2();
958  if (fMesonRecoMode < 2)
959  fHistoConvGammaPt[iCut]->Sumw2();
961  if (fHistoClusGammaPt[iCut]) fHistoClusGammaPt[iCut]->Sumw2();
962  if (fHistoClusGammaE[iCut]) fHistoClusGammaE[iCut]->Sumw2();
966  }
967  }
968  }
969 
970 
971  fHistoMotherInvMassPt[iCut] = new TH2F("ESD_Mother_InvMass_Pt", "ESD_Mother_InvMass_Pt; M_{inv} (GeV/c^{2}); p_{T,pair} (GeV/c)",
972  fMesonInvMassNBins, fMesonInvMassMin, fMesonInvMassMax, nBinsPt, arrPtBinning);
973  fESDList[iCut]->Add(fHistoMotherInvMassPt[iCut]);
974  fHistoMotherInvMassRejected[iCut] = new TH1F("ESD_Mother_InvMassRejected", "ESD_Mother_InvMassRejected; M_{inv} (GeV/c^{2})",
975  1200, 0, 1.2);
976  fESDList[iCut]->Add(fHistoMotherInvMassRejected[iCut]);
977 
978  if(!fDoLightOutput){
979  if (fMesonRecoMode == 1){
980  fHistoMotherMatchedInvMassPt[iCut] = new TH2F("ESD_MotherMatched_InvMass_Pt", "ESD_MotherMatched_InvMass_Pt; M_{inv} (GeV/c^{2}) matched conv e^{+/-}to cluster; p_{T,pair} (GeV/c)",
981  fMesonInvMassNBins, fMesonInvMassMin, fMesonInvMassMax, nBinsPt, arrPtBinning);
982  fESDList[iCut]->Add(fHistoMotherMatchedInvMassPt[iCut]);
983  }
984  }
985 
986  fHistoMotherBackInvMassPt[iCut] = new TH2F("ESD_Background_InvMass_Pt", "ESD_Background_InvMass_Pt; M_{inv, mxed}(GeV/c^{2}); p_{T,BG pair} (GeV/c)",
987  fMesonInvMassNBins, fMesonInvMassMin, fMesonInvMassMax, nBinsPt, arrPtBinning);
988  fESDList[iCut]->Add(fHistoMotherBackInvMassPt[iCut]);
989 
990  if (fIsMC > 1){
991  fHistoMotherInvMassPt[iCut]->Sumw2();
992  fHistoMotherBackInvMassPt[iCut]->Sumw2();
993  if(!fDoLightOutput && fMesonRecoMode == 1){
995  }
996  }
997 
998  if (fDoMesonQA > 0 ){
999  fHistoMotherMesonPtY[iCut] = new TH2F("ESD_MotherMeson_Pt_Y", "ESD_MotherMeson_Pt_Y; p_{T, meson cand} (GeV/c); y_{meson cand}",
1000  nBinsQAPt, arrQAPtBinning, 150, -1.5, 1.5);
1001  fESDList[iCut]->Add(fHistoMotherMesonPtY[iCut]);
1002  fHistoMotherMesonPtAlpha[iCut] = new TH2F("ESD_MotherMeson_Pt_Alpha", "ESD_MotherMeson_Pt_Alpha; p_{T, meson cand} (GeV/c); #alpha_{meson cand}",
1003  nBinsQAPt, arrQAPtBinning, 200, -1, 1);
1004  fESDList[iCut]->Add(fHistoMotherMesonPtAlpha[iCut]);
1005  fHistoMotherMesonPtOpenAngle[iCut] = new TH2F("ESD_MotherMeson_Pt_OpenAngle", "ESD_MotherMeson_Pt_OpenAngle; p_{T, meson cand} (GeV/c); #theta_{meson cand}",
1006  nBinsQAPt, arrQAPtBinning, 100, 0, 1);
1007  fESDList[iCut]->Add(fHistoMotherMesonPtOpenAngle[iCut]);
1008  fHistoMotherMesonConvPhotonEtaPhi[iCut] = new TH2F("ESD_MotherMesonConvPhoton_Eta_Phi", "ConvPhoton under meson peak; #phi_{#gamma_{conv}}(rad); #eta_{#gamma_{conv}}",
1009  600, 0, 2*TMath::Pi(), 200, -1, 1);
1010  fESDList[iCut]->Add(fHistoMotherMesonConvPhotonEtaPhi[iCut]);
1011  if (fIsMC > 1){
1012  fHistoMotherMesonPtY[iCut]->Sumw2();
1013  fHistoMotherMesonPtAlpha[iCut]->Sumw2();
1014  fHistoMotherMesonPtOpenAngle[iCut]->Sumw2();
1015  fHistoMotherMesonConvPhotonEtaPhi[iCut]->Sumw2();
1016  }
1017  }
1018  }
1019 
1020  InitBack(); // Init Background Handler
1021 
1022  if(fIsMC>0){
1023  // MC Histogramms
1024  fMCList = new TList*[fnCuts];
1025  // True Histogramms
1026  fTrueList = new TList*[fnCuts];
1027 
1028  if(!fDoLightOutput){
1029  fHistoMCHeaders = new TH1I*[fnCuts];
1030  if (fMesonRecoMode < 2){
1031  fHistoTrueConvGammaPt = new TH1F*[fnCuts];
1034  fHistoTruePrimaryConvGammaPt = new TH1F*[fnCuts];
1035  }
1036 
1038  fHistoTrueClusGammaPt = new TH1F*[fnCuts];
1039  fHistoTrueClusConvGammaPt = new TH1F*[fnCuts];
1040  fHistoTruePrimaryClusGammaPt = new TH1F*[fnCuts];
1046  }
1047  }
1048 
1049  fHistoMCMesonPt = new TH1F*[fnCuts];
1050  fHistoMCMesonWOWeightPt = new TH1F*[fnCuts];
1051  fHistoMCMesonInAccPt = new TH1F*[fnCuts];
1052  fHistoMCMesonWOWeightInAccPt = new TH1F*[fnCuts];
1053  if (fIsMC > 1){
1054  fHistoMCMesonWOEvtWeightPt = new TH1F*[fnCuts];
1056  }
1057 
1060  fHistoMultipleCountTrueMeson = new TH1F*[fnCuts];
1064  if (fMesonRecoMode == 1){
1066  }
1067 
1068  if (fDoMesonQA > 0){
1069  fHistoMCMesonPtY = new TH2F*[fnCuts];
1071  if (fIsMC == 2){
1073  }
1074 
1075  if (fIsMC < 2){
1076  if (fMesonRecoMode > 0){
1082  }
1083  if (fMesonRecoMode == 1){
1086  }
1087  if (fMesonRecoMode == 2){
1089  }
1091  }
1092  if(fDoMesonQA > 1){
1097  }
1098  fHistoTrueMesonPtY = new TH2F*[fnCuts];
1101  }
1102 
1103  for(Int_t iCut = 0; iCut<fnCuts;iCut++){
1104  TString cutstringEvent = ((AliConvEventCuts*)fEventCutArray->At(iCut))->GetCutNumber();
1105  TString cutstringConvGamma = "";
1106  TString cutstringCaloGamma = "";
1107  if (fMesonRecoMode < 2) cutstringConvGamma = ((AliConversionPhotonCuts*)fCutArray->At(iCut))->GetCutNumber();
1108  if (fMesonRecoMode > 0 || fEnableClusterCutsForTrigger) cutstringCaloGamma = ((AliCaloPhotonCuts*)fClusterCutArray->At(iCut))->GetCutNumber();
1109  TString cutstringMeson = ((AliConversionMesonCuts*)fMesonCutArray->At(iCut))->GetCutNumber();
1110 
1111  TString fullCutString = "";
1112  if (fMesonRecoMode == 0) fullCutString = Form("%i_%s_%s_%s",fMesonRecoMode, cutstringEvent.Data(), cutstringConvGamma.Data(), cutstringMeson.Data());
1113  else if (fMesonRecoMode == 1) fullCutString = Form("%i_%s_%s_%s_%s",fMesonRecoMode, cutstringEvent.Data(),cutstringConvGamma.Data(),cutstringCaloGamma.Data(), cutstringMeson.Data());
1114  else if (fMesonRecoMode == 2) fullCutString = Form("%i_%s_%s_%s",fMesonRecoMode, cutstringEvent.Data(), cutstringCaloGamma.Data(), cutstringMeson.Data());
1115 
1116  fMCList[iCut] = new TList();
1117  fMCList[iCut]->SetName(Form("%s MC histograms",fullCutString.Data()));
1118  fMCList[iCut]->SetOwner(kTRUE);
1119  fCutFolder[iCut]->Add(fMCList[iCut]);
1120  if(!fDoLightOutput){
1121  fHistoMCHeaders[iCut] = new TH1I("MC_Headers", "MC_Headers", 20, 0, 20);
1122  fMCList[iCut]->Add(fHistoMCHeaders[iCut]);
1123  }
1124 
1125  fHistoMCMesonPt[iCut] = new TH1F("MC_Meson_Pt", "MC_Meson_Pt; p_{T} (GeV/c)",
1126  (Int_t)((maxPt-minPt)/binWidthPt), minPt, maxPt);
1127  fHistoMCMesonPt[iCut]->Sumw2();
1128  fMCList[iCut]->Add(fHistoMCMesonPt[iCut]);
1129  fHistoMCMesonWOWeightPt[iCut] = new TH1F("MC_Meson_WOWeights_Pt", "MC_Meson_WOWeights_Pt; p_{T} (GeV/c)",
1130  (Int_t)((maxPt-minPt)/binWidthPt), minPt, maxPt);
1131  fMCList[iCut]->Add(fHistoMCMesonWOWeightPt[iCut]);
1132 
1133  fHistoMCMesonInAccPt[iCut] = new TH1F("MC_MesonInAcc_Pt", "MC_MesonInAcc_Pt; p_{T} (GeV/c)",
1134  (Int_t)((maxPt-minPt)/binWidthPt), minPt, maxPt);
1135  fHistoMCMesonInAccPt[iCut]->Sumw2();
1136  fMCList[iCut]->Add(fHistoMCMesonInAccPt[iCut]);
1137  fHistoMCMesonWOWeightInAccPt[iCut] = new TH1F("MC_MesonWOWeightInAcc_Pt", "MC_MesonWOWeightInAcc_Pt; p_{T} (GeV/c)",
1138  (Int_t)((maxPt-minPt)/binWidthPt), minPt, maxPt);
1139  fMCList[iCut]->Add(fHistoMCMesonWOWeightInAccPt[iCut]);
1140 
1141  if (fIsMC > 1){
1142  fHistoMCMesonWOWeightPt[iCut]->Sumw2();
1143  fHistoMCMesonWOWeightInAccPt[iCut]->Sumw2();
1144  fHistoMCMesonWOEvtWeightPt[iCut] = new TH1F("MC_Meson_WOEventWeights_Pt", "MC_Meson_WOEventWeights_Pt; p_{T} (GeV/c)",
1145  (Int_t)((maxPt-minPt)/binWidthPt), minPt, maxPt);
1146  fMCList[iCut]->Add(fHistoMCMesonWOEvtWeightPt[iCut]);
1147  fHistoMCMesonWOEvtWeightInAccPt[iCut] = new TH1F("MC_Meson_WOEventWeightsInAcc_Pt", "MC_Meson_WOEventWeightsInAcc_Pt; p_{T} (GeV/c)",
1148  (Int_t)((maxPt-minPt)/binWidthPt), minPt, maxPt);
1149  fMCList[iCut]->Add(fHistoMCMesonWOEvtWeightInAccPt[iCut]);
1150 
1151  if (fDoMesonQA > 0 && fIsMC == 2){
1152  fHistoMCMesonPtJetPt[iCut] = new TH2F("MC_Meson_Pt_JetPt", "MC_Meson_Pt_JetPt; p_{T} (GeV/c); p_{T,jet} (GeV/c)",
1153  nBinsQAPt, arrQAPtBinning, 200, 0, 200);
1154  fHistoMCMesonPtJetPt[iCut]->Sumw2();
1155  fMCList[iCut]->Add(fHistoMCMesonPtJetPt[iCut]);
1156  }
1157  }
1158 
1159  // book additional MC QA histograms
1160  if (fDoMesonQA > 0){
1161  fHistoMCMesonPtY[iCut] = new TH2F("MC_Meson_Pt_Y", "MC_Meson_Pt_Y; p_{T} (GeV/c); Y ",
1162  nBinsQAPt, arrQAPtBinning, 150, -1.5, 1.5);
1163  fHistoMCMesonPtY[iCut]->Sumw2();
1164  fMCList[iCut]->Add(fHistoMCMesonPtY[iCut]);
1165  fHistoMCMesonPtAlpha[iCut] = new TH2F("MC_Meson_Pt_Alpha", "MC_Meson_Pt_Alpha; p_{T} (GeV/c); #alpha",
1166  nBinsQAPt, arrQAPtBinning, 200, -1, 1);
1167  fMCList[iCut]->Add(fHistoMCMesonPtAlpha[iCut]);
1168 
1169  if (fIsMC == 2){
1170  fHistoMCMesonPtAlpha[iCut]->Sumw2();
1171  }
1172  }
1173 
1174  fTrueList[iCut] = new TList();
1175  fTrueList[iCut]->SetName(Form("%s True histograms",fullCutString.Data()));
1176  fTrueList[iCut]->SetOwner(kTRUE);
1177  fCutFolder[iCut]->Add(fTrueList[iCut]);
1178 
1179  if(!fDoLightOutput){
1180  if (fMesonRecoMode < 2){
1181  fHistoTrueConvGammaPt[iCut] = new TH1F("ESD_TrueConvGamma_Pt", "ESD_TrueConvGamma_Pt; p_{T,conv} (GeV/c); counts",
1182  (Int_t)((maxPt-minPt)/binWidthPt), minPt, maxPt);
1183  fTrueList[iCut]->Add(fHistoTrueConvGammaPt[iCut]);
1184  fHistoDoubleCountTrueConvGammaRPt[iCut] = new TH2F("ESD_TrueDoubleCountConvGamma_R_Pt", "ESD_TrueDoubleCountConvGamma_R_Pt; R (cm); p_{T,conv} (GeV/c)",
1185  800, 0, 200, (Int_t)((maxPt-minPt)/binWidthPt), minPt, maxPt);
1186  fTrueList[iCut]->Add(fHistoDoubleCountTrueConvGammaRPt[iCut]);
1187  fHistoMultipleCountTrueConvGamma[iCut] = new TH1F("ESD_TrueMultipleCountConvGamma", "ESD_TrueMultipleCountConvGamma", 10, 1, 11);
1188  fTrueList[iCut]->Add(fHistoMultipleCountTrueConvGamma[iCut]);
1189  fHistoTruePrimaryConvGammaPt[iCut] = new TH1F("ESD_TruePrimaryConvGamma_Pt", "ESD_TruePrimaryConvGamma_Pt;p_{T,conv} (GeV/c); counts", (Int_t)((maxPt-minPt)/binWidthPt), minPt, maxPt);
1190  fTrueList[iCut]->Add(fHistoTruePrimaryConvGammaPt[iCut]);
1191  }
1193  fHistoTrueClusGammaPt[iCut] = new TH1F("TrueClusGamma_Pt", "ESD_TrueClusGamma_Pt; p_{T,clus} (GeV/c); counts", nBinsClusterPt, arrClusPtBinning);
1194  fTrueList[iCut]->Add(fHistoTrueClusGammaPt[iCut]);
1195  fHistoTruePrimaryClusGammaPt[iCut] = new TH1F("TruePrimaryClusGamma_Pt", "ESD_TruePrimaryClusGamma_Pt; p_{T,clus} (GeV/c); counts", nBinsClusterPt, arrClusPtBinning);
1196  fTrueList[iCut]->Add(fHistoTruePrimaryClusGammaPt[iCut]);
1197  fHistoTrueClusConvGammaPt[iCut] = new TH1F("TrueClusConvGamma_Pt", "TrueClusConvGamma_Pt; p_{T,clus} (GeV/c); counts", nBinsClusterPt, arrClusPtBinning);
1198  fTrueList[iCut]->Add(fHistoTrueClusConvGammaPt[iCut]);
1199  fHistoTruePrimaryClusConvGammaPt[iCut] = new TH1F("TruePrimaryClusConvGamma_Pt", "ESD_TruePrimaryClusConvGamma_Pt; p_{T,clus} (GeV/c); counts", nBinsClusterPt, arrClusPtBinning);
1200  fTrueList[iCut]->Add(fHistoTruePrimaryClusConvGammaPt[iCut]);
1201  fHistoDoubleCountTrueClusterGammaPt[iCut] = new TH2F("TrueDoubleCountClusterGamma_Pt", "TrueDoubleCountClusterGamma_Pt; p_{T,clus} (GeV/c); counts",
1202  nBinsClusterPt, arrClusPtBinning, 2, 0, 2);
1204  fHistoMultipleCountTrueClusterGamma[iCut] = new TH1F("TrueMultipleCountClusterGamma", "TrueMultipleCountClusterGamma", 10, 1, 11);
1206  fHistoTrueClusConvGammaFullyPt[iCut] = new TH1F("TrueClusConvGammaFullyContained_Pt", "TrueClusConvGammaFullyContained_Pt; p_{T,clus} (GeV/c); counts",
1207  nBinsClusterPt, arrClusPtBinning);
1208  fTrueList[iCut]->Add(fHistoTrueClusConvGammaFullyPt[iCut]);
1209  fHistoTrueNLabelsInClusPt[iCut] = new TH2F("TrueNLabelsInClus_Pt", "TrueNLabelsInClus_Pt; p_{T,clus} (GeV/c); counts",
1210  100, -0.5, 99.5, nBinsClusterPt, arrClusPtBinning);
1211  fTrueList[iCut]->Add(fHistoTrueNLabelsInClusPt[iCut]);
1212  }
1213  }
1214 
1215  if (fIsMC > 1){
1216  if(!fDoLightOutput){
1217  if (fMesonRecoMode < 2){
1218  fHistoTrueConvGammaPt[iCut]->Sumw2();
1219  fHistoDoubleCountTrueConvGammaRPt[iCut]->Sumw2();
1220  fHistoMultipleCountTrueConvGamma[iCut]->Sumw2();
1221  }
1223  fHistoTruePrimaryConvGammaPt[iCut]->Sumw2();
1224  }
1226  fHistoTruePrimaryClusGammaPt[iCut]->Sumw2();
1227  fHistoTrueNLabelsInClusPt[iCut]->Sumw2();
1228  fHistoTrueClusConvGammaPt[iCut]->Sumw2();
1229  fHistoTruePrimaryClusConvGammaPt[iCut]->Sumw2();
1230  fHistoTrueClusGammaPt[iCut]->Sumw2();
1231  fHistoDoubleCountTrueClusterGammaPt[iCut]->Sumw2();
1232  fHistoMultipleCountTrueClusterGamma[iCut]->Sumw2();
1233  fHistoTrueClusConvGammaFullyPt[iCut]->Sumw2();
1234  }
1235  }
1236 
1237  }
1238 
1239  fHistoTrueMesonInvMassPt[iCut] = new TH2F("ESD_TrueMeson_InvMass_Pt", "ESD_TrueMeson_InvMass_Pt;M_{inv}(GeV/c^{2});p_{T}(GeV/c)",
1240  fMesonInvMassNBins, fMesonInvMassMin, fMesonInvMassMax, nBinsPt, arrPtBinning);
1241  fHistoTrueMesonInvMassPt[iCut]->Sumw2();
1242  fTrueList[iCut]->Add(fHistoTrueMesonInvMassPt[iCut]);
1243 
1244  fHistoDoubleCountTrueMesonInvMassPt[iCut] = new TH2F("ESD_TrueDoubleCountMeson_InvMass_Pt", "ESD_TrueDoubleCountMeson_InvMass_Pt;M_{inv}(GeV/c^{2});p_{T}(GeV/c)",
1245  fMesonInvMassNBins, fMesonInvMassMin, fMesonInvMassMax, nBinsPt, arrPtBinning);
1247  fHistoMultipleCountTrueMeson[iCut] = new TH1F("ESD_TrueMultipleCountMeson", "ESD_TrueMultipleCountMeson;# number of multiple counts;#", 10, 1, 11);
1248  fHistoMultipleCountTrueMeson[iCut]->Sumw2();
1249  fTrueList[iCut]->Add(fHistoMultipleCountTrueMeson[iCut]);
1250 
1251  fHistoTruePrimaryMesonInvMassPt[iCut] = new TH2F("ESD_TruePrimaryMeson_InvMass_Pt", "ESD_TruePrimaryMeson_InvMass_Pt;M_{inv}(GeV/c^{2});p_{T}(GeV/c)",
1252  fMesonInvMassNBins, fMesonInvMassMin, fMesonInvMassMax, nBinsPt, arrPtBinning);
1253  fHistoTruePrimaryMesonInvMassPt[iCut]->Sumw2();
1254  fTrueList[iCut]->Add(fHistoTruePrimaryMesonInvMassPt[iCut]);
1255 
1256  fHistoTruePrimaryMesonW0WeightingInvMassPt[iCut] = new TH2F("ESD_TruePrimaryMesonW0Weights_InvMass_Pt",
1257  "ESD_TruePrimaryMesonW0Weights_InvMass_Pt;M_{inv}(GeV/c^{2});p_{T}(GeV/c)",
1258  fMesonInvMassNBins, fMesonInvMassMin, fMesonInvMassMax, nBinsPt, arrPtBinning);
1260 
1261  fProfileTruePrimaryMesonWeightsInvMassPt[iCut] = new TProfile2D("ESD_TruePrimaryMesonWeights_InvMass_Pt",
1262  "ESD_TruePrimaryMesonWeights_InvMass_Pt;M_{inv}(GeV/c^{2});p_{T}(GeV/c)",
1263  fMesonInvMassNBins, fMesonInvMassMin, fMesonInvMassMax, nBinsPt, arrPtBinning);
1266 
1267  if (fMesonRecoMode == 1){
1268  fHistoTrueMesonMatchedInvMassPt[iCut] = new TH2F("ESD_TrueMeson_Matched_InvMass_Pt", "ESD_TrueMeson_Matched_InvMass_Pt;M_{inv}(GeV/c^{2});p_{T}(GeV/c)",
1269  fMesonInvMassNBins, fMesonInvMassMin, fMesonInvMassMax, nBinsPt, arrPtBinning);
1270  fHistoTrueMesonMatchedInvMassPt[iCut]->Sumw2();
1271  fTrueList[iCut]->Add(fHistoTrueMesonMatchedInvMassPt[iCut]);
1272  }
1273 
1274  if (fIsMC > 1){
1276  }
1277 
1278  if (fDoMesonQA > 0){
1279  if (fIsMC < 2){
1280  if (fMesonRecoMode > 0){
1281  fHistoTrueMesonCaloPhotonInvMassPt[iCut] = new TH2F( "ESD_TrueMesonCaloPhoton_InvMass_Pt",
1282  "ESD_TrueMesonCaloPhoton_InvMass_Pt; M_{inv}(GeV/c^{2}) #gamma #gamma;p_{T}(GeV/c) ",
1283  fMesonInvMassNBins, fMesonInvMassMin, fMesonInvMassMax, nBinsPt, arrPtBinning);
1285  fHistoTrueMesonCaloConvertedPhotonInvMassPt[iCut] = new TH2F( "ESD_TrueMesonCaloConvertedPhoton_InvMass_Pt",
1286  "ESD_TrueMesonCaloConvertedPhoton_InvMass_Pt;M_{inv}(GeV/c^{2}) #gamma #gamma_{conv};p_{T}(GeV/c)",
1287  fMesonInvMassNBins, fMesonInvMassMin, fMesonInvMassMax, nBinsPt, arrPtBinning);
1289  fHistoTrueMesonCaloElectronInvMassPt[iCut] = new TH2F("ESD_TrueMesonCaloElectron_InvMass_Pt",
1290  "ESD_TrueMesonCaloElectron_InvMass_Pt; M_{inv}(GeV/c^{2}) #gamma e^{#pm}; ; p_{T}(GeV/c)",
1291  fMesonInvMassNBins, fMesonInvMassMin, fMesonInvMassMax, nBinsPt, arrPtBinning);
1293  fHistoTrueMesonCaloMergedClusterInvMassPt[iCut] = new TH2F("ESD_TrueMesonCaloMergedCluster_InvMass_Pt",
1294  "ESD_TrueMesonCaloMergedCluster_InvMass_Pt; M_{inv}(GeV/c^{2}) #gamma merged cluster; p_{T}(GeV/c)",
1295  fMesonInvMassNBins, fMesonInvMassMin, fMesonInvMassMax, nBinsPt, arrPtBinning);
1297  fHistoTrueMesonCaloMergedClusterPartConvInvMassPt[iCut] = new TH2F( "ESD_TrueMesonCaloMergedClusterPartConv_InvMass_Pt",
1298  "ESD_TrueMesonCaloMergedClusterPartConv_InvMass_Pt; M_{inv}(GeV/c^{2}) #gamma merged cluster, part conv; p_{T}(GeV/c)",
1299  fMesonInvMassNBins, fMesonInvMassMin, fMesonInvMassMax, nBinsPt, arrPtBinning);
1301  }
1302  if (fMesonRecoMode == 1){
1303  fHistoTrueMesonCaloConvertedPhotonMatchedInvMassPt[iCut] = new TH2F("ESD_TrueMesonCaloConvertedPhotonMatched_InvMass_Pt",
1304  "ESD_TrueMesonCaloConvertedPhotonMatched_InvMass_Pt;M_{inv}(GeV/c^{2}) #gamma #gamma_{conv,matched};p_{T}(GeV/c)",
1305  fMesonInvMassNBins, fMesonInvMassMin, fMesonInvMassMax, nBinsPt, arrPtBinning);
1307  fHistoTrueMotherMesonConvPhotonEtaPhi[iCut] = new TH2F("ESD_TrueMotherMesonConvPhoton_Eta_Phi", "conv photons for true ; #phi_{#gamma_{conv}}(rad);#eta_{#gamma_{conv}}",
1308  600, 0, 2*TMath::Pi(), 200, -1, 1);
1310  }
1311  if (fMesonRecoMode == 2){
1312  fHistoTrueMesonCaloMixedPhotonConvPhotonInvMassPt[iCut] = new TH2F("ESD_TrueMesonCaloMixedPhotonConvertedPhoton_InvMass_Pt",
1313  "ESD_TrueMesonCaloMixedPhotonConvertedPhoton_InvMass_Pt;M_{inv}(GeV/c^{2}) #gamma #gamma_{conv,matched};p_{T}(GeV/c)",
1314  fMesonInvMassNBins, fMesonInvMassMin, fMesonInvMassMax, nBinsPt, arrPtBinning);
1316  }
1317  fHistoTruePrimaryMesonMCPtResolPt[iCut] = new TH2F("ESD_TruePrimaryMeson_MCPt_ResolPt",
1318  "ESD_TruePrimaryMeson_ResolPt_MCPt; p_{T,MC}(GeV/c); (p_{T,rec}-p_{T,MC})/p_{T,MC}()",
1319  500, 0.03, 25, 1000, -1., 1.);
1320  fHistoTruePrimaryMesonMCPtResolPt[iCut]->Sumw2();
1322  fTrueList[iCut]->Add(fHistoTruePrimaryMesonMCPtResolPt[iCut]);
1323  }
1324  if(fDoMesonQA>1){
1325  fHistoTrueBckGGInvMassPt[iCut] = new TH2F("ESD_TrueBckGG_InvMass_Pt",
1326  "ESD_TrueBckGG_InvMass_Pt; M_{inv} (GeV/c^{2}) #gamma #gamma no signal; #pair p_{T}(GeV/c)",
1327  fMesonInvMassNBins, fMesonInvMassMin, fMesonInvMassMax, nBinsPt, arrPtBinning);
1328  fTrueList[iCut]->Add(fHistoTrueBckGGInvMassPt[iCut]);
1329  fHistoTrueBckFullMesonContainedInOneClusterInvMassPt[iCut] = new TH2F( "ESD_TrueBckFullMesonContained_InvMass_Pt",
1330  "ESD_TrueBckFullMesonContained_InvMass_Pt; M_{inv} (GeV/c^{2}) #gamma #gamma, calo gamma with full pi0; #pair p_{T}(GeV/c)",
1331  fMesonInvMassNBins, fMesonInvMassMin, fMesonInvMassMax, nBinsPt, arrPtBinning);
1333  fHistoTrueBckAsymEClustersInvMassPt[iCut] = new TH2F("ESD_TrueBckAsymEClus_InvMass_Pt",
1334  "ESD_TrueBckAsymEClus_InvMass_Pt; M_{inv} (GeV/c^{2}) #gamma #gamma, calo gamma >70% of pi0 energy; #pair p_{T}(GeV/c)",
1335  fMesonInvMassNBins, fMesonInvMassMin, fMesonInvMassMax, nBinsPt, arrPtBinning);
1337  fHistoTrueBckContInvMassPt[iCut] = new TH2F("ESD_TrueBckCont_InvMass_Pt",
1338  "ESD_TrueBckCont_InvMass_Pt; M_{inv} (GeV/c^{2}) contamination; #pair p_{T}(GeV/c)",
1339  fMesonInvMassNBins, fMesonInvMassMin, fMesonInvMassMax, nBinsPt, arrPtBinning);
1340  fTrueList[iCut]->Add(fHistoTrueBckContInvMassPt[iCut]);
1341  }
1342  fHistoTrueMesonPtY[iCut] = new TH2F("ESD_TrueMeson_Pt_Y", "ESD_TrueMeson_Pt_Y; p_{T}(GeV/c); Y",
1343  nBinsQAPt, arrQAPtBinning, 150, -1.5, 1.5);
1344  fTrueList[iCut]->Add(fHistoTrueMesonPtY[iCut]);
1345  fHistoTrueMesonPtAlpha[iCut] = new TH2F("ESD_TrueMeson_Pt_Alpha", "ESD_TrueMeson_Pt_Alpha; p_{T}(GeV/c); #alpha",
1346  nBinsQAPt, arrQAPtBinning, 200, -1, 1);
1347  fTrueList[iCut]->Add(fHistoTrueMesonPtAlpha[iCut]);
1348  fHistoTrueMesonPtOpenAngle[iCut] = new TH2F("ESD_TrueMeson_Pt_OpenAngle", "ESD_TrueMeson_Pt_OpenAngle; p_{T}(GeV/c); #theta",
1349  nBinsQAPt, arrQAPtBinning, 100, 0, 1);
1350  fTrueList[iCut]->Add(fHistoTrueMesonPtOpenAngle[iCut]);
1351 
1352  if (fIsMC > 1){
1353  if(fDoMesonQA>1){
1354  fHistoTrueBckGGInvMassPt[iCut]->Sumw2();
1356  fHistoTrueBckAsymEClustersInvMassPt[iCut]->Sumw2();
1357  fHistoTrueBckContInvMassPt[iCut]->Sumw2();
1358  }
1359  fHistoTrueMesonPtY[iCut]->Sumw2();
1360  fHistoTrueMesonPtAlpha[iCut]->Sumw2();
1361  fHistoTrueMesonPtOpenAngle[iCut]->Sumw2();
1363  }
1364  }
1365  }
1366  }
1367 
1371 
1375 
1376  fVectorRecTrueMesons.clear();
1377 
1378  if(fV0Reader)
1379  if((AliConversionPhotonCuts*)fV0Reader->GetConversionCuts())
1380  if(((AliConversionPhotonCuts*)fV0Reader->GetConversionCuts())->GetCutHistograms())
1381  fOutputContainer->Add(((AliConversionPhotonCuts*)fV0Reader->GetConversionCuts())->GetCutHistograms());
1382  if(fV0Reader)
1383  if((AliConvEventCuts*)fV0Reader->GetEventCuts())
1384  if(((AliConvEventCuts*)fV0Reader->GetEventCuts())->GetCutHistograms())
1385  fOutputContainer->Add(((AliConvEventCuts*)fV0Reader->GetEventCuts())->GetCutHistograms());
1386 
1387  if( fV0Reader && fV0Reader->GetProduceV0FindingEfficiency())
1388  if (fV0Reader->GetV0FindingEfficiencyHistograms())
1389  fOutputContainer->Add(fV0Reader->GetV0FindingEfficiencyHistograms());
1390 
1391  for(Int_t iMatcherTask = 0; iMatcherTask < 3; iMatcherTask++){
1392  AliCaloTrackMatcher* temp = (AliCaloTrackMatcher*) (AliAnalysisManager::GetAnalysisManager()->GetTask(Form("CaloTrackMatcher_%i",iMatcherTask)));
1393  if(temp) fOutputContainer->Add(temp->GetCaloTrackMatcherHistograms());
1394  }
1395 
1396  for(Int_t iCut = 0; iCut<fnCuts;iCut++){
1397  if(!((AliConvEventCuts*)fEventCutArray->At(iCut))) continue;
1398  if(((AliConvEventCuts*)fEventCutArray->At(iCut))->GetCutHistograms()){
1399  fCutFolder[iCut]->Add(((AliConvEventCuts*)fEventCutArray->At(iCut))->GetCutHistograms());
1400  }
1401  if (fMesonRecoMode < 2){
1402  if(!((AliConversionPhotonCuts*)fCutArray->At(iCut))) continue;
1403  if(((AliConversionPhotonCuts*)fCutArray->At(iCut))->GetCutHistograms()){
1404  fCutFolder[iCut]->Add(((AliConversionPhotonCuts*)fCutArray->At(iCut))->GetCutHistograms());
1405  }
1406  }
1408  if(!((AliCaloPhotonCuts*)fClusterCutArray->At(iCut))) continue;
1409  if(((AliCaloPhotonCuts*)fClusterCutArray->At(iCut))->GetCutHistograms()){
1410  fCutFolder[iCut]->Add(((AliCaloPhotonCuts*)fClusterCutArray->At(iCut))->GetCutHistograms());
1411  }
1412  if(fSetPlotHistsExtQA){
1413  if(!((AliCaloPhotonCuts*)fClusterCutArray->At(iCut))) continue;
1414  if(((AliCaloPhotonCuts*)fClusterCutArray->At(iCut))->GetExtQAHistograms()){
1415  fCutFolder[iCut]->Add(((AliCaloPhotonCuts*)fClusterCutArray->At(iCut))->GetExtQAHistograms());
1416  }
1417  }
1418  }
1419  if(!((AliConversionMesonCuts*)fMesonCutArray->At(iCut))) continue;
1420  if(((AliConversionMesonCuts*)fMesonCutArray->At(iCut))->GetCutHistograms()){
1421  fCutFolder[iCut]->Add(((AliConversionMesonCuts*)fMesonCutArray->At(iCut))->GetCutHistograms());
1422  }
1423  }
1424 
1425  if (fIsMC > 0){
1426  fTreeBrokenFiles = new TTree("BrokenFiles", "BrokenFiles");
1427  fTreeBrokenFiles->Branch("fileName",&fFileNameBroken);
1429  }
1430 
1431  PostData(1, fOutputContainer);
1432 }
1433 //_____________________________________________________________________________
1435 
1436  for(Int_t iCut = 0; iCut<fnCuts;iCut++){
1437  if (((AliConvEventCuts*)fEventCutArray->At(iCut))->GetPeriodEnum() == AliConvEventCuts::kNoPeriod && ((AliConvEventCuts*)fV0Reader->GetEventCuts())->GetPeriodEnum() != AliConvEventCuts::kNoPeriod){
1438  ((AliConvEventCuts*)fEventCutArray->At(iCut))->SetPeriodEnumExplicit(((AliConvEventCuts*)fV0Reader->GetEventCuts())->GetPeriodEnum());
1439  } else if (((AliConvEventCuts*)fEventCutArray->At(iCut))->GetPeriodEnum() == AliConvEventCuts::kNoPeriod ){
1440  ((AliConvEventCuts*)fEventCutArray->At(iCut))->SetPeriodEnum(fV0Reader->GetPeriodName());
1441  }
1442 
1443  if(fIsHeavyIon == 2){
1444  if(!((AliConvEventCuts*)fEventCutArray->At(iCut))->GetDoEtaShift()){
1445  fProfileEtaShift[iCut]->Fill(0.,(((AliConvEventCuts*)fEventCutArray->At(iCut))->GetEtaShift()));
1446  continue; // No Eta Shift requested, continue
1447  }
1448  if(((AliConvEventCuts*)fEventCutArray->At(iCut))->GetEtaShift() == 0.0){ // Eta Shift requested but not set, get shift automatically
1449  ((AliConvEventCuts*)fEventCutArray->At(iCut))->GetCorrectEtaShiftFromPeriod();
1450  fProfileEtaShift[iCut]->Fill(0.,(((AliConvEventCuts*)fEventCutArray->At(iCut))->GetEtaShift()));
1451  ((AliConvEventCuts*)fEventCutArray->At(iCut))->DoEtaShift(kFALSE); // Eta Shift Set, make sure that it is called only once
1452  continue;
1453  } else{
1454  printf(" Gamma Conversion Task %s :: Eta Shift Manually Set to %f \n\n",
1455  (((AliConvEventCuts*)fEventCutArray->At(iCut))->GetCutNumber()).Data(),((AliConvEventCuts*)fEventCutArray->At(iCut))->GetEtaShift());
1456  fProfileEtaShift[iCut]->Fill(0.,(((AliConvEventCuts*)fEventCutArray->At(iCut))->GetEtaShift()));
1457  ((AliConvEventCuts*)fEventCutArray->At(iCut))->DoEtaShift(kFALSE); // Eta Shift Set, make sure that it is called only once
1458  }
1459  }
1460  }
1461 
1462  return kTRUE;
1463 }
1464 //_____________________________________________________________________________
1466  //
1467  // Called for each event
1468  //
1469  fInputEvent = InputEvent();
1470 
1471  if(fIsMC > 0) fMCEvent = MCEvent();
1472 
1473  Int_t eventQuality = ((AliConvEventCuts*)fV0Reader->GetEventCuts())->GetEventQuality();
1474  if(fInputEvent->IsIncompleteDAQ()==kTRUE) eventQuality = 2;// incomplete event
1475  // Event Not Accepted due to MC event missing or wrong trigger for V0ReaderV1 or because it is incomplete => abort processing of this event/file
1476  if(eventQuality == 2 || eventQuality == 3){
1477  // write out name of broken file for first event
1478  if (fIsMC > 0){
1479  if (fInputEvent->IsA()==AliESDEvent::Class()){
1480  if (((AliESDEvent*)fInputEvent)->GetEventNumberInFile() == 0){
1481  fFileNameBroken = new TObjString(Form("%s",((TString)fV0Reader->GetCurrentFileName()).Data()));
1482  if (fTreeBrokenFiles) fTreeBrokenFiles->Fill();
1483  delete fFileNameBroken;
1484  }
1485  }
1486  }
1487 
1488  for(Int_t iCut = 0; iCut<fnCuts; iCut++){
1489  fHistoNEvents[iCut]->Fill(eventQuality);
1490  if (fIsMC > 1) fHistoNEventsWOWeight[iCut]->Fill(eventQuality);
1491  }
1492  return;
1493  }
1494 
1495 // if(fInputEvent->IsA()==AliAODEvent::Class()){
1496 // fInputEvent->InitMagneticField();
1497 // }
1498 
1499  fReaderGammas = fV0Reader->GetReconstructedGammas(); // Gammas from default Cut
1500 
1501  // ------------------- BeginEvent ----------------------------
1502  AliEventplane *EventPlane = fInputEvent->GetEventplane();
1503  if(fIsHeavyIon ==1)fEventPlaneAngle = EventPlane->GetEventplane("V0",fInputEvent,2);
1504  else fEventPlaneAngle=0.0;
1505 
1506  if(fIsMC>0 && fInputEvent->IsA()==AliAODEvent::Class() && !(fV0Reader->AreAODsRelabeled()) && fMesonRecoMode < 2){
1507  RelabelAODPhotonCandidates(kTRUE);// In case of AODMC relabeling MC
1508  fV0Reader->RelabelAODs(kTRUE);
1509  }
1510 
1511  for(Int_t iCut = 0; iCut<fnCuts; iCut++){
1512 
1513  fiCut = iCut;
1514  Bool_t isRunningEMCALrelAna = kFALSE;
1516  if (((AliCaloPhotonCuts*)fClusterCutArray->At(fiCut))->GetClusterType() == 1) isRunningEMCALrelAna = kTRUE;
1517  }
1518 
1519  Int_t eventNotAccepted = ((AliConvEventCuts*)fEventCutArray->At(iCut))->IsEventAcceptedByCut(fV0Reader->GetEventCuts(),fInputEvent,fMCEvent,fIsHeavyIon,isRunningEMCALrelAna);
1520 
1521  if(fIsMC==2){
1522  Float_t xsection = -1.;
1523  Float_t ntrials = -1.;
1524  ((AliConvEventCuts*)fEventCutArray->At(iCut))->GetXSectionAndNTrials(fMCEvent,xsection,ntrials,fInputEvent);
1525  if((xsection==-1.) || (ntrials==-1.)) AliFatal(Form("ERROR: GetXSectionAndNTrials returned invalid xsection/ntrials, periodName from V0Reader: '%s'",fV0Reader->GetPeriodName().Data()));
1526  fProfileJetJetXSection[iCut]->Fill(0.,xsection);
1527  fHistoJetJetNTrials[iCut]->Fill("#sum{NTrials}",ntrials);
1528  }
1529 
1530  if(fIsMC>0){
1531  fWeightJetJetMC = 1;
1532  // cout << fMCEvent << endl;
1533  Float_t pthard = -1;
1534  Bool_t isMCJet = ((AliConvEventCuts*)fEventCutArray->At(iCut))->IsJetJetMCEventAccepted( fMCEvent, fWeightJetJetMC, pthard, fInputEvent );
1535  if (fIsMC == 3){
1536  Double_t weightMult = ((AliConvEventCuts*)fEventCutArray->At(iCut))->GetWeightForMultiplicity(fV0Reader->GetNumberOfPrimaryTracks());
1537  fWeightJetJetMC = fWeightJetJetMC*weightMult;
1538  }
1539  if(fIsMC==1) fWeightJetJetMC = 1;
1540  if (!isMCJet){
1541  fHistoNEvents[iCut]->Fill(10,fWeightJetJetMC);
1542  if (fIsMC>1) fHistoNEventsWOWeight[iCut]->Fill(10);
1543  continue;
1544  }
1545  }
1546 
1547  Bool_t triggered = kTRUE;
1548 
1549  if(eventNotAccepted!= 0){
1550  // cout << "event rejected due to wrong trigger: " <<eventNotAccepted << endl;
1551  fHistoNEvents[iCut]->Fill(eventNotAccepted, fWeightJetJetMC); // Check Centrality, PileUp, SDD and V0AND --> Not Accepted => eventQuality = 1
1552  if (fIsMC>1) fHistoNEventsWOWeight[iCut]->Fill(eventNotAccepted);
1553  if (eventNotAccepted==3 && fIsMC > 0){
1554  triggered = kFALSE;
1555  }else {
1556  continue;
1557  }
1558  }
1559 
1560  if(eventQuality != 0){// Event Not Accepted
1561  //cout << "event rejected due to: " <<eventQuality << endl;
1562  fHistoNEvents[iCut]->Fill(eventQuality, fWeightJetJetMC);
1563  if (fIsMC>1) fHistoNEventsWOWeight[iCut]->Fill(eventQuality);
1564  continue;
1565  }
1566 
1567  if(fMesonRecoMode < 2){
1568  if(((AliConversionPhotonCuts*)fCutArray->At(fiCut))->GetDoElecDeDxPostCalibration()){
1569  if(!((AliConversionPhotonCuts*)fCutArray->At(fiCut))->LoadElecDeDxPostCalibration(fInputEvent->GetRunNumber())){
1570  AliFatal(Form("ERROR: LoadElecDeDxPostCalibration returned kFALSE for %d despite being requested!",fInputEvent->GetRunNumber()));
1571  }
1572  }
1573  }
1574 
1575  if (triggered==kTRUE){
1576  fHistoNEvents[iCut]->Fill(eventQuality, fWeightJetJetMC); // Should be 0 here
1577  if (fIsMC>1) fHistoNEventsWOWeight[iCut]->Fill(eventQuality); // Should be 0 here
1578 
1579  fHistoNGoodESDTracks[iCut]->Fill(fV0Reader->GetNumberOfPrimaryTracks(), fWeightJetJetMC);
1580  fHistoVertexZ[iCut]->Fill(fInputEvent->GetPrimaryVertex()->GetZ(), fWeightJetJetMC);
1581  if(!fDoLightOutput){
1582  fHistoVertexX[iCut]->Fill(fInputEvent->GetPrimaryVertex()->GetX(), fWeightJetJetMC);
1583  fHistoVertexY[iCut]->Fill(fInputEvent->GetPrimaryVertex()->GetY(), fWeightJetJetMC);
1584  fHistoSPDClusterTrackletBackground[iCut]->Fill(fInputEvent->GetMultiplicity()->GetNumberOfTracklets(),(fInputEvent->GetNumberOfITSClusters(0)+fInputEvent->GetNumberOfITSClusters(1)),fWeightJetJetMC);
1585  if(((AliConvEventCuts*)fEventCutArray->At(iCut))->IsHeavyIon() == 2) fHistoNV0Tracks[iCut]->Fill(fInputEvent->GetVZEROData()->GetMTotV0A(), fWeightJetJetMC);
1586  else fHistoNV0Tracks[iCut]->Fill(fInputEvent->GetVZEROData()->GetMTotV0A()+fInputEvent->GetVZEROData()->GetMTotV0C(), fWeightJetJetMC);
1587  }
1588  }
1589 
1590  if(fIsMC>0){
1591  // Process MC Particle
1592  if(((AliConvEventCuts*)fEventCutArray->At(iCut))->GetSignalRejection() != 0){
1593  if(fInputEvent->IsA()==AliESDEvent::Class()){
1594  ((AliConvEventCuts*)fEventCutArray->At(iCut))->GetNotRejectedParticles(((AliConvEventCuts*)fEventCutArray->At(iCut))->GetSignalRejection(),
1595  ((AliConvEventCuts*)fEventCutArray->At(iCut))->GetAcceptedHeader(),
1596  fMCEvent);
1597  }
1598  else if(fInputEvent->IsA()==AliAODEvent::Class()){
1599  ((AliConvEventCuts*)fEventCutArray->At(iCut))->GetNotRejectedParticles(((AliConvEventCuts*)fEventCutArray->At(iCut))->GetSignalRejection(),
1600  ((AliConvEventCuts*)fEventCutArray->At(iCut))->GetAcceptedHeader(),
1601  fInputEvent);
1602  }
1603 
1604  if(((AliConvEventCuts*)fEventCutArray->At(iCut))->GetAcceptedHeader()){
1605  for(Int_t i = 0;i<(((AliConvEventCuts*)fEventCutArray->At(iCut))->GetAcceptedHeader())->GetEntries();i++){
1606  TString nameBin= fHistoMCHeaders[iCut]->GetXaxis()->GetBinLabel(i+1);
1607  if (nameBin.CompareTo("")== 0){
1608  TString nameHeader = ((TObjString*)((TList*)((AliConvEventCuts*)fEventCutArray->At(iCut))
1609  ->GetAcceptedHeader())->At(i))->GetString();
1610  fHistoMCHeaders[iCut]->GetXaxis()->SetBinLabel(i+1,nameHeader.Data());
1611  }
1612  }
1613  }
1614  }
1615  }
1616  if(fIsMC>0){
1617  if(fInputEvent->IsA()==AliESDEvent::Class())
1619  if(fInputEvent->IsA()==AliAODEvent::Class())
1621  }
1622 
1623  if (triggered==kFALSE) continue;
1624 
1625  // it is in the loop to have the same conversion cut string (used also for MC stuff that should be same for V0 and Cluster)
1626  if (fMesonRecoMode > 0 || fEnableClusterCutsForTrigger) // process calo clusters
1627  ProcessClusters();
1628  if (fMesonRecoMode < 2) // Process this cuts gammas
1630 
1631  if (fMesonRecoMode < 2) fHistoNGammaConvCandidates[iCut]->Fill(fGammaCandidates->GetEntries(),fWeightJetJetMC);
1633 
1634  if (fMesonRecoMode < 2){
1635  if(((AliConversionMesonCuts*)fMesonCutArray->At(iCut))->UseMCPSmearing() && fIsMC>0){
1636  fUnsmearedPx = new Double_t[fGammaCandidates->GetEntries()]; // Store unsmeared Momenta
1637  fUnsmearedPy = new Double_t[fGammaCandidates->GetEntries()];
1638  fUnsmearedPz = new Double_t[fGammaCandidates->GetEntries()];
1639  fUnsmearedE = new Double_t[fGammaCandidates->GetEntries()];
1640 
1641  for(Int_t gamma=0;gamma<fGammaCandidates->GetEntries();gamma++){ // Smear the AODPhotons in MC
1642  fUnsmearedPx[gamma] = ((AliAODConversionPhoton*)fGammaCandidates->At(gamma))->Px();
1643  fUnsmearedPy[gamma] = ((AliAODConversionPhoton*)fGammaCandidates->At(gamma))->Py();
1644  fUnsmearedPz[gamma] = ((AliAODConversionPhoton*)fGammaCandidates->At(gamma))->Pz();
1645  fUnsmearedE[gamma] = ((AliAODConversionPhoton*)fGammaCandidates->At(gamma))->E();
1646  ((AliConversionMesonCuts*)fMesonCutArray->At(iCut))->SmearParticle(dynamic_cast<AliAODConversionPhoton*>(fGammaCandidates->At(gamma)));
1647  }
1648  }
1649  }
1650 
1651  // check gamma gamma pairs and veto if necessary
1652  if (!((AliConversionMesonCuts*)fMesonCutArray->At(iCut))->GetAcceptMassFlag())
1653  SetPhotonVeto();
1654 
1655  CalculateMesonCandidates(); // Combine Gammas from conversion and from calo
1656 
1657  if(((AliConversionMesonCuts*)fMesonCutArray->At(iCut))->DoBGCalculation()){
1658  if(((AliConversionMesonCuts*)fMesonCutArray->At(iCut))->BackgroundHandlerType() == 0){
1659  CalculateBackground(); // Combinatorial Background
1660  UpdateEventByEventData(); // Store Event for mixed Events
1661  } else{
1662  CalculateBackgroundRP(); // Combinatorial Background
1663  fBGHandlerRP[iCut]->AddEvent(fGammaCandidates,fInputEvent); // Store Event for mixed Events
1664  fBGClusHandlerRP[iCut]->AddEvent(fClusterCandidates,fInputEvent); // Store Event for mixed Events
1665  }
1666  }
1667 
1668  if (fMesonRecoMode < 2){
1669  if(((AliConversionMesonCuts*)fMesonCutArray->At(iCut))->UseMCPSmearing() && fIsMC>0){
1670  for(Int_t gamma=0;gamma<fGammaCandidates->GetEntries();gamma++){ // Smear the AODPhotons in MC
1671  ((AliAODConversionPhoton*)fGammaCandidates->At(gamma))->SetPx(fUnsmearedPx[gamma]); // Reset Unsmeared Momenta
1672  ((AliAODConversionPhoton*)fGammaCandidates->At(gamma))->SetPy(fUnsmearedPy[gamma]);
1673  ((AliAODConversionPhoton*)fGammaCandidates->At(gamma))->SetPz(fUnsmearedPz[gamma]);
1674  ((AliAODConversionPhoton*)fGammaCandidates->At(gamma))->SetE(fUnsmearedE[gamma]);
1675  }
1676  delete[] fUnsmearedPx; fUnsmearedPx = 0x0;
1677  delete[] fUnsmearedPy; fUnsmearedPy = 0x0;
1678  delete[] fUnsmearedPz; fUnsmearedPz = 0x0;
1679  delete[] fUnsmearedE;fUnsmearedE = 0x0;
1680  }
1681  }
1682 
1683  if(fIsMC>0){
1684  fVectorRecTrueMesons.clear();
1691  }
1692 
1693  fGammaCandidates->Clear(); // delete this cuts good gammas
1694  fClusterCandidates->Clear(); // delete cluster candidates
1695  }
1696 
1697  if(fIsMC>0 && fInputEvent->IsA()==AliAODEvent::Class() && !(fV0Reader->AreAODsRelabeled()) && fMesonRecoMode < 2){
1698  RelabelAODPhotonCandidates(kFALSE); // Back to ESDMC Label
1699  fV0Reader->RelabelAODs(kFALSE);
1700  }
1701 
1702  PostData(1, fOutputContainer);
1703 }
1704 
1705 //________________________________________________________________________
1707  Int_t nclus = 0;
1708  TClonesArray * arrClustersProcess = NULL;
1709  if(!fCorrTaskSetting.CompareTo("")){
1710  nclus = fInputEvent->GetNumberOfCaloClusters();
1711  } else {
1712  arrClustersProcess = dynamic_cast<TClonesArray*>(fInputEvent->FindListObject(Form("%sClustersBranch",fCorrTaskSetting.Data())));
1713  if(!arrClustersProcess)
1714  AliFatal(Form("%sClustersBranch was not found in AliAnalysisTaskHeavyNeutralMesonToGG! Check the correction framework settings!",fCorrTaskSetting.Data()));
1715  nclus = arrClustersProcess->GetEntries();
1716  }
1717 
1718  vector<AliAODConversionPhoton*> vectorCurrentClusters;
1719  vector<Int_t> vectorRejectCluster;
1720  vector<Double_t> vectorPhotonWeight;
1721  vector<Double_t> vectorClusterM02;
1722 
1723  // cout << nclus << endl;
1724 
1725  if(nclus == 0) return;
1726 
1727  // plotting histograms on cell/tower level, only if extendedMatchAndQA > 1
1728  ((AliCaloPhotonCuts*)fClusterCutArray->At(fiCut))->FillHistogramsExtendedQA(fInputEvent,fIsMC);
1729 
1730  // match tracks to clusters
1731  if(fDoPrimaryTrackMatching) ((AliCaloPhotonCuts*)fClusterCutArray->At(fiCut))->MatchTracksToClusters(fInputEvent,fWeightJetJetMC,kFALSE);
1732 
1733  // vertex
1734  Double_t vertex[3] = {0,0,0};
1735  InputEvent()->GetPrimaryVertex()->GetXYZ(vertex);
1736 
1737  // Loop over EMCal clusters
1738  for(Int_t i = 0; i < nclus; i++){
1739  AliVCluster* clus = NULL;
1740  if(fInputEvent->IsA()==AliESDEvent::Class()){
1741  if(arrClustersProcess)
1742  clus = new AliESDCaloCluster(*(AliESDCaloCluster*)arrClustersProcess->At(i));
1743  else
1744  clus = new AliESDCaloCluster(*(AliESDCaloCluster*)fInputEvent->GetCaloCluster(i));
1745  } else if(fInputEvent->IsA()==AliAODEvent::Class()){
1746  if(arrClustersProcess)
1747  clus = new AliAODCaloCluster(*(AliAODCaloCluster*)arrClustersProcess->At(i));
1748  else
1749  clus = new AliAODCaloCluster(*(AliAODCaloCluster*)fInputEvent->GetCaloCluster(i));
1750  }
1751 
1752  if (!clus) continue;
1753  if(!((AliCaloPhotonCuts*)fClusterCutArray->At(fiCut))->ClusterIsSelected(clus,fInputEvent,fMCEvent,fIsMC,fWeightJetJetMC,i)){
1754  delete clus;
1755  continue;
1756  }
1757 
1758  // TLorentzvector with cluster
1759  TLorentzVector clusterVector;
1760  clus->GetMomentum(clusterVector,vertex);
1761 
1762  TLorentzVector* tmpvec = new TLorentzVector();
1763  tmpvec->SetPxPyPzE(clusterVector.Px(),clusterVector.Py(),clusterVector.Pz(),clusterVector.E());
1764 
1765  // convert to AODConversionPhoton
1766  AliAODConversionPhoton *PhotonCandidate = new AliAODConversionPhoton(tmpvec);
1767  if(!PhotonCandidate){ delete clus; delete tmpvec; continue;}
1768 
1769  // Flag Photon as CaloPhoton
1770  PhotonCandidate->SetIsCaloPhoton(((AliCaloPhotonCuts*)fClusterCutArray->At(fiCut))->GetClusterType());
1771  PhotonCandidate->SetCaloClusterRef((Long_t)i);
1772  PhotonCandidate->SetLeadingCellID(((AliCaloPhotonCuts*)fClusterCutArray->At(fiCut))->FindLargestCellInCluster(clus,fInputEvent));
1773 
1774  // get MC label
1775  if(fIsMC>0){
1776  PhotonCandidate->SetNCaloPhotonMCLabels(clus->GetNLabels());
1777  Int_t* mclabelsCluster = clus->GetLabels();
1778  PhotonCandidate->SetNCaloPhotonMCLabels(clus->GetNLabels());
1779  // cout << clus->GetNLabels() << endl;
1780  if (clus->GetNLabels()>0){
1781  for (Int_t k =0; k<(Int_t)clus->GetNLabels(); k++){
1782  if (k<50)PhotonCandidate->SetCaloPhotonMCLabel(k,mclabelsCluster[k]);
1783  // Int_t pdgCode = fMCEvent->Particle(mclabelsCluster[k])->GetPdgCode();
1784  // cout << "label " << k << "\t" << mclabelsCluster[k] << " pdg code: " << pdgCode << endl;
1785  }
1786  }
1787  }
1788 
1789  fIsFromDesiredHeader = kTRUE;
1791  //TString periodName = fV0Reader->GetPeriodName();
1792  // test whether largest contribution to cluster orginates in added signals
1793  if (fIsMC>0 && ((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetSignalRejection() > 0){
1794  if (((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsParticleFromBGEvent(PhotonCandidate->GetCaloPhotonMCLabel(0), fMCEvent, fInputEvent) == 0){
1795  fIsFromDesiredHeader = kFALSE;
1796  }
1797  if (clus->GetNLabels()>1){
1798  Int_t* mclabelsCluster = clus->GetLabels();
1799  for (Int_t l = 1; l < (Int_t)clus->GetNLabels(); l++ ){
1800  if (((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsParticleFromBGEvent(mclabelsCluster[l], fMCEvent, fInputEvent) == 0){
1802  }
1803  }
1804  }
1805  }
1807  fHistoClusAllHeadersGammaPt[fiCut]->Fill(PhotonCandidate->Pt(),fWeightJetJetMC);
1809  fHistoClusRejectedHeadersGammaPt[fiCut]->Fill(PhotonCandidate->Pt(),fWeightJetJetMC);
1811  fHistoClusOverlapHeadersGammaPt[fiCut]->Fill(PhotonCandidate->Pt(),fWeightJetJetMC);
1812 
1814 
1815  vectorCurrentClusters.push_back(PhotonCandidate);
1816  vectorPhotonWeight.push_back(fWeightJetJetMC);
1817  vectorClusterM02.push_back(clus->GetM02());
1818  }else{
1819  delete PhotonCandidate;
1820  }
1821 
1822  delete clus;
1823  delete tmpvec;
1824  }
1825 
1826  Bool_t rejected = kFALSE;
1827  // run conversion recovery in addition
1828  if (((AliCaloPhotonCuts*)fClusterCutArray->At(fiCut))->GetIsConversionRecovery()){
1829  rejected = ((AliCaloPhotonCuts*)fClusterCutArray->At(fiCut))->CheckForReconstructedConversionPairs(vectorCurrentClusters,vectorRejectCluster);
1830  }
1831 
1832  for (Int_t iter = 0; iter < (Int_t)vectorCurrentClusters.size();iter++){
1833  if (!((AliCaloPhotonCuts*)fClusterCutArray->At(fiCut))->CheckVectorForIndexAndAdd(vectorRejectCluster, iter,kFALSE)){
1834  fHistoClusGammaPt[fiCut]->Fill(vectorCurrentClusters.at(iter)->Pt(), vectorPhotonWeight.at(iter));
1835  fHistoClusGammaE[fiCut]->Fill(vectorCurrentClusters.at(iter)->E(), vectorPhotonWeight.at(iter));
1836  if(fIsMC> 0){
1837  if(fInputEvent->IsA()==AliESDEvent::Class()){
1838  ProcessTrueClusterCandidates(vectorCurrentClusters.at(iter),vectorClusterM02.at(iter));
1839  } else {
1840  ProcessTrueClusterCandidatesAOD(vectorCurrentClusters.at(iter),vectorClusterM02.at(iter));
1841  }
1842  }
1843  fClusterCandidates->Add(vectorCurrentClusters.at(iter));
1844  }
1845  }
1846  vectorRejectCluster.clear();
1847  vectorPhotonWeight.clear();
1848  vectorClusterM02.clear();
1849 }
1850 
1851 //________________________________________________________________________
1852 void AliAnalysisTaskHeavyNeutralMesonToGG::ProcessTrueClusterCandidates(AliAODConversionPhoton *TruePhotonCandidate, Float_t clusM02){
1853 
1854  TParticle *Photon = NULL;
1855  if (TruePhotonCandidate->GetIsCaloPhoton() == 0) AliFatal("CaloPhotonFlag has not been set task will abort");
1856  if (TruePhotonCandidate->GetCaloPhotonMCLabel(0) < 0) return;
1857  if(!fDoLightOutput) fHistoTrueNLabelsInClusPt[fiCut]->Fill(TruePhotonCandidate->GetNCaloPhotonMCLabels(),TruePhotonCandidate->Pt(),fWeightJetJetMC);
1858 
1859  if (TruePhotonCandidate->GetNCaloPhotonMCLabels()>0)Photon = fMCEvent->Particle(TruePhotonCandidate->GetCaloPhotonMCLabel(0));
1860  else return;
1861 
1862  if(Photon == NULL){
1863  // cout << "no photon" << endl;
1864  return;
1865  }
1866 
1867  TruePhotonCandidate->SetCaloPhotonMCFlags(fMCEvent, fEnableSortForClusMC);
1868 
1869  const AliVVertex* primVtxMC = fMCEvent->GetPrimaryVertex();
1870  Double_t mcProdVtxX = primVtxMC->GetX();
1871  Double_t mcProdVtxY = primVtxMC->GetY();
1872  Double_t mcProdVtxZ = primVtxMC->GetZ();
1873  Bool_t isPrimary = ((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsConversionPrimaryESD( fMCEvent, TruePhotonCandidate->GetCaloPhotonMCLabel(0), mcProdVtxX, mcProdVtxY, mcProdVtxZ);
1874 
1875  // to get primary distrinction right put mother of conversion electron as particle to check
1876  if (TruePhotonCandidate->IsLargestComponentElectron() && TruePhotonCandidate->IsConversion())
1877  isPrimary = ((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsConversionPrimaryESD( fMCEvent, Photon->GetMother(0), mcProdVtxX, mcProdVtxY, mcProdVtxZ);
1878 
1879 
1880  // Fill histograms for inclusive gamma corrections
1881  // --> a) all clusters with leading real or converted photons
1882  if (TruePhotonCandidate->IsLargestComponentPhoton() || (TruePhotonCandidate->IsLargestComponentElectron() && TruePhotonCandidate->IsConversion()) ){
1883  fHistoTrueClusGammaPt[fiCut]->Fill(TruePhotonCandidate->Pt(),fWeightJetJetMC);
1884  if(!fDoLightOutput){
1885  // how many of those clusters are from converted photons
1886  if (TruePhotonCandidate->IsLargestComponentElectron() && TruePhotonCandidate->IsConversion()){
1887  fHistoTrueClusConvGammaPt[fiCut]->Fill(TruePhotonCandidate->Pt(), fWeightJetJetMC);
1888  }
1889  // --> b) which of these are primary
1890  if(isPrimary){
1891  fHistoTruePrimaryClusGammaPt[fiCut]->Fill(TruePhotonCandidate->Pt(),fWeightJetJetMC);
1892  // --> c) which are from conversions? Just additonal information
1893  if (TruePhotonCandidate->IsLargestComponentElectron() && TruePhotonCandidate->IsConversion() && (Photon->GetMother(0)>-1)){
1894  fHistoTruePrimaryClusConvGammaPt[fiCut]->Fill(TruePhotonCandidate->Pt(), fWeightJetJetMC);
1895  }
1896  // --> d) how do the secondaries look like
1897  }
1898  }
1899  }
1900 
1901  // Some additional QA
1902  if (fDoClusterQA > 0){
1903  // how many of the converted photons are fully contained in the cluster
1904  if (TruePhotonCandidate->IsLargestComponentElectron() && TruePhotonCandidate->IsConversion() && TruePhotonCandidate->IsConversionFullyContained())
1905  fHistoTrueClusConvGammaFullyPt[fiCut]->Fill(TruePhotonCandidate->Pt(), fWeightJetJetMC);
1906  }
1907 
1908  // Check if we are double counting photons
1909  Int_t motherLab = Photon->GetMother(0);
1910  if (motherLab > -1){
1911  if (TruePhotonCandidate->IsLargestComponentPhoton()){
1913  fHistoDoubleCountTrueClusterGammaPt[fiCut]->Fill(TruePhotonCandidate->Pt(),(Double_t)0,fWeightJetJetMC);
1915  }
1916  }
1917  Int_t grandMotherLab = fMCEvent->Particle(motherLab)->GetMother(0);
1918  if (grandMotherLab > -1){
1919  if (TruePhotonCandidate->IsLargestComponentElectron() && TruePhotonCandidate->IsConversion()){
1921  fHistoDoubleCountTrueClusterGammaPt[fiCut]->Fill(TruePhotonCandidate->Pt(),(Double_t)1,fWeightJetJetMC);
1923  }
1924  }
1925  }
1926  }
1927  return;
1928 }
1929 
1930 
1931 //________________________________________________________________________
1932 void AliAnalysisTaskHeavyNeutralMesonToGG::ProcessTrueClusterCandidatesAOD(AliAODConversionPhoton *TruePhotonCandidate, Float_t clusM02){
1933 
1934  AliAODMCParticle *Photon = NULL;
1935  TClonesArray *AODMCTrackArray = dynamic_cast<TClonesArray*>(fInputEvent->FindListObject(AliAODMCParticle::StdBranchName()));
1936  if (AODMCTrackArray){
1937  if (TruePhotonCandidate->GetIsCaloPhoton() == 0) AliFatal("CaloPhotonFlag has not been set task will abort");
1938  if (TruePhotonCandidate->GetNCaloPhotonMCLabels()>0) Photon = (AliAODMCParticle*) AODMCTrackArray->At(TruePhotonCandidate->GetCaloPhotonMCLabel(0));
1939  else return;
1940  }else {
1941  AliInfo("AODMCTrackArray could not be loaded");
1942  return;
1943  }
1944 
1945  if(Photon == NULL){
1946  // cout << "no photon" << endl;
1947  return;
1948  }
1949 
1950  TruePhotonCandidate->SetCaloPhotonMCFlagsAOD(fInputEvent, fEnableSortForClusMC);
1951  if(!fDoLightOutput) fHistoTrueNLabelsInClusPt[fiCut]->Fill(TruePhotonCandidate->GetNCaloPhotonMCLabels(),TruePhotonCandidate->Pt(),fWeightJetJetMC);
1952 
1953  const AliVVertex* primVtxMC = fMCEvent->GetPrimaryVertex();
1954  Double_t mcProdVtxX = primVtxMC->GetX();
1955  Double_t mcProdVtxY = primVtxMC->GetY();
1956  Double_t mcProdVtxZ = primVtxMC->GetZ();
1957  Bool_t isPrimary = ((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsConversionPrimaryAOD(fInputEvent, Photon, mcProdVtxX, mcProdVtxY, mcProdVtxZ);
1958 
1959  // to get primary distrinction right put mother of conversion electron as particle to check
1960  if (TruePhotonCandidate->IsLargestComponentElectron() && TruePhotonCandidate->IsConversion()){
1961  if (Photon->GetMother()> -1){
1962  AliAODMCParticle *Mother = (AliAODMCParticle*) AODMCTrackArray->At(Photon->GetMother());
1963  isPrimary = ((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsConversionPrimaryAOD( fInputEvent, Mother, mcProdVtxX, mcProdVtxY, mcProdVtxZ);
1964  }
1965  }
1966 
1967  // True Photon
1968  if (TruePhotonCandidate->IsLargestComponentPhoton() || (TruePhotonCandidate->IsLargestComponentElectron() && TruePhotonCandidate->IsConversion()) ){
1969  fHistoTrueClusGammaPt[fiCut]->Fill(TruePhotonCandidate->Pt(),fWeightJetJetMC);
1970  if(!fDoLightOutput){
1971  if (TruePhotonCandidate->IsLargestComponentElectron() && TruePhotonCandidate->IsConversion()){
1972  fHistoTrueClusConvGammaPt[fiCut]->Fill(TruePhotonCandidate->Pt(), fWeightJetJetMC);
1973  }
1974  if(isPrimary){
1975  fHistoTruePrimaryClusGammaPt[fiCut]->Fill(TruePhotonCandidate->Pt(),fWeightJetJetMC);
1976  if (TruePhotonCandidate->IsLargestComponentElectron() && TruePhotonCandidate->IsConversion()){
1977  fHistoTruePrimaryClusConvGammaPt[fiCut]->Fill(TruePhotonCandidate->Pt(), fWeightJetJetMC);
1978  }
1979  }
1980  }
1981  }
1982  if (fDoClusterQA > 0){
1983  if (TruePhotonCandidate->IsLargestComponentElectron() && TruePhotonCandidate->IsConversion() && TruePhotonCandidate->IsConversionFullyContained())
1984  fHistoTrueClusConvGammaFullyPt[fiCut]->Fill(TruePhotonCandidate->Pt(), fWeightJetJetMC);
1985  }
1986  Int_t motherLab = Photon->GetMother();
1987  if (motherLab > -1){
1988  if (TruePhotonCandidate->IsLargestComponentPhoton()){
1990  fHistoDoubleCountTrueClusterGammaPt[fiCut]->Fill(TruePhotonCandidate->Pt(),(Double_t)0,fWeightJetJetMC);
1992  }
1993  }
1994  Int_t grandMotherLab = ((AliAODMCParticle*) AODMCTrackArray->At(motherLab))->GetMother();
1995  if (grandMotherLab > -1){
1996  if (TruePhotonCandidate->IsLargestComponentElectron() && TruePhotonCandidate->IsConversion()){
1998  fHistoDoubleCountTrueClusterGammaPt[fiCut]->Fill(TruePhotonCandidate->Pt(),(Double_t)1,fWeightJetJetMC);
2000  }
2001  }
2002  }
2003  }
2004 
2005  return;
2006 }
2007 
2008 //________________________________________________________________________
2010 
2011  Int_t nV0 = 0;
2012  TList *GammaCandidatesStepOne = new TList();
2013  TList *GammaCandidatesStepTwo = new TList();
2014  // Loop over Photon Candidates allocated by ReaderV1
2015  for(Int_t i = 0; i < fReaderGammas->GetEntriesFast(); i++){
2016  AliAODConversionPhoton* PhotonCandidate = (AliAODConversionPhoton*) fReaderGammas->At(i);
2017  if(!PhotonCandidate) continue;
2018  fIsFromDesiredHeader = kTRUE;
2019  if(fIsMC>0 && ((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetSignalRejection() != 0){
2020  Int_t isPosFromMBHeader = ((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsParticleFromBGEvent(PhotonCandidate->GetMCLabelPositive(), fMCEvent, fInputEvent);
2021  if(isPosFromMBHeader == 0 && ((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetSignalRejection() != 3) continue;
2022  Int_t isNegFromMBHeader = ((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsParticleFromBGEvent(PhotonCandidate->GetMCLabelNegative(), fMCEvent, fInputEvent);
2023  if(isNegFromMBHeader == 0 && ((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetSignalRejection() != 3) continue;
2024  if( (isNegFromMBHeader+isPosFromMBHeader) != 4) fIsFromDesiredHeader = kFALSE;
2025  }
2026 
2027  if(!((AliConversionPhotonCuts*)fCutArray->At(fiCut))->PhotonIsSelected(PhotonCandidate,fInputEvent)) continue;
2028  if(!((AliConversionPhotonCuts*)fCutArray->At(fiCut))->InPlaneOutOfPlaneCut(PhotonCandidate->GetPhotonPhi(),fEventPlaneAngle)) continue;
2029  if(!((AliConversionPhotonCuts*)fCutArray->At(fiCut))->UseElecSharingCut() &&
2030  !((AliConversionPhotonCuts*)fCutArray->At(fiCut))->UseToCloseV0sCut()){
2031  fGammaCandidates->Add(PhotonCandidate); // if no second loop is required add to events good gammas
2032 
2034  if(!fDoLightOutput) fHistoConvGammaPt[fiCut]->Fill(PhotonCandidate->Pt(),fWeightJetJetMC);
2035  }
2036  if(fIsMC>0){
2037  if(fInputEvent->IsA()==AliESDEvent::Class())
2038  ProcessTruePhotonCandidates(PhotonCandidate);
2039  if(fInputEvent->IsA()==AliAODEvent::Class())
2040  ProcessTruePhotonCandidatesAOD(PhotonCandidate);
2041  }
2042  }else if(((AliConversionPhotonCuts*)fCutArray->At(fiCut))->UseElecSharingCut()){ // if Shared Electron cut is enabled, Fill array, add to step one
2043  ((AliConversionPhotonCuts*)fCutArray->At(fiCut))->FillElectonLabelArray(PhotonCandidate,nV0);
2044  nV0++;
2045  GammaCandidatesStepOne->Add(PhotonCandidate);
2046  }else if(!((AliConversionPhotonCuts*)fCutArray->At(fiCut))->UseElecSharingCut() &&
2047  ((AliConversionPhotonCuts*)fCutArray->At(fiCut))->UseToCloseV0sCut()){ // shared electron is disabled, step one not needed -> step two
2048  GammaCandidatesStepTwo->Add(PhotonCandidate);
2049  }
2050  }
2051 
2052  if(((AliConversionPhotonCuts*)fCutArray->At(fiCut))->UseElecSharingCut()){
2053  for(Int_t i = 0;i<GammaCandidatesStepOne->GetEntries();i++){
2054  AliAODConversionPhoton *PhotonCandidate= (AliAODConversionPhoton*) GammaCandidatesStepOne->At(i);
2055  if(!PhotonCandidate) continue;
2056  fIsFromDesiredHeader = kTRUE;
2057  if(fMCEvent && ((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetSignalRejection() != 0){
2058  Int_t isPosFromMBHeader = ((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsParticleFromBGEvent(PhotonCandidate->GetMCLabelPositive(), fMCEvent, fInputEvent);
2059  Int_t isNegFromMBHeader = ((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsParticleFromBGEvent(PhotonCandidate->GetMCLabelNegative(), fMCEvent, fInputEvent);
2060  if( (isNegFromMBHeader+isPosFromMBHeader) != 4) fIsFromDesiredHeader = kFALSE;
2061  }
2062  if(!((AliConversionPhotonCuts*)fCutArray->At(fiCut))->RejectSharedElectronV0s(PhotonCandidate,i,GammaCandidatesStepOne->GetEntries())) continue;
2063  if(!((AliConversionPhotonCuts*)fCutArray->At(fiCut))->UseToCloseV0sCut()){ // To Colse v0s cut diabled, step two not needed
2064  fGammaCandidates->Add(PhotonCandidate);
2066  if(!fDoLightOutput) fHistoConvGammaPt[fiCut]->Fill(PhotonCandidate->Pt(),fWeightJetJetMC);
2067  }
2068  if(fIsMC>0){
2069  if(fInputEvent->IsA()==AliESDEvent::Class())
2070  ProcessTruePhotonCandidates(PhotonCandidate);
2071  if(fInputEvent->IsA()==AliAODEvent::Class())
2072  ProcessTruePhotonCandidatesAOD(PhotonCandidate);
2073  }
2074  } else GammaCandidatesStepTwo->Add(PhotonCandidate); // Close v0s cut enabled -> add to list two
2075  }
2076  }
2077 
2078  if(((AliConversionPhotonCuts*)fCutArray->At(fiCut))->UseToCloseV0sCut()){
2079  for(Int_t i = 0;i<GammaCandidatesStepTwo->GetEntries();i++){
2080  AliAODConversionPhoton* PhotonCandidate = (AliAODConversionPhoton*) GammaCandidatesStepTwo->At(i);
2081  if(!PhotonCandidate) continue;
2082  fIsFromDesiredHeader = kTRUE;
2083  if(fMCEvent && ((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetSignalRejection() != 0){
2084  Int_t isPosFromMBHeader = ((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsParticleFromBGEvent(PhotonCandidate->GetMCLabelPositive(), fMCEvent, fInputEvent);
2085  Int_t isNegFromMBHeader = ((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsParticleFromBGEvent(PhotonCandidate->GetMCLabelNegative(), fMCEvent, fInputEvent);
2086  if( (isNegFromMBHeader+isPosFromMBHeader) != 4) fIsFromDesiredHeader = kFALSE;
2087  }
2088  if(!((AliConversionPhotonCuts*)fCutArray->At(fiCut))->RejectToCloseV0s(PhotonCandidate,GammaCandidatesStepTwo,i)) continue;
2089  fGammaCandidates->Add(PhotonCandidate); // Add gamma to current cut TList
2091  if(!fDoLightOutput) fHistoConvGammaPt[fiCut]->Fill(PhotonCandidate->Pt(),fWeightJetJetMC);
2092  }
2093  if(fIsMC>0){
2094  if(fInputEvent->IsA()==AliESDEvent::Class())
2095  ProcessTruePhotonCandidates(PhotonCandidate);
2096  if(fInputEvent->IsA()==AliAODEvent::Class())
2097  ProcessTruePhotonCandidatesAOD(PhotonCandidate);
2098  }
2099  }
2100  }
2101 
2102  delete GammaCandidatesStepOne;
2103  GammaCandidatesStepOne = 0x0;
2104  delete GammaCandidatesStepTwo;
2105  GammaCandidatesStepTwo = 0x0;
2106 
2107 }
2108 
2109 //________________________________________________________________________
2110 void AliAnalysisTaskHeavyNeutralMesonToGG::ProcessTruePhotonCandidatesAOD(AliAODConversionPhoton *TruePhotonCandidate){
2111 
2112  const AliVVertex* primVtxMC = fMCEvent->GetPrimaryVertex();
2113  Double_t mcProdVtxX = primVtxMC->GetX();
2114  Double_t mcProdVtxY = primVtxMC->GetY();
2115  Double_t mcProdVtxZ = primVtxMC->GetZ();
2116 
2117  TClonesArray *AODMCTrackArray = dynamic_cast<TClonesArray*>(fInputEvent->FindListObject(AliAODMCParticle::StdBranchName()));
2118  if (AODMCTrackArray == NULL) return;
2119  AliAODMCParticle *posDaughter = (AliAODMCParticle*) AODMCTrackArray->At(TruePhotonCandidate->GetMCLabelPositive());
2120  AliAODMCParticle *negDaughter = (AliAODMCParticle*) AODMCTrackArray->At(TruePhotonCandidate->GetMCLabelNegative());
2121 
2122  if(posDaughter == NULL || negDaughter == NULL) return; // One particle does not exist
2123  Int_t pdgCode[2] = {TMath::Abs(posDaughter->GetPdgCode()),TMath::Abs(negDaughter->GetPdgCode())};
2124 
2125  if(posDaughter->GetMother() != negDaughter->GetMother()){
2126  return;
2127  } else if(posDaughter->GetMother() == -1){
2128  return;
2129  }
2130 
2131  if(pdgCode[0]!=11 || pdgCode[1]!=11){
2132  return; //One Particle is not a electron
2133  }
2134 
2135  if(posDaughter->GetPdgCode()==negDaughter->GetPdgCode()){
2136  return; // Same Charge
2137  }
2138 
2139  AliAODMCParticle *Photon = (AliAODMCParticle*) AODMCTrackArray->At(posDaughter->GetMother());
2140  if(Photon->GetPdgCode() != 22){
2141  return; // Mother is no Photon
2142  }
2143 
2144  if(((posDaughter->GetMCProcessCode())) != 5 || ((negDaughter->GetMCProcessCode())) != 5){
2145  return;// check if the daughters come from a conversion
2146  }
2147  // STILL A BUG IN ALIROOT >>8 HAS TPO BE REMOVED AFTER FIX
2148 
2149  // True Photon
2151  if(!fDoLightOutput) fHistoTrueConvGammaPt[fiCut]->Fill(TruePhotonCandidate->Pt(),fWeightJetJetMC);
2152  if (CheckVectorForDoubleCount(fVectorDoubleCountTrueConvGammas,posDaughter->GetMother())){
2153  if(!fDoLightOutput) fHistoDoubleCountTrueConvGammaRPt[fiCut]->Fill(TruePhotonCandidate->GetConversionRadius(),TruePhotonCandidate->Pt(),fWeightJetJetMC);
2154  FillMultipleCountMap(fMapMultipleCountTrueConvGammas,posDaughter->GetMother());
2155  }
2156  }
2157 
2158  Bool_t isPrimary = ((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsConversionPrimaryAOD(fInputEvent, Photon, mcProdVtxX, mcProdVtxY, mcProdVtxZ);
2159  if(isPrimary){
2160  // Count just primary MC Gammas as true --> For Ratio esdtruegamma / mcconvgamma
2162  if(!fDoLightOutput){
2163  fHistoTruePrimaryConvGammaPt[fiCut]->Fill(TruePhotonCandidate->Pt(),fWeightJetJetMC);
2164  }
2165  }
2166  // (Not Filled for i6, Extra Signal Gamma (parambox) are secondary)
2167  }
2168  TruePhotonCandidate->SetIsTrueConvertedPhoton();
2169  return;
2170 }
2171 
2172 //________________________________________________________________________
2173 void AliAnalysisTaskHeavyNeutralMesonToGG::ProcessTruePhotonCandidates(AliAODConversionPhoton *TruePhotonCandidate){
2174 
2175  const AliVVertex* primVtxMC = fMCEvent->GetPrimaryVertex();
2176  Double_t mcProdVtxX = primVtxMC->GetX();
2177  Double_t mcProdVtxY = primVtxMC->GetY();
2178  Double_t mcProdVtxZ = primVtxMC->GetZ();
2179 
2180  // Process True Photons
2181  TParticle *posDaughter = TruePhotonCandidate->GetPositiveMCDaughter(fMCEvent);
2182  TParticle *negDaughter = TruePhotonCandidate->GetNegativeMCDaughter(fMCEvent);
2183 
2184  if(posDaughter == NULL || negDaughter == NULL) return; // One particle does not exist
2185  Int_t pdgCode[2] = {TMath::Abs(posDaughter->GetPdgCode()),TMath::Abs(negDaughter->GetPdgCode())};
2186  if(posDaughter->GetMother(0) != negDaughter->GetMother(0)){
2187  return;
2188  }
2189  else if(posDaughter->GetMother(0) == -1){
2190  return;
2191  }
2192 
2193  if(pdgCode[0]!=11 || pdgCode[1]!=11) return; //One Particle is not a electron
2194 
2195  if(posDaughter->GetPdgCode()==negDaughter->GetPdgCode()) return; // Same Charge
2196 
2197  TParticle *Photon = TruePhotonCandidate->GetMCParticle(fMCEvent);
2198 
2199  if(Photon->GetPdgCode() != 22){
2200  return; // Mother is no Photon
2201  }
2202 
2203  if(posDaughter->GetUniqueID() != 5 || negDaughter->GetUniqueID() !=5) return;// check if the daughters come from a conversion
2204 
2205  // True Photon
2207  if(!fDoLightOutput) fHistoTrueConvGammaPt[fiCut]->Fill(TruePhotonCandidate->Pt(),fWeightJetJetMC);
2208  if (CheckVectorForDoubleCount(fVectorDoubleCountTrueConvGammas,posDaughter->GetMother(0))){
2209  if(!fDoLightOutput) fHistoDoubleCountTrueConvGammaRPt[fiCut]->Fill(TruePhotonCandidate->GetConversionRadius(),TruePhotonCandidate->Pt(),fWeightJetJetMC);
2210  FillMultipleCountMap(fMapMultipleCountTrueConvGammas,posDaughter->GetMother(0));
2211  }
2212  }
2213  Bool_t isPrimary = ((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsConversionPrimaryESD( fMCEvent, posDaughter->GetMother(0), mcProdVtxX, mcProdVtxY, mcProdVtxZ);
2214  if(isPrimary){
2215  // Count just primary MC Gammas as true --> For Ratio esdtruegamma / mcconvgamma
2217  if(!fDoLightOutput){
2218  fHistoTruePrimaryConvGammaPt[fiCut]->Fill(TruePhotonCandidate->Pt(),fWeightJetJetMC);
2219  }
2220  }
2221  // (Not Filled for i6, Extra Signal Gamma (parambox) are secondary)
2222  }
2223  TruePhotonCandidate->SetIsTrueConvertedPhoton();
2224  return;
2225 }
2226 //________________________________________________________________________
2228 
2229  const AliVVertex* primVtxMC = fMCEvent->GetPrimaryVertex();
2230  Double_t mcProdVtxX = primVtxMC->GetX();
2231  Double_t mcProdVtxY = primVtxMC->GetY();
2232  Double_t mcProdVtxZ = primVtxMC->GetZ();
2233 
2234  TClonesArray *AODMCTrackArray = dynamic_cast<TClonesArray*>(fInputEvent->FindListObject(AliAODMCParticle::StdBranchName()));
2235  if (AODMCTrackArray == NULL) return;
2236 
2237  // Loop over all primary MC particle
2238  for(Long_t i = 0; i < AODMCTrackArray->GetEntriesFast(); i++) {
2239 
2240  AliAODMCParticle* particle = static_cast<AliAODMCParticle*>(AODMCTrackArray->At(i));
2241  if (!particle) continue;
2242 
2243  Bool_t isPrimary = ((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsConversionPrimaryAOD(fInputEvent, particle, mcProdVtxX, mcProdVtxY, mcProdVtxZ);
2244  if (isPrimary) {
2245 
2246  Int_t isMCFromMBHeader = -1;
2247  if(((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetSignalRejection() != 0){
2248  isMCFromMBHeader = ((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsParticleFromBGEvent(i, fMCEvent, fInputEvent);
2249  if(isMCFromMBHeader == 0 && ((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetSignalRejection() != 3) continue;
2250  }
2251 
2252  if (fMesonRecoMode < 2){
2253  if(!((AliConversionPhotonCuts*)fCutArray->At(fiCut))->InPlaneOutOfPlaneCut(particle->Phi(),fEventPlaneAngle,kFALSE)) continue;
2254  }
2255 
2256  Double_t mesonY = 1.e30;
2257  Double_t ratio = 0;
2258  if (particle->E() != TMath::Abs(particle->Pz())){
2259  ratio = (particle->E()+particle->Pz()) / (particle->E()-particle->Pz());
2260  }
2261  if( !(ratio <= 0) ){
2262  mesonY = particle->Y()-((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetEtaShift();
2263  }
2264 
2265  // check neutral mesons
2266  if(((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->MesonIsSelectedAODMC(particle,AODMCTrackArray,((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetEtaShift())){
2267  AliAODMCParticle* daughter0 = static_cast<AliAODMCParticle*>(AODMCTrackArray->At(particle->GetDaughterLabel(0)));
2268  AliAODMCParticle* daughter1 = static_cast<AliAODMCParticle*>(AODMCTrackArray->At(particle->GetDaughterLabel(1)));
2269  Float_t weighted= 1;
2270  if(((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsParticleFromBGEvent(i, fMCEvent, fInputEvent)){
2271  if (particle->Pt()>0.005){
2272  weighted= ((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetWeightForMeson(i, 0x0, fInputEvent);
2273  }
2274  }
2275  Double_t alpha = -10;
2276 
2277  if(particle->GetPdgCode() == fMesonPDG){
2278  alpha = (daughter0->E() - daughter1->E())/(daughter0->E() + daughter1->E());
2279  if (fHistoMCMesonPt[fiCut])fHistoMCMesonPt[fiCut]->Fill(particle->Pt(),weighted*fWeightJetJetMC); // All MC Meson
2280  if (fHistoMCMesonWOWeightPt[fiCut]) fHistoMCMesonWOWeightPt[fiCut]->Fill(particle->Pt(),fWeightJetJetMC);
2281  if (fIsMC > 1 && fHistoMCMesonWOEvtWeightPt[fiCut])fHistoMCMesonWOEvtWeightPt[fiCut]->Fill(particle->Pt());
2282  if (fDoMesonQA > 0){
2283  if (fHistoMCMesonPtY[fiCut])fHistoMCMesonPtY[fiCut]->Fill(particle->Pt(),mesonY,weighted*fWeightJetJetMC);
2284  if (fHistoMCMesonPtAlpha[fiCut])fHistoMCMesonPtAlpha[fiCut]->Fill(particle->Pt(),alpha,fWeightJetJetMC);
2285  if (fIsMC == 2 && fHistoMCMesonPtJetPt[fiCut])fHistoMCMesonPtJetPt[fiCut]->Fill(particle->Pt(),((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetMaxPtJet(),fWeightJetJetMC);
2286  }
2287 
2288  // Check the acceptance for both gammas
2289  if (fMesonRecoMode == 0){
2290  if( ((AliConversionPhotonCuts*)fCutArray->At(fiCut))->PhotonIsSelectedAODMC(daughter0,AODMCTrackArray,kFALSE) &&
2291  ((AliConversionPhotonCuts*)fCutArray->At(fiCut))->PhotonIsSelectedAODMC(daughter1,AODMCTrackArray,kFALSE) &&
2292  ((AliConversionPhotonCuts*)fCutArray->At(fiCut))->InPlaneOutOfPlaneCut(daughter0->Phi(),fEventPlaneAngle,kFALSE) &&
2293  ((AliConversionPhotonCuts*)fCutArray->At(fiCut))->InPlaneOutOfPlaneCut(daughter1->Phi(),fEventPlaneAngle,kFALSE)){
2294  // check acceptance of clusters as well, true if one of them points into the Calo acceptance
2295  if (fHistoMCMesonInAccPt[fiCut])fHistoMCMesonInAccPt[fiCut]->Fill(particle->Pt(),weighted*fWeightJetJetMC); // MC Meson with gamma in acc
2296  if (fHistoMCMesonWOWeightInAccPt[fiCut])fHistoMCMesonWOWeightInAccPt[fiCut]->Fill(particle->Pt(),fWeightJetJetMC); // MC Meson with gamma in acc wo weight
2297  if(fIsMC > 1 && fHistoMCMesonWOEvtWeightInAccPt[fiCut])fHistoMCMesonWOEvtWeightInAccPt[fiCut]->Fill(particle->Pt()); // MC Meson with gamma in acc wo any weight
2298  }
2299  } else if (fMesonRecoMode == 1){
2300  if( ((AliConversionPhotonCuts*)fCutArray->At(fiCut))->PhotonIsSelectedAODMC(daughter0,AODMCTrackArray,kFALSE) &&
2301  ((AliConversionPhotonCuts*)fCutArray->At(fiCut))->PhotonIsSelectedAODMC(daughter1,AODMCTrackArray,kFALSE) &&
2302  ((AliConversionPhotonCuts*)fCutArray->At(fiCut))->InPlaneOutOfPlaneCut(daughter0->Phi(),fEventPlaneAngle,kFALSE) &&
2303  ((AliConversionPhotonCuts*)fCutArray->At(fiCut))->InPlaneOutOfPlaneCut(daughter1->Phi(),fEventPlaneAngle,kFALSE)){
2304  // check acceptance of clusters as well, true if one of them points into the Calo acceptance
2305  if( ((AliCaloPhotonCuts*)fClusterCutArray->At(fiCut))->ClusterIsSelectedAODMC(daughter0,AODMCTrackArray) ||
2306  ((AliCaloPhotonCuts*)fClusterCutArray->At(fiCut))->ClusterIsSelectedAODMC(daughter1,AODMCTrackArray) ){
2307  if (fHistoMCMesonInAccPt[fiCut])fHistoMCMesonInAccPt[fiCut]->Fill(particle->Pt(),weighted*fWeightJetJetMC); // MC Meson with gamma in acc
2308  if (fHistoMCMesonWOWeightInAccPt[fiCut])fHistoMCMesonWOWeightInAccPt[fiCut]->Fill(particle->Pt(),fWeightJetJetMC); // MC Meson with gamma in acc wo weight
2309  if(fIsMC > 1 && fHistoMCMesonWOEvtWeightInAccPt[fiCut])fHistoMCMesonWOEvtWeightInAccPt[fiCut]->Fill(particle->Pt()); // MC Meson with gamma in acc wo any weight
2310  }
2311  }
2312  } else if (fMesonRecoMode == 2){
2313  if (((AliCaloPhotonCuts*)fClusterCutArray->At(fiCut))->ClusterIsSelectedAODMC(daughter0,AODMCTrackArray) &&
2314  ((AliCaloPhotonCuts*)fClusterCutArray->At(fiCut))->ClusterIsSelectedAODMC(daughter1,AODMCTrackArray) ){
2315  if (fHistoMCMesonInAccPt[fiCut])fHistoMCMesonInAccPt[fiCut]->Fill(particle->Pt(),weighted*fWeightJetJetMC); // MC Meson with gamma in acc
2316  if (fHistoMCMesonWOWeightInAccPt[fiCut])fHistoMCMesonWOWeightInAccPt[fiCut]->Fill(particle->Pt(),fWeightJetJetMC); // MC Meson with gamma in acc wo weight
2317  if(fIsMC > 1 && fHistoMCMesonWOEvtWeightInAccPt[fiCut])fHistoMCMesonWOEvtWeightInAccPt[fiCut]->Fill(particle->Pt()); // MC Meson with gamma in acc wo any weight
2318  }
2319  }
2320  }
2321  }
2322  }
2323  }
2324 }
2325 
2326 //________________________________________________________________________
2328 
2329  const AliVVertex* primVtxMC = fMCEvent->GetPrimaryVertex();
2330  Double_t mcProdVtxX = primVtxMC->GetX();
2331  Double_t mcProdVtxY = primVtxMC->GetY();
2332  Double_t mcProdVtxZ = primVtxMC->GetZ();
2333  // cout << mcProdVtxX <<"\t" << mcProdVtxY << "\t" << mcProdVtxZ << endl;
2334 
2335  // Loop over all primary MC particle
2336  for(Long_t i = 0; i < fMCEvent->GetNumberOfTracks(); i++) {
2337  if (((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsConversionPrimaryESD( fMCEvent, i, mcProdVtxX, mcProdVtxY, mcProdVtxZ)){
2338  // fill primary histograms
2339  TParticle* particle = (TParticle *)fMCEvent->Particle(i);
2340  if (!particle) continue;
2341 
2342  Int_t isMCFromMBHeader = -1;
2343  if(((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetSignalRejection() != 0){
2344  isMCFromMBHeader = ((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsParticleFromBGEvent(i, fMCEvent, fInputEvent);
2345  if(isMCFromMBHeader == 0 && ((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetSignalRejection() != 3) continue;
2346  }
2347 
2348  if (fMesonRecoMode < 2){
2349  if(!((AliConversionPhotonCuts*)fCutArray->At(fiCut))->InPlaneOutOfPlaneCut(particle->Phi(),fEventPlaneAngle,kFALSE)) continue;
2350  }
2351 
2352  // Fill histograms for other particles
2353  Double_t mesonY = 1.e30;
2354  Double_t ratio = 0;
2355  if (particle->Energy() != TMath::Abs(particle->Pz())){
2356  ratio = (particle->Energy()+particle->Pz()) / (particle->Energy()-particle->Pz());
2357  }
2358  if( !(ratio <= 0) ){
2359  mesonY = particle->Y()-((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetEtaShift();
2360  }
2361 
2362  // check neutral mesons
2363  if(((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->MesonIsSelectedMC(particle,fMCEvent,((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetEtaShift())){
2364  TParticle* daughter0 = (TParticle*)fMCEvent->Particle(particle->GetFirstDaughter());
2365  TParticle* daughter1 = (TParticle*)fMCEvent->Particle(particle->GetLastDaughter());
2366 
2367  Float_t weighted= 1;
2368  if(((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsParticleFromBGEvent(i, fMCEvent, fInputEvent)){
2369  if (particle->Pt()>0.005){
2370  weighted= ((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetWeightForMeson(i, fMCEvent, fInputEvent);
2371  }
2372  }
2373 
2374  Double_t alpha = -10;
2375  if(particle->GetPdgCode() == fMesonPDG){
2376  alpha = (daughter0->Energy() - daughter1->Energy())/(daughter0->Energy() + daughter1->Energy());
2377  if (fHistoMCMesonPt[fiCut])fHistoMCMesonPt[fiCut]->Fill(particle->Pt(),weighted*fWeightJetJetMC); // All MC Meson
2378  if (fHistoMCMesonWOWeightPt[fiCut]) fHistoMCMesonWOWeightPt[fiCut]->Fill(particle->Pt(),fWeightJetJetMC);
2379  if (fIsMC > 1 && fHistoMCMesonWOEvtWeightPt[fiCut])fHistoMCMesonWOEvtWeightPt[fiCut]->Fill(particle->Pt());
2380  if (fDoMesonQA > 0){
2381  if (fHistoMCMesonPtY[fiCut])fHistoMCMesonPtY[fiCut]->Fill(particle->Pt(),mesonY,weighted*fWeightJetJetMC); // All MC Meson
2382  if (fHistoMCMesonPtAlpha[fiCut])fHistoMCMesonPtAlpha[fiCut]->Fill(particle->Pt(),alpha,fWeightJetJetMC); // All MC Meson
2383  if (fIsMC == 2 && fHistoMCMesonPtJetPt[fiCut])fHistoMCMesonPtJetPt[fiCut]->Fill(particle->Pt(),((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetMaxPtJet(),fWeightJetJetMC);
2384  }
2385 
2386  // Check the acceptance for both gammas & whether they are counted as primaries as well
2387  Bool_t kDaughter0IsPrim = ((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsConversionPrimaryESD( fMCEvent, particle->GetFirstDaughter(), mcProdVtxX, mcProdVtxY, mcProdVtxZ);
2388  Bool_t kDaughter1IsPrim = ((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsConversionPrimaryESD( fMCEvent, particle->GetLastDaughter(), mcProdVtxX, mcProdVtxY, mcProdVtxZ);
2389  if( kDaughter0IsPrim && kDaughter1IsPrim ){
2390  if (fMesonRecoMode == 0){
2391  if( ((AliConversionPhotonCuts*)fCutArray->At(fiCut))->PhotonIsSelectedMC(daughter0,fMCEvent,kFALSE) &&
2392  ((AliConversionPhotonCuts*)fCutArray->At(fiCut))->PhotonIsSelectedMC(daughter1,fMCEvent,kFALSE) &&
2393  ((AliConversionPhotonCuts*)fCutArray->At(fiCut))->InPlaneOutOfPlaneCut(daughter0->Phi(),fEventPlaneAngle,kFALSE) &&
2394  ((AliConversionPhotonCuts*)fCutArray->At(fiCut))->InPlaneOutOfPlaneCut(daughter1->Phi(),fEventPlaneAngle,kFALSE)){
2395  // check acceptance of clusters as well, true if one of them points into the Calo acceptance
2396  if (fHistoMCMesonInAccPt[fiCut])fHistoMCMesonInAccPt[fiCut]->Fill(particle->Pt(),weighted*fWeightJetJetMC); // MC Meson with gamma in acc
2397  if (fHistoMCMesonWOWeightInAccPt[fiCut])fHistoMCMesonWOWeightInAccPt[fiCut]->Fill(particle->Pt(),fWeightJetJetMC); // MC Meson with gamma in acc wo weighting
2398  if(fIsMC > 1 && fHistoMCMesonWOEvtWeightInAccPt[fiCut])fHistoMCMesonWOEvtWeightInAccPt[fiCut]->Fill(particle->Pt()); // MC Meson with gamma in acc wo any weight
2399  }
2400  } else if (fMesonRecoMode == 1){
2401  if (((AliConversionPhotonCuts*)fCutArray->At(fiCut))->PhotonIsSelectedMC(daughter0,fMCEvent,kFALSE) &&
2402  ((AliConversionPhotonCuts*)fCutArray->At(fiCut))->PhotonIsSelectedMC(daughter1,fMCEvent,kFALSE) &&
2403  ((AliConversionPhotonCuts*)fCutArray->At(fiCut))->InPlaneOutOfPlaneCut(daughter0->Phi(),fEventPlaneAngle,kFALSE) &&
2404  ((AliConversionPhotonCuts*)fCutArray->At(fiCut))->InPlaneOutOfPlaneCut(daughter1->Phi(),fEventPlaneAngle,kFALSE)){
2405  // check acceptance of clusters as well, true if one of them points into the Calo acceptance
2406  if (((AliCaloPhotonCuts*)fClusterCutArray->At(fiCut))->ClusterIsSelectedMC(daughter0,fMCEvent) ||
2407  ((AliCaloPhotonCuts*)fClusterCutArray->At(fiCut))->ClusterIsSelectedMC(daughter1,fMCEvent) ){
2408  if (fHistoMCMesonInAccPt[fiCut])fHistoMCMesonInAccPt[fiCut]->Fill(particle->Pt(),weighted*fWeightJetJetMC); // MC Meson with gamma in acc
2409  if (fHistoMCMesonWOWeightInAccPt[fiCut])fHistoMCMesonWOWeightInAccPt[fiCut]->Fill(particle->Pt(),fWeightJetJetMC); // MC Meson with gamma in acc wo weighting
2410  if(fIsMC > 1 && fHistoMCMesonWOEvtWeightInAccPt[fiCut])fHistoMCMesonWOEvtWeightInAccPt[fiCut]->Fill(particle->Pt()); // MC Meson with gamma in acc wo any weight
2411  }
2412  }
2413  } else if (fMesonRecoMode == 2){
2414  if (((AliCaloPhotonCuts*)fClusterCutArray->At(fiCut))->ClusterIsSelectedMC(daughter0,fMCEvent) &&
2415  ((AliCaloPhotonCuts*)fClusterCutArray->At(fiCut))->ClusterIsSelectedMC(daughter1,fMCEvent) ){
2416  if (fHistoMCMesonInAccPt[fiCut])fHistoMCMesonInAccPt[fiCut]->Fill(particle->Pt(),weighted*fWeightJetJetMC); // MC Meson with gamma in acc
2417  if (fHistoMCMesonWOWeightInAccPt[fiCut])fHistoMCMesonWOWeightInAccPt[fiCut]->Fill(particle->Pt(),fWeightJetJetMC); // MC Meson with gamma in acc wo weighting
2418  if(fIsMC > 1 && fHistoMCMesonWOEvtWeightInAccPt[fiCut])fHistoMCMesonWOEvtWeightInAccPt[fiCut]->Fill(particle->Pt()); // MC Meson with gamma in acc wo any weight
2419  }
2420  }
2421  }
2422  }
2423  }
2424  }
2425  }
2426 }
2427 
2428 
2429 //________________________________________________________________________
2430 // function to reject photons in specific invariant mass window
2431 //________________________________________________________________________
2433  TClonesArray * arrClustersMesonCand = NULL;
2434  if(fCorrTaskSetting.CompareTo("") && fMesonRecoMode > 0)
2435  arrClustersMesonCand = dynamic_cast<TClonesArray*>(fInputEvent->FindListObject(Form("%sClustersBranch",fCorrTaskSetting.Data())));
2436 
2437  if (fMesonRecoMode == 0){
2438  if(fGammaCandidates->GetEntries()>1){
2439  for(Int_t firstGammaIndex=0;firstGammaIndex<fGammaCandidates->GetEntries()-1;firstGammaIndex++){
2440  AliAODConversionPhoton *gamma0=dynamic_cast<AliAODConversionPhoton*>(fGammaCandidates->At(firstGammaIndex));
2441  if (gamma0==NULL) continue;
2442  for(Int_t secondGammaIndex=firstGammaIndex+1;secondGammaIndex<fGammaCandidates->GetEntries();secondGammaIndex++){
2443  AliAODConversionPhoton *gamma1=dynamic_cast<AliAODConversionPhoton*>(fGammaCandidates->At(secondGammaIndex));
2444  //Check for same Electron ID
2445  if (gamma1==NULL) continue;
2446  if(gamma0->GetTrackLabelPositive() == gamma1->GetTrackLabelPositive() ||
2447  gamma0->GetTrackLabelNegative() == gamma1->GetTrackLabelNegative() ||
2448  gamma0->GetTrackLabelNegative() == gamma1->GetTrackLabelPositive() ||
2449  gamma0->GetTrackLabelPositive() == gamma1->GetTrackLabelNegative() ) continue;
2450 
2451  AliAODConversionMother *mesonCand = new AliAODConversionMother(gamma0,gamma1);
2452  mesonCand->SetLabels(firstGammaIndex,secondGammaIndex);
2453  mesonCand->CalculateDistanceOfClossetApproachToPrimVtx(fInputEvent->GetPrimaryVertex());
2454 
2455  if((((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->MesonIsSelected(mesonCand,kTRUE,((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetEtaShift()))){
2456  if (!(((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->MesonIsSelectedByMassCut(mesonCand, 0))){
2457  gamma0->SetUseForMesonPair(kFALSE);
2458  gamma1->SetUseForMesonPair(kFALSE);
2459  fHistoMotherInvMassRejected[fiCut]->Fill(mesonCand->M());
2460  }
2461  }
2462  delete mesonCand;
2463  mesonCand=0x0;
2464  }
2465  }
2466  }
2467  } else if (fMesonRecoMode == 1){
2468  if(fGammaCandidates->GetEntries()>0){
2469  for(Int_t firstGammaIndex=0;firstGammaIndex<fGammaCandidates->GetEntries();firstGammaIndex++){
2470  AliAODConversionPhoton *gamma0=dynamic_cast<AliAODConversionPhoton*>(fGammaCandidates->At(firstGammaIndex));
2471  if (gamma0==NULL) continue;
2472 
2473  for(Int_t secondGammaIndex=0;secondGammaIndex<fClusterCandidates->GetEntries();secondGammaIndex++){
2474  Bool_t matched = kFALSE;
2475  AliAODConversionPhoton *gamma1=dynamic_cast<AliAODConversionPhoton*>(fClusterCandidates->At(secondGammaIndex));
2476  if (gamma1==NULL) continue;
2477 
2478  if (gamma1->GetIsCaloPhoton() > 0){
2479  AliVCluster* cluster = NULL;
2480  if(fInputEvent->IsA()==AliESDEvent::Class()){
2481  if(arrClustersMesonCand)
2482  cluster = new AliESDCaloCluster(*(AliESDCaloCluster*)arrClustersMesonCand->At(gamma1->GetCaloClusterRef()));
2483  else
2484  cluster = fInputEvent->GetCaloCluster(gamma1->GetCaloClusterRef());
2485  } else if(fInputEvent->IsA()==AliAODEvent::Class()){
2486  if(arrClustersMesonCand)
2487  cluster = new AliAODCaloCluster(*(AliAODCaloCluster*)arrClustersMesonCand->At(gamma1->GetCaloClusterRef()));
2488  else
2489  cluster = fInputEvent->GetCaloCluster(gamma1->GetCaloClusterRef());
2490  }
2491 
2492  matched = ((AliCaloPhotonCuts*)fClusterCutArray->At(fiCut))->MatchConvPhotonToCluster(gamma0,cluster, fInputEvent, fWeightJetJetMC);
2493  if(arrClustersMesonCand) delete cluster;
2494  }
2495 
2496  AliAODConversionMother *mesonCand = new AliAODConversionMother(gamma0,gamma1);
2497  mesonCand->SetLabels(firstGammaIndex,secondGammaIndex);
2498 
2499  if((((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->MesonIsSelected(mesonCand,kTRUE,((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetEtaShift()))){
2500  if (!matched){
2501  if (!(((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->MesonIsSelectedByMassCut(mesonCand, 0))){
2502  gamma0->SetUseForMesonPair(kFALSE);
2503  gamma1->SetUseForMesonPair(kFALSE);
2504  fHistoMotherInvMassRejected[fiCut]->Fill(mesonCand->M());
2505  }
2506  }
2507  }
2508  delete mesonCand;
2509  mesonCand=0x0;
2510  }
2511  }
2512  }
2513  } else if (fMesonRecoMode == 2){
2514  if(fClusterCandidates->GetEntries()>0){
2515  for(Int_t firstGammaIndex=0;firstGammaIndex<fClusterCandidates->GetEntries();firstGammaIndex++){
2516  AliAODConversionPhoton *gamma0=dynamic_cast<AliAODConversionPhoton*>(fClusterCandidates->At(firstGammaIndex));
2517  if (gamma0==NULL) continue;
2518  for(Int_t secondGammaIndex=firstGammaIndex+1;secondGammaIndex<fClusterCandidates->GetEntries();secondGammaIndex++){
2519  AliAODConversionPhoton *gamma1=dynamic_cast<AliAODConversionPhoton*>(fClusterCandidates->At(secondGammaIndex));
2520  if (gamma1==NULL) continue;
2521 
2522  AliAODConversionMother *mesonCand = new AliAODConversionMother(gamma0,gamma1);
2523  mesonCand->SetLabels(firstGammaIndex,secondGammaIndex);
2524 
2525  if((((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->MesonIsSelected(mesonCand, kTRUE, ((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetEtaShift(), gamma0->GetLeadingCellID(), gamma1->GetLeadingCellID(), gamma0->GetIsCaloPhoton(), gamma1->GetIsCaloPhoton()))){
2526  if (!(((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->MesonIsSelectedByMassCut(mesonCand, 0))){
2527  gamma0->SetUseForMesonPair(kFALSE);
2528  gamma1->SetUseForMesonPair(kFALSE);
2529  fHistoMotherInvMassRejected[fiCut]->Fill(mesonCand->M());
2530  }
2531  }
2532  delete mesonCand;
2533  mesonCand=0x0;
2534  }
2535  }
2536  }
2537  }
2538 }
2539 
2540 
2541 //________________________________________________________________________
2543  TClonesArray * arrClustersMesonCand = NULL;
2544  if(fCorrTaskSetting.CompareTo("") && fMesonRecoMode > 0)
2545  arrClustersMesonCand = dynamic_cast<TClonesArray*>(fInputEvent->FindListObject(Form("%sClustersBranch",fCorrTaskSetting.Data())));
2546 
2547  if (fMesonRecoMode == 0){
2548  if(fGammaCandidates->GetEntries()>1){
2549  for(Int_t firstGammaIndex=0;firstGammaIndex<fGammaCandidates->GetEntries()-1;firstGammaIndex++){
2550  AliAODConversionPhoton *gamma0=dynamic_cast<AliAODConversionPhoton*>(fGammaCandidates->At(firstGammaIndex));
2551  if (gamma0==NULL) continue;
2552  for(Int_t secondGammaIndex=firstGammaIndex+1;secondGammaIndex<fGammaCandidates->GetEntries();secondGammaIndex++){
2553  AliAODConversionPhoton *gamma1=dynamic_cast<AliAODConversionPhoton*>(fGammaCandidates->At(secondGammaIndex));
2554  //Check for same Electron ID
2555  if (gamma1==NULL) continue;
2556  if(gamma0->GetTrackLabelPositive() == gamma1->GetTrackLabelPositive() ||
2557  gamma0->GetTrackLabelNegative() == gamma1->GetTrackLabelNegative() ||
2558  gamma0->GetTrackLabelNegative() == gamma1->GetTrackLabelPositive() ||
2559  gamma0->GetTrackLabelPositive() == gamma1->GetTrackLabelNegative() ) continue;
2560 
2561 
2562  AliAODConversionMother *mesonCand = new AliAODConversionMother(gamma0,gamma1);
2563  mesonCand->SetLabels(firstGammaIndex,secondGammaIndex);
2564  mesonCand->CalculateDistanceOfClossetApproachToPrimVtx(fInputEvent->GetPrimaryVertex());
2565 
2566  if((((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->MesonIsSelected(mesonCand,kTRUE,((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetEtaShift()))){
2567  if (!gamma0->GetUseForMesonPair() || !gamma1->GetUseForMesonPair()){
2568  fHistoMotherInvMassRejected[fiCut]->Fill(mesonCand->M());
2569  delete mesonCand;
2570  mesonCand=0x0;
2571  continue;
2572  }
2573  if (fHistoMotherInvMassPt[fiCut]) fHistoMotherInvMassPt[fiCut]->Fill(mesonCand->M(),mesonCand->Pt(), fWeightJetJetMC);
2574 
2575  if (fDoMesonQA > 0){
2576  if (mesonCand->M() > fMesonInvMassWindow[0] && mesonCand->M() < fMesonInvMassWindow[1]){
2577  if (fHistoMotherMesonPtY[fiCut]) fHistoMotherMesonPtY[fiCut]->Fill(mesonCand->Pt(),mesonCand->Rapidity()-((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetEtaShift(),fWeightJetJetMC);
2578  if (fHistoMotherMesonPtAlpha[fiCut]) fHistoMotherMesonPtAlpha[fiCut]->Fill(mesonCand->Pt(),mesonCand->GetAlpha(),fWeightJetJetMC);
2579  if (fHistoMotherMesonPtOpenAngle[fiCut]) fHistoMotherMesonPtOpenAngle[fiCut]->Fill(mesonCand->Pt(),mesonCand->GetOpeningAngle(),fWeightJetJetMC);
2580  }
2581  }
2582  if(fDoTHnSparse && ((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->DoBGCalculation()){
2583  Int_t zbin = 0;
2584  Int_t mbin = 0;
2585  if(((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->BackgroundHandlerType() == 0){
2586  zbin = fBGHandler[fiCut]->GetZBinIndex(fInputEvent->GetPrimaryVertex()->GetZ());
2587  if(((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->UseTrackMultiplicity()){
2588  mbin = fBGHandler[fiCut]->GetMultiplicityBinIndex(fV0Reader->GetNumberOfPrimaryTracks());
2589  }else {
2590  mbin = fBGHandler[fiCut]->GetMultiplicityBinIndex(fGammaCandidates->GetEntries());
2591  }
2592  }else{
2593  zbin = fBGHandlerRP[fiCut]->GetZBinIndex(fInputEvent->GetPrimaryVertex()->GetZ());
2594  if(((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->UseTrackMultiplicity()){
2595  mbin = fBGHandlerRP[fiCut]->GetMultiplicityBinIndex(fV0Reader->GetNumberOfPrimaryTracks());
2596  }else {
2597  mbin = fBGHandlerRP[fiCut]->GetMultiplicityBinIndex(fGammaCandidates->GetEntries());
2598  }
2599  }
2600 
2601  Double_t sparesFill[4] = {mesonCand->M(),mesonCand->Pt(),(Double_t)zbin,(Double_t)mbin};
2602  if (fSparseMotherInvMassPtZM[fiCut]) fSparseMotherInvMassPtZM[fiCut]->Fill(sparesFill,1);
2603  }
2604 
2605 
2606  if( fIsMC > 0 ){
2607  if(fInputEvent->IsA()==AliESDEvent::Class())
2608  ProcessTrueMesonCandidatesConv(mesonCand, gamma0, gamma1);
2609  if(fInputEvent->IsA()==AliAODEvent::Class())
2610  ProcessTrueMesonCandidatesConvAOD(mesonCand, gamma0, gamma1);
2611  }
2612  }
2613  delete mesonCand;
2614  mesonCand=0x0;
2615  }
2616  }
2617  }
2618  } else if (fMesonRecoMode == 1){
2619  if(fGammaCandidates->GetEntries()>0){
2620  for(Int_t firstGammaIndex=0;firstGammaIndex<fGammaCandidates->GetEntries();firstGammaIndex++){
2621  AliAODConversionPhoton *gamma0=dynamic_cast<AliAODConversionPhoton*>(fGammaCandidates->At(firstGammaIndex));
2622  if (gamma0==NULL) continue;
2623 
2624  for(Int_t secondGammaIndex=0;secondGammaIndex<fClusterCandidates->GetEntries();secondGammaIndex++){
2625  Bool_t matched = kFALSE;
2626  AliAODConversionPhoton *gamma1=dynamic_cast<AliAODConversionPhoton*>(fClusterCandidates->At(secondGammaIndex));
2627  if (gamma1==NULL) continue;
2628 
2629  if (gamma1->GetIsCaloPhoton() > 0){
2630  AliVCluster* cluster = NULL;
2631  if(fInputEvent->IsA()==AliESDEvent::Class()){
2632  if(arrClustersMesonCand)
2633  cluster = new AliESDCaloCluster(*(AliESDCaloCluster*)arrClustersMesonCand->At(gamma1->GetCaloClusterRef()));
2634  else
2635  cluster = fInputEvent->GetCaloCluster(gamma1->GetCaloClusterRef());
2636  } else if(fInputEvent->IsA()==AliAODEvent::Class()){
2637  if(arrClustersMesonCand)
2638  cluster = new AliAODCaloCluster(*(AliAODCaloCluster*)arrClustersMesonCand->At(gamma1->GetCaloClusterRef()));
2639  else
2640  cluster = fInputEvent->GetCaloCluster(gamma1->GetCaloClusterRef());
2641  }
2642 
2643  matched = ((AliCaloPhotonCuts*)fClusterCutArray->At(fiCut))->MatchConvPhotonToCluster(gamma0,cluster, fInputEvent, fWeightJetJetMC);
2644  if(arrClustersMesonCand) delete cluster;
2645  }
2646 
2647  AliAODConversionMother *mesonCand = new AliAODConversionMother(gamma0,gamma1);
2648  mesonCand->SetLabels(firstGammaIndex,secondGammaIndex);
2649 
2650  if((((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->MesonIsSelected(mesonCand,kTRUE,((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetEtaShift()))){
2651  if (matched){
2652  if(!fDoLightOutput && fMesonRecoMode == 1){
2653  if(fHistoMotherMatchedInvMassPt[fiCut]) fHistoMotherMatchedInvMassPt[fiCut]->Fill(mesonCand->M(),mesonCand->Pt(),fWeightJetJetMC);
2654  }
2655  }else {
2656  if (!gamma0->GetUseForMesonPair() || !gamma1->GetUseForMesonPair()){
2657  fHistoMotherInvMassRejected[fiCut]->Fill(mesonCand->M());
2658  delete mesonCand;
2659  mesonCand=0x0;
2660  continue;
2661  }
2662  if (fHistoMotherInvMassPt[fiCut]) fHistoMotherInvMassPt[fiCut]->Fill(mesonCand->M(),mesonCand->Pt(),fWeightJetJetMC);
2663  }
2664 
2665  // fill new histograms
2666  if (!matched){
2667  if (fDoMesonQA > 0){
2668  if (mesonCand->M() > fMesonInvMassWindow[0] && mesonCand->M() < fMesonInvMassWindow[1]){
2669  if (fHistoMotherMesonPtY[fiCut]) fHistoMotherMesonPtY[fiCut]->Fill(mesonCand->Pt(),mesonCand->Rapidity()-((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetEtaShift(),fWeightJetJetMC);
2670  if (fHistoMotherMesonPtAlpha[fiCut]) fHistoMotherMesonPtAlpha[fiCut]->Fill(mesonCand->Pt(),mesonCand->GetAlpha(),fWeightJetJetMC);
2671  if (fHistoMotherMesonPtOpenAngle[fiCut]) fHistoMotherMesonPtOpenAngle[fiCut]->Fill(mesonCand->Pt(),mesonCand->GetOpeningAngle(),fWeightJetJetMC);
2672  if (fHistoMotherMesonConvPhotonEtaPhi[fiCut]) fHistoMotherMesonConvPhotonEtaPhi[fiCut]->Fill(gamma0->GetPhotonPhi(), gamma0->GetPhotonEta(),fWeightJetJetMC);
2673  }
2674  }
2675  if(fDoTHnSparse && ((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->DoBGCalculation()){
2676  Int_t zbin = 0;
2677  Int_t mbin = 0;
2678 
2679  if(((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->BackgroundHandlerType() == 0){
2680  zbin = fBGHandler[fiCut]->GetZBinIndex(fInputEvent->GetPrimaryVertex()->GetZ());
2681  if(((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->UseTrackMultiplicity()){
2682  mbin = fBGHandler[fiCut]->GetMultiplicityBinIndex(fV0Reader->GetNumberOfPrimaryTracks());
2683  }else {
2684  mbin = fBGHandler[fiCut]->GetMultiplicityBinIndex(fGammaCandidates->GetEntries());
2685  }
2686  }else{
2687  zbin = fBGHandlerRP[fiCut]->GetZBinIndex(fInputEvent->GetPrimaryVertex()->GetZ());
2688  if(((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->UseTrackMultiplicity()){
2689  mbin = fBGHandlerRP[fiCut]->GetMultiplicityBinIndex(fV0Reader->GetNumberOfPrimaryTracks());
2690  }else {
2691  mbin = fBGHandlerRP[fiCut]->GetMultiplicityBinIndex(fGammaCandidates->GetEntries());
2692  }
2693  }
2694  Double_t sparesFill[4] = {mesonCand->M(),mesonCand->Pt(),(Double_t)zbin,(Double_t)mbin};
2695  if (fSparseMotherInvMassPtZM[fiCut]) fSparseMotherInvMassPtZM[fiCut]->Fill(sparesFill,1);
2696  }
2697  }
2698 
2699  if(fIsMC>0){
2700  if(fInputEvent->IsA()==AliESDEvent::Class())
2701  ProcessTrueMesonCandidatesConvCalo(mesonCand, gamma0, gamma1, matched);
2702  if(fInputEvent->IsA()==AliAODEvent::Class())
2703  ProcessTrueMesonCandidatesConvCaloAOD(mesonCand, gamma0, gamma1, matched);
2704  }
2705  }
2706  delete mesonCand;
2707  mesonCand=0x0;
2708  }
2709  }
2710  }
2711  } else if (fMesonRecoMode == 2){
2712  if(fClusterCandidates->GetEntries()>0){
2713  for(Int_t firstGammaIndex=0;firstGammaIndex<fClusterCandidates->GetEntries();firstGammaIndex++){
2714  AliAODConversionPhoton *gamma0=dynamic_cast<AliAODConversionPhoton*>(fClusterCandidates->At(firstGammaIndex));
2715  if (gamma0==NULL) continue;
2716  for(Int_t secondGammaIndex=firstGammaIndex+1;secondGammaIndex<fClusterCandidates->GetEntries();secondGammaIndex++){
2717  AliAODConversionPhoton *gamma1=dynamic_cast<AliAODConversionPhoton*>(fClusterCandidates->At(secondGammaIndex));
2718  if (gamma1==NULL) continue;
2719 
2720  AliAODConversionMother *mesonCand = new AliAODConversionMother(gamma0,gamma1);
2721  mesonCand->SetLabels(firstGammaIndex,secondGammaIndex);
2722 
2723  if((((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->MesonIsSelected(mesonCand, kTRUE, ((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetEtaShift(), gamma0->GetLeadingCellID(), gamma1->GetLeadingCellID(), gamma0->GetIsCaloPhoton(), gamma1->GetIsCaloPhoton()))){
2724  if (!gamma0->GetUseForMesonPair() || !gamma1->GetUseForMesonPair()){
2725  fHistoMotherInvMassRejected[fiCut]->Fill(mesonCand->M());
2726  delete mesonCand;
2727  mesonCand=0x0;
2728  continue;
2729  }
2730  if (fHistoMotherInvMassPt[fiCut]) fHistoMotherInvMassPt[fiCut]->Fill(mesonCand->M(),mesonCand->Pt(), fWeightJetJetMC);
2731  // fill new histograms
2732  if (fDoMesonQA > 0 && fDoMesonQA < 3){
2733  if ( mesonCand->M() > fMesonInvMassWindow[0] && mesonCand->M() < fMesonInvMassWindow[1]){
2734  fHistoMotherMesonPtY[fiCut]->Fill(mesonCand->Pt(),mesonCand->Rapidity()-((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetEtaShift(), fWeightJetJetMC);
2735  fHistoMotherMesonPtAlpha[fiCut]->Fill(mesonCand->Pt(),mesonCand->GetAlpha(), fWeightJetJetMC);
2736  fHistoMotherMesonPtOpenAngle[fiCut]->Fill(mesonCand->Pt(),mesonCand->GetOpeningAngle(), fWeightJetJetMC);
2737  }
2738  }
2739  if(fDoTHnSparse && ((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->DoBGCalculation()){
2740  Int_t zbin = 0;
2741  Int_t mbin = 0;
2742 
2743  if(((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->BackgroundHandlerType() == 0){
2744  zbin = fBGHandler[fiCut]->GetZBinIndex(fInputEvent->GetPrimaryVertex()->GetZ());
2745  if(((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->UseTrackMultiplicity()){
2746  mbin = fBGHandler[fiCut]->GetMultiplicityBinIndex(fV0Reader->GetNumberOfPrimaryTracks());
2747  } else {
2749  }
2750  }
2751  Double_t sparesFill[4] = {mesonCand->M(),mesonCand->Pt(),(Double_t)zbin,(Double_t)mbin};
2752  if (fSparseMotherInvMassPtZM[fiCut]) fSparseMotherInvMassPtZM[fiCut]->Fill(sparesFill,1);
2753  }
2754 
2755  if(fIsMC> 0){
2756  if(fInputEvent->IsA()==AliESDEvent::Class())
2757  ProcessTrueMesonCandidatesCalo(mesonCand,gamma0,gamma1);
2758  if(fInputEvent->IsA()==AliAODEvent::Class())
2759  ProcessTrueMesonCandidatesCaloAOD(mesonCand,gamma0,gamma1);
2760  }
2761  }
2762  delete mesonCand;
2763  mesonCand=0x0;
2764  }
2765  }
2766  }
2767  }
2768 }
2769 
2770 //______________________________________________________________________
2772  AliAODConversionPhoton *TrueGammaCandidate0,
2773  AliAODConversionPhoton *TrueGammaCandidate1,
2774  Bool_t matched )
2775 {
2776  // obtain MC vertex
2777  const AliVVertex* primVtxMC = fMCEvent->GetPrimaryVertex();
2778  Double_t mcProdVtxX = primVtxMC->GetX();
2779  Double_t mcProdVtxY = primVtxMC->GetY();
2780  Double_t mcProdVtxZ = primVtxMC->GetZ();
2781 
2782  // Process True Mesons
2783  if(TrueGammaCandidate0->GetV0Index()<fInputEvent->GetNumberOfV0s()){
2784  Bool_t isTrueMeson = kFALSE;
2785  Int_t gamma0MCLabel = -1;
2786  Int_t gamma0MotherLabel = -1;
2787  if (TrueGammaCandidate0->IsTrueConvertedPhoton()){
2788  gamma0MCLabel = TrueGammaCandidate0->GetMCParticleLabel(fMCEvent);
2789  if(gamma0MCLabel>-1){
2790  TParticle * gammaMC0 = (TParticle*)fMCEvent->Particle(gamma0MCLabel);
2791  gamma0MotherLabel=gammaMC0->GetFirstMother();
2792  }
2793  }
2794  if (TrueGammaCandidate1->GetIsCaloPhoton() == 0) AliFatal("CaloPhotonFlag has not been set. Aborting");
2795 
2796  Int_t gamma1MCLabel = TrueGammaCandidate1->GetCaloPhotonMCLabel(0); // get most probable MC label
2797  Int_t gamma1MotherLabel = -1;
2798  // check if
2799 
2800  TParticle * gammaMC1 = 0x0;
2801  if(gamma1MCLabel != -1){ // Gamma is Combinatorial; MC Particles don't belong to the same Mother
2802  // Daughters Gamma 1
2803  gammaMC1 = (TParticle*)fMCEvent->Particle(gamma1MCLabel);
2804  if (TrueGammaCandidate1->IsLargestComponentPhoton() || TrueGammaCandidate1->IsLargestComponentElectron()){ // largest component is electro magnetic
2805  // get mother of interest (pi0 or eta)
2806  if (TrueGammaCandidate1->IsLargestComponentPhoton()){ // for photons its the direct mother
2807  gamma1MotherLabel=gammaMC1->GetMother(0);
2808  }else if (TrueGammaCandidate1->IsLargestComponentElectron()){ // for electrons its either the direct mother or for conversions the grandmother
2809  if (TrueGammaCandidate1->IsConversion() && gammaMC1->GetMother(0)>-1) gamma1MotherLabel=fMCEvent->Particle(gammaMC1->GetMother(0))->GetMother(0);
2810  else gamma1MotherLabel=gammaMC1->GetMother(0);
2811  }
2812  }
2813  }
2814 
2815  if(gamma0MotherLabel>=0 && gamma0MotherLabel==gamma1MotherLabel){
2816  if(((TParticle*)fMCEvent->Particle(gamma1MotherLabel))->GetPdgCode() == fMesonPDG){
2817  isTrueMeson=kTRUE;
2818  }
2819  }
2820 
2821  if(isTrueMeson ){// True Pion or Eta
2822  if (!matched){
2823  if (fHistoTrueMesonInvMassPt[fiCut]) fHistoTrueMesonInvMassPt[fiCut]->Fill(mesonCand->M(),mesonCand->Pt(),fWeightJetJetMC);
2824  }else{
2825  if (fHistoTrueMesonMatchedInvMassPt[fiCut]) fHistoTrueMesonMatchedInvMassPt[fiCut]->Fill(mesonCand->M(),mesonCand->Pt(),fWeightJetJetMC);
2826  }
2827  if (fDoMesonQA > 0 && fIsMC < 2){
2828  if (TrueGammaCandidate1->IsLargestComponentPhoton() && !matched){
2829  if (fHistoTrueMesonCaloPhotonInvMassPt[fiCut]) fHistoTrueMesonCaloPhotonInvMassPt[fiCut]->Fill(mesonCand->M(),mesonCand->Pt());
2830  }
2831  if (TrueGammaCandidate1->IsLargestComponentElectron() && !matched){
2832  if (fHistoTrueMesonCaloElectronInvMassPt[fiCut]) fHistoTrueMesonCaloElectronInvMassPt[fiCut]->Fill(mesonCand->M(),mesonCand->Pt());
2833  }
2834  if (TrueGammaCandidate1->IsLargestComponentElectron() && TrueGammaCandidate1->IsConversion() ){
2835  if ( !matched){
2836  fHistoTrueMesonCaloConvertedPhotonInvMassPt[fiCut]->Fill(mesonCand->M(),mesonCand->Pt());
2837  }
2838 
2839  if ((TrueGammaCandidate0->GetMCLabelPositive() == gamma1MCLabel || TrueGammaCandidate0->GetMCLabelNegative() == gamma1MCLabel) && isTrueMeson){
2841  }
2842  }
2843  if ((TrueGammaCandidate1->IsMerged() || TrueGammaCandidate1->IsMergedPartConv() || TrueGammaCandidate1->IsDalitzMerged()) && !matched ){
2845  }
2846  if (TrueGammaCandidate1->IsMergedPartConv() && !matched){
2848  }
2849  }
2850  if (!matched){
2851  if (fDoMesonQA > 0){
2852  if (isTrueMeson){
2853  if (mesonCand->M() > fMesonInvMassWindow[0] && mesonCand->M() < fMesonInvMassWindow[1] ){
2854  if (fHistoTrueMesonPtY[fiCut]) fHistoTrueMesonPtY[fiCut]->Fill(mesonCand->Pt(),mesonCand->Rapidity()-((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetEtaShift(),fWeightJetJetMC);
2855  if (fHistoTrueMesonPtAlpha[fiCut]) fHistoTrueMesonPtAlpha[fiCut]->Fill(mesonCand->Pt(),mesonCand->GetAlpha(),fWeightJetJetMC);
2856  if (fHistoTrueMesonPtOpenAngle[fiCut]) fHistoTrueMesonPtOpenAngle[fiCut]->Fill(mesonCand->Pt(),mesonCand->GetOpeningAngle(),fWeightJetJetMC);
2857  if (fHistoTrueMotherMesonConvPhotonEtaPhi[fiCut])fHistoTrueMotherMesonConvPhotonEtaPhi[fiCut]->Fill(TrueGammaCandidate0->GetPhotonPhi(), TrueGammaCandidate0->GetPhotonEta(),fWeightJetJetMC);
2858  }
2859  }
2860  }
2861  Bool_t isPrimary = ((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsConversionPrimaryESD( fMCEvent, gamma0MotherLabel, mcProdVtxX, mcProdVtxY, mcProdVtxZ);
2862  if (isPrimary) {
2863  Float_t weighted= 1;
2864  if(((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsParticleFromBGEvent(gamma1MotherLabel, fMCEvent, fInputEvent)){
2865  if (((TParticle*)fMCEvent->Particle(gamma1MotherLabel))->Pt()>0.005){
2866  weighted= ((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetWeightForMeson(gamma1MotherLabel, fMCEvent, fInputEvent);
2867  // cout << "rec \t " <<gamma1MotherLabel << "\t" << weighted << endl;
2868  }
2869  }
2870  if (isTrueMeson){
2871  if (fHistoTruePrimaryMesonInvMassPt[fiCut]) fHistoTruePrimaryMesonInvMassPt[fiCut]->Fill(mesonCand->M(),mesonCand->Pt(),weighted*fWeightJetJetMC);
2873  if (fProfileTruePrimaryMesonWeightsInvMassPt[fiCut]) fProfileTruePrimaryMesonWeightsInvMassPt[fiCut]->Fill(mesonCand->M(),mesonCand->Pt(),weighted*fWeightJetJetMC);
2875  if (fHistoDoubleCountTrueMesonInvMassPt[fiCut]) fHistoDoubleCountTrueMesonInvMassPt[fiCut]->Fill(mesonCand->M(),mesonCand->Pt(),weighted*fWeightJetJetMC);
2877  }
2878  }
2879 
2880  if (fDoMesonQA > 0 && fIsMC < 2){
2881  if(isTrueMeson){ // Only primary pi0 for resolution
2882  if (fHistoTruePrimaryMesonMCPtResolPt[fiCut]) fHistoTruePrimaryMesonMCPtResolPt[fiCut]->Fill(((TParticle*)fMCEvent->Particle(gamma1MotherLabel))->Pt(),(mesonCand->Pt()-((TParticle*)fMCEvent->Particle(gamma1MotherLabel))->Pt())/((TParticle*)fMCEvent->Particle(gamma1MotherLabel))->Pt(),weighted*fWeightJetJetMC);
2883  }
2884  }
2885  }
2886  }
2887  }else if(!isTrueMeson ){ // Background
2888  if (fDoMesonQA > 1){
2889  if(gamma0MotherLabel>-1 && gamma1MotherLabel>-1){ // Both Tracks are Photons and have a mother but not Meson or Eta
2890  if (fHistoTrueBckGGInvMassPt[fiCut])fHistoTrueBckGGInvMassPt[fiCut]->Fill(mesonCand->M(),mesonCand->Pt(), fWeightJetJetMC);
2891  if(
2892  ( ((TParticle*)fMCEvent->Particle(gamma1MotherLabel))->GetPdgCode() == fMesonPDG
2893  && (TrueGammaCandidate1->IsMerged() || TrueGammaCandidate1->IsMergedPartConv()))
2894  ){
2896  }else if( TrueGammaCandidate1->E()/mesonCand->E() > 0.7 ){
2897  if (fHistoTrueBckAsymEClustersInvMassPt[fiCut]) fHistoTrueBckAsymEClustersInvMassPt[fiCut]->Fill(mesonCand->M(),mesonCand->Pt(), fWeightJetJetMC);
2898  }
2899  }else { // No photon or without mother
2900  if (fHistoTrueBckContInvMassPt[fiCut]) fHistoTrueBckContInvMassPt[fiCut]->Fill(mesonCand->M(),mesonCand->Pt(), fWeightJetJetMC);
2901  }
2902  }
2903  }
2904  if (isTrueMeson && !matched){
2905  fVectorRecTrueMesons.push_back(gamma0MotherLabel);
2906  }
2907  }
2908 }
2909 
2910 //______________________________________________________________________
2912  AliAODConversionPhoton *TrueGammaCandidate0,
2913  AliAODConversionPhoton *TrueGammaCandidate1,
2914  Bool_t matched )
2915 {
2916  const AliVVertex* primVtxMC = fMCEvent->GetPrimaryVertex();
2917  Double_t mcProdVtxX = primVtxMC->GetX();
2918  Double_t mcProdVtxY = primVtxMC->GetY();
2919  Double_t mcProdVtxZ = primVtxMC->GetZ();
2920 
2921  // Process True Mesons
2922  TClonesArray *AODMCTrackArray = dynamic_cast<TClonesArray*>(fInputEvent->FindListObject(AliAODMCParticle::StdBranchName()));
2923  if (AODMCTrackArray == NULL) return;
2924  Bool_t isTrueMeson = kFALSE;
2925 
2926  AliAODMCParticle *positiveMC = static_cast<AliAODMCParticle*>(AODMCTrackArray->At(TrueGammaCandidate0->GetMCLabelPositive()));
2927  AliAODMCParticle *negativeMC = static_cast<AliAODMCParticle*>(AODMCTrackArray->At(TrueGammaCandidate0->GetMCLabelNegative()));
2928 
2929  Int_t gamma0MCLabel = -1;
2930  Int_t gamma0MotherLabel = -1;
2931  if(!positiveMC||!negativeMC)
2932  return;
2933 
2934  if (TrueGammaCandidate0->IsTrueConvertedPhoton()){
2935  gamma0MCLabel = positiveMC->GetMother();
2936  AliAODMCParticle * gammaMC0 = static_cast<AliAODMCParticle*>(AODMCTrackArray->At(gamma0MCLabel));
2937  gamma0MotherLabel=gammaMC0->GetMother();
2938  }
2939 
2940  if (TrueGammaCandidate1->GetIsCaloPhoton() == 0) AliFatal("CaloPhotonFlag has not been set. Aborting");
2941  Int_t gamma1MCLabel = TrueGammaCandidate1->GetCaloPhotonMCLabel(0); // get most probable MC label
2942  Int_t gamma1MotherLabel = -1;
2943  // check if
2944 
2945  AliAODMCParticle * gammaMC1 = 0x0;
2946  if(gamma1MCLabel != -1){ // Gamma is Combinatorial; MC Particles don't belong to the same Mother
2947  // Daughters Gamma 1
2948  gammaMC1 = static_cast<AliAODMCParticle*>(AODMCTrackArray->At(gamma1MCLabel));
2949  if (TrueGammaCandidate1->IsLargestComponentPhoton() || TrueGammaCandidate1->IsLargestComponentElectron()){ // largest component is electro magnetic
2950  // get mother of interest (pi0 or eta)
2951  if (TrueGammaCandidate1->IsLargestComponentPhoton()){ // for photons its the direct mother
2952  gamma1MotherLabel=gammaMC1->GetMother();
2953  }else if (TrueGammaCandidate1->IsLargestComponentElectron()){ // for electrons its either the direct mother or for conversions the grandmother
2954  if (TrueGammaCandidate1->IsConversion()){
2955  AliAODMCParticle * gammaGrandMotherMC1 = static_cast<AliAODMCParticle*>(AODMCTrackArray->At(gammaMC1->GetMother()));
2956  gamma1MotherLabel=gammaGrandMotherMC1->GetMother();
2957  }else gamma1MotherLabel=gammaMC1->GetMother();
2958  }
2959  }
2960  }
2961 
2962  if(gamma0MotherLabel>=0 && gamma0MotherLabel==gamma1MotherLabel){
2963  if(((AliAODMCParticle*)AODMCTrackArray->At(gamma1MotherLabel))->GetPdgCode() == fMesonPDG){
2964  isTrueMeson=kTRUE;
2965  }
2966  }
2967 
2968  if(isTrueMeson ){// True Pion or Eta
2969  if (!matched){
2970  if (fHistoTrueMesonInvMassPt[fiCut]) fHistoTrueMesonInvMassPt[fiCut]->Fill(mesonCand->M(),mesonCand->Pt(),fWeightJetJetMC);
2971 
2972  }else{
2973  if (fHistoTrueMesonMatchedInvMassPt[fiCut]) fHistoTrueMesonMatchedInvMassPt[fiCut]->Fill(mesonCand->M(),mesonCand->Pt(),fWeightJetJetMC);
2974  }
2975  if (fDoMesonQA > 0 && fIsMC < 2){
2976  if (TrueGammaCandidate1->IsLargestComponentPhoton() && !matched){
2977  if (fHistoTrueMesonCaloPhotonInvMassPt[fiCut]) fHistoTrueMesonCaloPhotonInvMassPt[fiCut]->Fill(mesonCand->M(),mesonCand->Pt());
2978  }
2979  if (TrueGammaCandidate1->IsLargestComponentElectron() && !matched) {
2980  if (fHistoTrueMesonCaloElectronInvMassPt[fiCut]) fHistoTrueMesonCaloElectronInvMassPt[fiCut]->Fill(mesonCand->M(),mesonCand->Pt());
2981  }
2982  if (TrueGammaCandidate1->IsLargestComponentElectron() && TrueGammaCandidate1->IsConversion()){
2983  if (!matched)fHistoTrueMesonCaloConvertedPhotonInvMassPt[fiCut]->Fill(mesonCand->M(),mesonCand->Pt());
2984  if ((TrueGammaCandidate0->GetMCLabelPositive() == gamma1MCLabel || TrueGammaCandidate0->GetMCLabelNegative() == gamma1MCLabel) && isTrueMeson)
2986  }
2987  if ((TrueGammaCandidate1->IsMerged() || TrueGammaCandidate1->IsMergedPartConv() || TrueGammaCandidate1->IsDalitzMerged()) && !matched ){
2989  }
2990  if (TrueGammaCandidate1->IsMergedPartConv() && !matched) {
2992  }
2993  }
2994 
2995  if ( !matched){
2996  if (fDoMesonQA > 0){
2997  if (mesonCand->M() > fMesonInvMassWindow[0] && mesonCand->M() < fMesonInvMassWindow[1] ){
2998  if (fHistoTrueMesonPtY[fiCut]) fHistoTrueMesonPtY[fiCut]->Fill(mesonCand->Pt(),mesonCand->Rapidity()-((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetEtaShift(),fWeightJetJetMC);
2999  if (fHistoTrueMesonPtAlpha[fiCut]) fHistoTrueMesonPtAlpha[fiCut]->Fill(mesonCand->Pt(),mesonCand->GetAlpha(),fWeightJetJetMC);
3000  if (fHistoTrueMesonPtOpenAngle[fiCut]) fHistoTrueMesonPtOpenAngle[fiCut]->Fill(mesonCand->Pt(),mesonCand->GetOpeningAngle(),fWeightJetJetMC);
3001  if (fHistoTrueMotherMesonConvPhotonEtaPhi[fiCut]) fHistoTrueMotherMesonConvPhotonEtaPhi[fiCut]->Fill(TrueGammaCandidate0->GetPhotonPhi(), TrueGammaCandidate0->GetPhotonEta(),fWeightJetJetMC);
3002  }
3003  }
3004 
3005  Bool_t isPrimary = ((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsConversionPrimaryAOD(fInputEvent, static_cast<AliAODMCParticle*>(AODMCTrackArray->At(gamma0MotherLabel)), mcProdVtxX, mcProdVtxY, mcProdVtxZ);
3006  if(isPrimary){
3007  // Only primary pi0 for efficiency calculation
3008  Float_t weighted= 1;
3009  if(((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsParticleFromBGEvent(gamma1MotherLabel, 0x0, fInputEvent)){
3010  if (static_cast<AliAODMCParticle*>(AODMCTrackArray->At(gamma1MotherLabel))->Pt()>0.005){
3011  weighted= ((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetWeightForMeson(gamma1MotherLabel, 0x0, fInputEvent);
3012  // cout << "rec \t " <<gamma1MotherLabel << "\t" << weighted << endl;
3013  }
3014  }
3015  if (isTrueMeson){
3016  if (fHistoTruePrimaryMesonInvMassPt[fiCut]) fHistoTruePrimaryMesonInvMassPt[fiCut]->Fill(mesonCand->M(),mesonCand->Pt(),weighted*fWeightJetJetMC);
3018  if (fProfileTruePrimaryMesonWeightsInvMassPt[fiCut]) fProfileTruePrimaryMesonWeightsInvMassPt[fiCut]->Fill(mesonCand->M(),mesonCand->Pt(),weighted*fWeightJetJetMC);
3020  if (fHistoDoubleCountTrueMesonInvMassPt[fiCut]) fHistoDoubleCountTrueMesonInvMassPt[fiCut]->Fill(mesonCand->M(),mesonCand->Pt(),weighted*fWeightJetJetMC);
3022  }
3023  }
3024  if (fDoMesonQA > 0 && fIsMC < 2){
3025  if(isTrueMeson){ // Only primary pi0 for resolution
3026  if (fHistoTruePrimaryMesonMCPtResolPt[fiCut]) fHistoTruePrimaryMesonMCPtResolPt[fiCut]->Fill(static_cast<AliAODMCParticle*>(AODMCTrackArray->At(gamma1MotherLabel))->Pt(),
3027  (mesonCand->Pt()-static_cast<AliAODMCParticle*>(AODMCTrackArray->At(gamma1MotherLabel))->Pt())/static_cast<AliAODMCParticle*>(AODMCTrackArray->At(gamma1MotherLabel))->Pt(),weighted*fWeightJetJetMC);
3028 
3029  }
3030  }
3031  }
3032  }
3033  }else if(!isTrueMeson ) { // Background
3034  if (fDoMesonQA > 1){
3035  if(gamma0MotherLabel>-1 && gamma1MotherLabel>-1){ // Both Tracks are Photons and have a mother but not Meson or Eta
3036  if (fHistoTrueBckGGInvMassPt[fiCut]) fHistoTrueBckGGInvMassPt[fiCut]->Fill(mesonCand->M(),mesonCand->Pt());
3037 
3038  if( (((AliAODMCParticle*)AODMCTrackArray->At(gamma1MotherLabel))->GetPdgCode() == fMesonPDG
3039  && ((TrueGammaCandidate1->IsMerged() || TrueGammaCandidate1->IsMergedPartConv())))
3040  ){
3042  }else if( TrueGammaCandidate1->E()/mesonCand->E() > 0.7 ){
3043  if (fHistoTrueBckAsymEClustersInvMassPt[fiCut]) fHistoTrueBckAsymEClustersInvMassPt[fiCut]->Fill(mesonCand->M(),mesonCand->Pt(), fWeightJetJetMC);
3044  }
3045  }else { // No photon or without mother
3046  if (fHistoTrueBckContInvMassPt[fiCut]) fHistoTrueBckContInvMassPt[fiCut]->Fill(mesonCand->M(),mesonCand->Pt());
3047  }
3048  }
3049  }
3050  if (isTrueMeson && !matched){
3051  fVectorRecTrueMesons.push_back(gamma0MotherLabel);
3052  }
3053 }
3054 
3055 
3056 //______________________________________________________________________
3058  AliAODConversionPhoton *TrueGammaCandidate0,
3059  AliAODConversionPhoton *TrueGammaCandidate1)
3060 {
3061  // Process True Mesons
3062  const AliVVertex* primVtxMC = fMCEvent->GetPrimaryVertex();
3063  Double_t mcProdVtxX = primVtxMC->GetX();
3064  Double_t mcProdVtxY = primVtxMC->GetY();
3065  Double_t mcProdVtxZ = primVtxMC->GetZ();
3066 
3067  Bool_t isTrueMeson = kFALSE;
3068  Int_t gamma0MCLabel = TrueGammaCandidate0->GetCaloPhotonMCLabel(0); // get most probable MC label
3069  Int_t gamma0MotherLabel = -1;
3070 
3071  TParticle * gammaMC0 = 0x0;
3072  if(gamma0MCLabel != -1){ // Gamma is Combinatorial; MC Particles don't belong to the same Mother
3073  gammaMC0 = (TParticle*)fMCEvent->Particle(gamma0MCLabel);
3074  if (TrueGammaCandidate0->IsLargestComponentPhoton() || TrueGammaCandidate0->IsLargestComponentElectron()){ // largest component is electro magnetic
3075  // get mother of interest (pi0 or eta)
3076  if (TrueGammaCandidate0->IsLargestComponentPhoton()){ // for photons its the direct mother
3077  gamma0MotherLabel=gammaMC0->GetMother(0);
3078  } else if (TrueGammaCandidate0->IsLargestComponentElectron()){ // for electrons its either the direct mother or for conversions the grandmother
3079  if (TrueGammaCandidate0->IsConversion() && (gammaMC0->GetMother(0) > -1)){
3080  gamma0MotherLabel=fMCEvent->Particle(gammaMC0->GetMother(0))->GetMother(0);
3081  } else {
3082  gamma0MotherLabel=gammaMC0->GetMother(0);
3083  }
3084  }
3085  }
3086  }
3087  if (TrueGammaCandidate1->GetIsCaloPhoton() == 0) AliFatal("CaloPhotonFlag has not been set. Aborting");
3088 
3089  Int_t gamma1MCLabel = TrueGammaCandidate1->GetCaloPhotonMCLabel(0); // get most probable MC label
3090  Int_t gamma1MotherLabel = -1;
3091  // check if
3092 
3093  TParticle * gammaMC1 = 0x0;
3094  if(gamma1MCLabel != -1){ // Gamma is Combinatorial; MC Particles don't belong to the same Mother
3095  // Daughters Gamma 1
3096  gammaMC1 = (TParticle*)fMCEvent->Particle(gamma1MCLabel);
3097  if (TrueGammaCandidate1->IsLargestComponentPhoton() || TrueGammaCandidate1->IsLargestComponentElectron()){ // largest component is electro magnetic
3098  // get mother of interest (pi0 or eta)
3099  if (TrueGammaCandidate1->IsLargestComponentPhoton()){ // for photons its the direct mother
3100  gamma1MotherLabel = gammaMC1->GetMother(0);
3101  } else if (TrueGammaCandidate1->IsLargestComponentElectron()){ // for electrons its either the direct mother or for conversions the grandmother
3102  if (TrueGammaCandidate1->IsConversion()){
3103  if(gammaMC1->GetMother(0) > -1) gamma1MotherLabel = fMCEvent->Particle(gammaMC1->GetMother(0))->GetMother(0);
3104  } else {
3105  gamma1MotherLabel=gammaMC1->GetMother(0);
3106  }
3107  }
3108  }
3109  }
3110 
3111  if(gamma0MotherLabel>=0 && gamma0MotherLabel==gamma1MotherLabel){
3112  if(((TParticle*)fMCEvent->Particle(gamma1MotherLabel))->GetPdgCode() == fMesonPDG){
3113  isTrueMeson=kTRUE;
3114  }
3115  }
3116 
3117  if(isTrueMeson ){// True Meson
3118  if (fHistoTrueMesonInvMassPt[fiCut]) fHistoTrueMesonInvMassPt[fiCut]->Fill(mesonCand->M(),mesonCand->Pt(), fWeightJetJetMC);
3119  if (fDoMesonQA > 0 && fDoMesonQA < 3 && fIsMC < 2){
3120  // both gammas are real gammas
3121  if (TrueGammaCandidate0->IsLargestComponentPhoton() && TrueGammaCandidate1->IsLargestComponentPhoton()) {
3122  if (fHistoTrueMesonCaloPhotonInvMassPt[fiCut]) fHistoTrueMesonCaloPhotonInvMassPt[fiCut]->Fill(mesonCand->M(),mesonCand->Pt());
3123  }
3124  // both particles are electrons
3125  if (TrueGammaCandidate0->IsLargestComponentElectron() && TrueGammaCandidate1->IsLargestComponentElectron() ) {
3126  if (fHistoTrueMesonCaloElectronInvMassPt[fiCut]) fHistoTrueMesonCaloElectronInvMassPt[fiCut]->Fill(mesonCand->M(),mesonCand->Pt());
3127  }
3128  // both particles are converted electrons
3129  if ((TrueGammaCandidate0->IsLargestComponentElectron() && TrueGammaCandidate0->IsConversion()) && (TrueGammaCandidate1->IsLargestComponentElectron() && TrueGammaCandidate1->IsConversion()) ){
3130  fHistoTrueMesonCaloConvertedPhotonInvMassPt[fiCut]->Fill(mesonCand->M(),mesonCand->Pt());
3131  }
3132  // 1 gamma is converted the other one is normals
3133  if ( (TrueGammaCandidate0->IsLargestComponentPhoton() && TrueGammaCandidate1->IsLargestComponentElectron() && TrueGammaCandidate1->IsConversion()) ||
3134  (TrueGammaCandidate1->IsLargestComponentPhoton() && TrueGammaCandidate0->IsLargestComponentElectron() && TrueGammaCandidate0->IsConversion())
3135  ) {
3136  fHistoTrueMesonCaloMixedPhotonConvPhotonInvMassPt[fiCut]->Fill(mesonCand->M(),mesonCand->Pt());
3137  }
3138 
3139  // at least one of the photon is merged
3140  if (TrueGammaCandidate0->IsMerged() || TrueGammaCandidate0->IsMergedPartConv() || TrueGammaCandidate0->IsDalitzMerged() || TrueGammaCandidate1->IsMerged() || TrueGammaCandidate1->IsMergedPartConv() || TrueGammaCandidate1->IsDalitzMerged() ){
3141  if (fHistoTrueMesonCaloMergedClusterInvMassPt[fiCut]) fHistoTrueMesonCaloMergedClusterInvMassPt[fiCut]->Fill(mesonCand->M(),mesonCand->Pt());
3142  }
3143  // at least one of the photon is merged and part conv
3144  if (TrueGammaCandidate1->IsMergedPartConv() || TrueGammaCandidate0->IsMergedPartConv()) {
3146  }
3147  }
3148 
3149  if (fDoMesonQA > 0 && fDoMesonQA < 3){
3150  if ( mesonCand->M() > fMesonInvMassWindow[0] && mesonCand->M() < fMesonInvMassWindow[1] ){
3151  if (fHistoTrueMesonPtY[fiCut]) fHistoTrueMesonPtY[fiCut]->Fill(mesonCand->Pt(),mesonCand->Rapidity()-((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetEtaShift(), fWeightJetJetMC);
3152  if (fHistoTrueMesonPtAlpha[fiCut]) fHistoTrueMesonPtAlpha[fiCut]->Fill(mesonCand->Pt(),mesonCand->GetAlpha(), fWeightJetJetMC);
3153  if (fHistoTrueMesonPtOpenAngle[fiCut]) fHistoTrueMesonPtOpenAngle[fiCut]->Fill(mesonCand->Pt(),mesonCand->GetOpeningAngle(), fWeightJetJetMC);
3154  }
3155 
3156  }
3157  Bool_t isPrimary = ((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsConversionPrimaryESD( fMCEvent, gamma0MotherLabel, mcProdVtxX, mcProdVtxY, mcProdVtxZ);
3158  if(isPrimary){ // Only primary pi0 for efficiency calculation
3159  // filling primary histograms
3160  Float_t weighted= 1;
3161  if(((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsParticleFromBGEvent(gamma1MotherLabel, fMCEvent, fInputEvent)){
3162  if (((TParticle*)fMCEvent->Particle(gamma1MotherLabel))->Pt()>0.005){
3163  weighted= ((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetWeightForMeson(gamma1MotherLabel, fMCEvent, fInputEvent);
3164  // cout << "rec \t " <<gamma1MotherLabel << "\t" << weighted << endl;
3165  }
3166  }
3167  if (isTrueMeson){
3168  if (fHistoTruePrimaryMesonInvMassPt[fiCut]) fHistoTruePrimaryMesonInvMassPt[fiCut]->Fill(mesonCand->M(),mesonCand->Pt(),weighted* fWeightJetJetMC);
3170  if (fProfileTruePrimaryMesonWeightsInvMassPt[fiCut]) fProfileTruePrimaryMesonWeightsInvMassPt[fiCut]->Fill(mesonCand->M(),mesonCand->Pt(),weighted* fWeightJetJetMC);
3171  if (CheckVectorForDoubleCount(fVectorDoubleCountTrueMesons,gamma0MotherLabel) && fHistoDoubleCountTrueMesonInvMassPt[fiCut]) fHistoDoubleCountTrueMesonInvMassPt[fiCut]->Fill(mesonCand->M(),mesonCand->Pt(), weighted*fWeightJetJetMC);
3172  }
3173  if (fDoMesonQA > 0 && fDoMesonQA < 3 && fIsMC<2){
3174  if(isTrueMeson){ // Only primary pi0 for resolution
3175  if (fHistoTruePrimaryMesonMCPtResolPt[fiCut]) fHistoTruePrimaryMesonMCPtResolPt[fiCut]->Fill(((TParticle*)fMCEvent->Particle(gamma1MotherLabel))->Pt(),(mesonCand->Pt()-((TParticle*)fMCEvent->Particle(gamma1MotherLabel))->Pt())/((TParticle*)fMCEvent->Particle(gamma1MotherLabel))->Pt(),weighted* fWeightJetJetMC);
3176  }
3177  }
3178  }
3179  } else if(!isTrueMeson ){ // Background
3180  if (fDoMesonQA > 1 && fDoMesonQA < 3){
3181  if(gamma0MotherLabel>-1 && gamma1MotherLabel>-1){ // Both Tracks are Photons and have a mother but not Pi0 or Eta
3182  if (fHistoTrueBckGGInvMassPt[fiCut]) fHistoTrueBckGGInvMassPt[fiCut]->Fill(mesonCand->M(),mesonCand->Pt(), fWeightJetJetMC);
3183 
3184  if( (((TParticle*)fMCEvent->Particle(gamma0MotherLabel))->GetPdgCode() == fMesonPDG
3185  && (TrueGammaCandidate0->IsMerged() || TrueGammaCandidate0->IsMergedPartConv()))
3186  ||
3187  (((TParticle*)fMCEvent->Particle(gamma1MotherLabel))->GetPdgCode() == fMesonPDG
3188  && (TrueGammaCandidate1->IsMerged() || TrueGammaCandidate1->IsMergedPartConv()))
3189  ){
3191  }else if( (TrueGammaCandidate0->E()/mesonCand->E() > 0.7) || (TrueGammaCandidate1->E()/mesonCand->E() > 0.7) ){
3192  if (fHistoTrueBckAsymEClustersInvMassPt[fiCut]) fHistoTrueBckAsymEClustersInvMassPt[fiCut]->Fill(mesonCand->M(),mesonCand->Pt(), fWeightJetJetMC);
3193  }
3194  } else { // No photon or without mother
3195  if (fHistoTrueBckContInvMassPt[fiCut]) fHistoTrueBckContInvMassPt[fiCut]->Fill(mesonCand->M(),mesonCand->Pt(), fWeightJetJetMC);
3196  }
3197  }
3198  }
3199 }
3200 
3201 //______________________________________________________________________
3203  AliAODConversionPhoton *TrueGammaCandidate0,
3204  AliAODConversionPhoton *TrueGammaCandidate1)
3205 {
3206  const AliVVertex* primVtxMC = fMCEvent->GetPrimaryVertex();
3207  Double_t mcProdVtxX = primVtxMC->GetX();
3208  Double_t mcProdVtxY = primVtxMC->GetY();
3209  Double_t mcProdVtxZ = primVtxMC->GetZ();
3210 
3211  // Process True Mesons
3212  TClonesArray *AODMCTrackArray = dynamic_cast<TClonesArray*>(fInputEvent->FindListObject(AliAODMCParticle::StdBranchName()));
3213  if (AODMCTrackArray == NULL) return;
3214 
3215  Bool_t isTrueMeson = kFALSE;
3216  Int_t gamma0MCLabel = TrueGammaCandidate0->GetCaloPhotonMCLabel(0); // get most probable MC label
3217  Int_t gamma0MotherLabel = -1;
3218 
3219  // check if
3220  AliAODMCParticle * gammaMC0 = 0x0;
3221  if(gamma0MCLabel != -1){ // Gamma is Combinatorial; MC Particles don't belong to the same Mother
3222  // Daughters Gamma 0
3223  gammaMC0 = static_cast<AliAODMCParticle*>(AODMCTrackArray->At(gamma0MCLabel));
3224  if (TrueGammaCandidate0->IsLargestComponentPhoton() || TrueGammaCandidate0->IsLargestComponentElectron()){ // largest component is electro magnetic
3225  // get mother of interest (pi0 or eta)
3226  if (TrueGammaCandidate0->IsLargestComponentPhoton()){ // for photons its the direct mother
3227  gamma0MotherLabel=gammaMC0->GetMother();
3228  } else if (TrueGammaCandidate0->IsLargestComponentElectron()){ // for electrons its either the direct mother or for conversions the grandmother
3229  if (TrueGammaCandidate0->IsConversion()){
3230  AliAODMCParticle * gammaGrandMotherMC0 = static_cast<AliAODMCParticle*>(AODMCTrackArray->At(gammaMC0->GetMother()));
3231  gamma0MotherLabel=gammaGrandMotherMC0->GetMother();
3232  } else gamma0MotherLabel=gammaMC0->GetMother();
3233  }
3234  }
3235  }
3236 
3237  Int_t gamma1MCLabel = TrueGammaCandidate1->GetCaloPhotonMCLabel(0); // get most probable MC label
3238  Int_t gamma1MotherLabel = -1;
3239 
3240  // check if
3241  AliAODMCParticle *gammaMC1 = 0x0;
3242  if(gamma1MCLabel != -1){ // Gamma is Combinatorial; MC Particles don't belong to the same Mother
3243  // Daughters Gamma 1
3244  gammaMC1 = static_cast<AliAODMCParticle*>(AODMCTrackArray->At(gamma1MCLabel));
3245  if (TrueGammaCandidate1->IsLargestComponentPhoton() || TrueGammaCandidate1->IsLargestComponentElectron()){ // largest component is electro magnetic
3246  // get mother of interest (pi0 or eta)
3247  if (TrueGammaCandidate1->IsLargestComponentPhoton()){ // for photons its the direct mother
3248  gamma1MotherLabel=gammaMC1->GetMother();
3249  } else if (TrueGammaCandidate1->IsLargestComponentElectron()){ // for electrons its either the direct mother or for conversions the grandmother
3250  if (TrueGammaCandidate1->IsConversion()){
3251  AliAODMCParticle * gammaGrandMotherMC1 = static_cast<AliAODMCParticle*>(AODMCTrackArray->At(gammaMC1->GetMother()));
3252  gamma1MotherLabel=gammaGrandMotherMC1->GetMother();
3253  } else gamma1MotherLabel=gammaMC1->GetMother();
3254  }
3255  }
3256  }
3257 
3258  if(gamma0MotherLabel>=0 && gamma0MotherLabel==gamma1MotherLabel){
3259  if(((AliAODMCParticle*)AODMCTrackArray->At(gamma1MotherLabel))->GetPdgCode() == fMesonPDG){
3260  isTrueMeson=kTRUE;
3261  }
3262  }
3263 
3264  if(isTrueMeson ){// True Meson
3265  if (fHistoTrueMesonInvMassPt[fiCut]) fHistoTrueMesonInvMassPt[fiCut]->Fill(mesonCand->M(),mesonCand->Pt(), fWeightJetJetMC);
3266  if (fDoMesonQA > 0 && fDoMesonQA < 3 && fIsMC < 2){
3267  // both gammas are real gammas
3268  if (TrueGammaCandidate0->IsLargestComponentPhoton() && TrueGammaCandidate1->IsLargestComponentPhoton()) {
3269  if (fHistoTrueMesonCaloPhotonInvMassPt[fiCut]) fHistoTrueMesonCaloPhotonInvMassPt[fiCut]->Fill(mesonCand->M(),mesonCand->Pt());
3270  }
3271  // both particles are electrons
3272  if (TrueGammaCandidate0->IsLargestComponentElectron() && TrueGammaCandidate1->IsLargestComponentElectron() ) {
3273  if (fHistoTrueMesonCaloElectronInvMassPt[fiCut]) fHistoTrueMesonCaloElectronInvMassPt[fiCut]->Fill(mesonCand->M(),mesonCand->Pt());
3274  }
3275  // both particles are converted electrons
3276  if ((TrueGammaCandidate0->IsLargestComponentElectron() && TrueGammaCandidate0->IsConversion()) && (TrueGammaCandidate1->IsLargestComponentElectron() && TrueGammaCandidate1->IsConversion()) ){
3277  fHistoTrueMesonCaloConvertedPhotonInvMassPt[fiCut]->Fill(mesonCand->M(),mesonCand->Pt());
3278  }
3279  // 1 gamma is converted the other one is normals
3280  if ( (TrueGammaCandidate0->IsLargestComponentPhoton() && TrueGammaCandidate1->IsLargestComponentElectron() && TrueGammaCandidate1->IsConversion()) ||
3281  (TrueGammaCandidate1->IsLargestComponentPhoton() && TrueGammaCandidate0->IsLargestComponentElectron() && TrueGammaCandidate0->IsConversion())
3282  ) {
3283  fHistoTrueMesonCaloMixedPhotonConvPhotonInvMassPt[fiCut]->Fill(mesonCand->M(),mesonCand->Pt());
3284  }
3285 
3286  // at least one of the photon is merged
3287  if (TrueGammaCandidate0->IsMerged() || TrueGammaCandidate0->IsMergedPartConv() || TrueGammaCandidate0->IsDalitzMerged() || TrueGammaCandidate1->IsMerged() || TrueGammaCandidate1->IsMergedPartConv() || TrueGammaCandidate1->IsDalitzMerged() ){
3288  if (fHistoTrueMesonCaloMergedClusterInvMassPt[fiCut]) fHistoTrueMesonCaloMergedClusterInvMassPt[fiCut]->Fill(mesonCand->M(),mesonCand->Pt());
3289  }
3290  // at least one of the photon is merged and part conv
3291  if (TrueGammaCandidate1->IsMergedPartConv() || TrueGammaCandidate0->IsMergedPartConv()) {
3293  }
3294  }
3295 
3296  if (fDoMesonQA > 0 && fDoMesonQA < 3){
3297  if ( mesonCand->M() > fMesonInvMassWindow[0] && mesonCand->M() < fMesonInvMassWindow[1] ){
3298  if (fHistoTrueMesonPtY[fiCut]) fHistoTrueMesonPtY[fiCut]->Fill(mesonCand->Pt(),mesonCand->Rapidity()-((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetEtaShift(), fWeightJetJetMC);
3299  if (fHistoTrueMesonPtAlpha[fiCut]) fHistoTrueMesonPtAlpha[fiCut]->Fill(mesonCand->Pt(),mesonCand->GetAlpha(), fWeightJetJetMC);
3300  if (fHistoTrueMesonPtOpenAngle[fiCut]) fHistoTrueMesonPtOpenAngle[fiCut]->Fill(mesonCand->Pt(),mesonCand->GetOpeningAngle(), fWeightJetJetMC);
3301  }
3302  }
3303 
3304  Bool_t isPrimary = ((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsConversionPrimaryAOD(fInputEvent, static_cast<AliAODMCParticle*>(AODMCTrackArray->At(gamma0MotherLabel)), mcProdVtxX, mcProdVtxY, mcProdVtxZ);
3305  if(isPrimary){ // Only primary pi0 for efficiency calculation
3306  Float_t weighted= 1;
3307  if(((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsParticleFromBGEvent(gamma1MotherLabel, 0x0, fInputEvent)){
3308  if (static_cast<AliAODMCParticle*>(AODMCTrackArray->At(gamma1MotherLabel))->Pt()>0.005){
3309  weighted= ((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetWeightForMeson(gamma1MotherLabel, 0x0, fInputEvent);
3310  }
3311  }
3312 
3313  if (fHistoTruePrimaryMesonInvMassPt[fiCut]) fHistoTruePrimaryMesonInvMassPt[fiCut]->Fill(mesonCand->M(),mesonCand->Pt(),weighted* fWeightJetJetMC);
3315  if (fProfileTruePrimaryMesonWeightsInvMassPt[fiCut]) fProfileTruePrimaryMesonWeightsInvMassPt[fiCut]->Fill(mesonCand->M(),mesonCand->Pt(),weighted* fWeightJetJetMC);
3316  if (CheckVectorForDoubleCount(fVectorDoubleCountTrueMesons,gamma0MotherLabel) && fHistoDoubleCountTrueMesonInvMassPt[fiCut]) fHistoDoubleCountTrueMesonInvMassPt[fiCut]->Fill(mesonCand->M(),mesonCand->Pt(), weighted*fWeightJetJetMC);
3317  if (fDoMesonQA > 0 && fDoMesonQA < 3 && fIsMC<2){
3318  if (fHistoTruePrimaryMesonMCPtResolPt[fiCut]) fHistoTruePrimaryMesonMCPtResolPt[fiCut]->Fill(static_cast<AliAODMCParticle*>(AODMCTrackArray->At(gamma1MotherLabel))->Pt(),
3319  (mesonCand->Pt()-static_cast<AliAODMCParticle*>(AODMCTrackArray->At(gamma1MotherLabel))->Pt())/static_cast<AliAODMCParticle*>(AODMCTrackArray->At(gamma1MotherLabel))->Pt(),weighted* fWeightJetJetMC);
3320  }
3321  }
3322  } else if(!isTrueMeson ) { // Background
3323  if (fDoMesonQA > 1 && fDoMesonQA < 3){
3324  if(gamma0MotherLabel>-1 && gamma1MotherLabel>-1){ // Both Tracks are Photons and have a mother but not Pi0 or Eta
3325  if (fHistoTrueBckGGInvMassPt[fiCut]) fHistoTrueBckGGInvMassPt[fiCut]->Fill(mesonCand->M(),mesonCand->Pt(), fWeightJetJetMC);
3326 
3327  if( (((AliAODMCParticle*)AODMCTrackArray->At(gamma0MotherLabel))->GetPdgCode() == fMesonPDG
3328  && (TrueGammaCandidate0->IsMerged() || TrueGammaCandidate0->IsMergedPartConv()))
3329  ||
3330  (((AliAODMCParticle*)AODMCTrackArray->At(gamma1MotherLabel))->GetPdgCode() == fMesonPDG
3331  && (TrueGammaCandidate1->IsMerged() || TrueGammaCandidate1->IsMergedPartConv()))
3332  ){
3334  }else if( (TrueGammaCandidate0->E()/mesonCand->E() > 0.7) || (TrueGammaCandidate1->E()/mesonCand->E() > 0.7) ){
3335  if (fHistoTrueBckAsymEClustersInvMassPt[fiCut]) fHistoTrueBckAsymEClustersInvMassPt[fiCut]->Fill(mesonCand->M(),mesonCand->Pt(), fWeightJetJetMC);
3336  }
3337  } else { // No photon or without mother
3338  if (fHistoTrueBckContInvMassPt[fiCut]) fHistoTrueBckContInvMassPt[fiCut]->Fill(mesonCand->M(),mesonCand->Pt(), fWeightJetJetMC);
3339  }
3340  }
3341  }
3342 }
3343 
3344 
3345 //______________________________________________________________________
3347  AliAODConversionPhoton *TrueGammaCandidate0,
3348  AliAODConversionPhoton *TrueGammaCandidate1)
3349 {
3350  // Process True Mesons
3351  const AliVVertex* primVtxMC = fMCEvent->GetPrimaryVertex();
3352  Double_t mcProdVtxX = primVtxMC->GetX();
3353  Double_t mcProdVtxY = primVtxMC->GetY();
3354  Double_t mcProdVtxZ = primVtxMC->GetZ();
3355 
3356  if(TrueGammaCandidate0->GetV0Index()<fInputEvent->GetNumberOfV0s()){
3357  Bool_t isTrueMeson = kFALSE;
3358  Bool_t isTrueMesonDalitz = kFALSE;
3359  Bool_t gamma0DalitzCand = kFALSE;
3360  Bool_t gamma1DalitzCand = kFALSE;
3361  Int_t gamma0MCLabel = TrueGammaCandidate0->GetMCParticleLabel(fMCEvent);
3362  Int_t gamma0MotherLabel = -1;
3363  if(gamma0MCLabel != -1){ // Gamma is Combinatorial; MC Particles don't belong to the same Mother
3364  // Daughters Gamma 0
3365  TParticle * negativeMC = (TParticle*)TrueGammaCandidate0->GetNegativeMCDaughter(fMCEvent);
3366  TParticle * positiveMC = (TParticle*)TrueGammaCandidate0->GetPositiveMCDaughter(fMCEvent);
3367  TParticle * gammaMC0 = (TParticle*)fMCEvent->Particle(gamma0MCLabel);
3368  if(TMath::Abs(negativeMC->GetPdgCode())==11 && TMath::Abs(positiveMC->GetPdgCode())==11){ // Electrons ...
3369  if(negativeMC->GetUniqueID() == 5 && positiveMC->GetUniqueID() ==5){ // ... From Conversion ...
3370  if(gammaMC0->GetPdgCode() == 22){ // ... with Gamma Mother
3371  gamma0MotherLabel=gammaMC0->GetFirstMother();
3372  }
3373  }
3374  if(gammaMC0->GetPdgCode() ==fMesonPDG){ // Dalitz candidate
3375  gamma0DalitzCand = kTRUE;
3376  gamma0MotherLabel=-fMesonPDG;
3377  }
3378  }
3379  }
3380  if(TrueGammaCandidate1->GetV0Index()<fInputEvent->GetNumberOfV0s()){
3381  Int_t gamma1MCLabel = TrueGammaCandidate1->GetMCParticleLabel(fMCEvent);
3382  Int_t gamma1MotherLabel = -1;
3383  if(gamma1MCLabel != -1){ // Gamma is Combinatorial; MC Particles don't belong to the same Mother
3384  // Daughters Gamma 1
3385  TParticle * negativeMC = (TParticle*)TrueGammaCandidate1->GetNegativeMCDaughter(fMCEvent);
3386  TParticle * positiveMC = (TParticle*)TrueGammaCandidate1->GetPositiveMCDaughter(fMCEvent);
3387  TParticle * gammaMC1 = (TParticle*)fMCEvent->Particle(gamma1MCLabel);
3388  if(TMath::Abs(negativeMC->GetPdgCode())==11 && TMath::Abs(positiveMC->GetPdgCode())==11){ // Electrons ...
3389  if(negativeMC->GetUniqueID() == 5 && positiveMC->GetUniqueID() ==5){ // ... From Conversion ...
3390  if(gammaMC1->GetPdgCode() == 22){ // ... with Gamma Mother
3391  gamma1MotherLabel=gammaMC1->GetFirstMother();
3392  }
3393  }
3394  if(gammaMC1->GetPdgCode() ==fMesonPDG ){ // Dalitz candidate
3395  gamma1DalitzCand = kTRUE;
3396  gamma1MotherLabel=-fMesonPDG;
3397  }
3398  }
3399  }
3400  if(gamma0MotherLabel>=0 && gamma0MotherLabel==gamma1MotherLabel){
3401  if(((TParticle*)fMCEvent->Particle(gamma1MotherLabel))->GetPdgCode() == fMesonPDG){
3402  isTrueMeson=kTRUE;
3406  }
3407  }
3408  }
3409 
3410  //Identify Dalitz candidate
3411  if (gamma1DalitzCand || gamma0DalitzCand){
3412  if (gamma0DalitzCand && gamma0MCLabel >=0 && gamma0MCLabel==gamma1MotherLabel){
3413  if (gamma0MotherLabel == -fMesonPDG) isTrueMesonDalitz = kTRUE;
3414  }
3415  if (gamma1DalitzCand && gamma1MCLabel >=0 && gamma1MCLabel==gamma0MotherLabel){
3416  if (gamma1MotherLabel == -fMesonPDG) isTrueMesonDalitz = kTRUE;
3417  }
3418  }
3419 
3420 
3421  if(isTrueMeson ){// True Meosn
3422  if (fHistoTrueMesonInvMassPt[fiCut])fHistoTrueMesonInvMassPt[fiCut]->Fill(mesonCand->M(),mesonCand->Pt(),fWeightJetJetMC);
3423  if (fDoMesonQA > 0){
3424  if ( mesonCand->M() > fMesonInvMassWindow[0] && mesonCand->M() < fMesonInvMassWindow[1]){
3425  if(fIsMC < 2){
3426  if (fHistoTrueMesonPtY[fiCut]) fHistoTrueMesonPtY[fiCut]->Fill(mesonCand->Pt(),mesonCand->Rapidity()-((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetEtaShift());
3427  if (fHistoTrueMesonPtOpenAngle[fiCut]) fHistoTrueMesonPtOpenAngle[fiCut]->Fill(mesonCand->Pt(),mesonCand->GetOpeningAngle());
3428  if (fHistoTrueMesonPtAlpha[fiCut]) fHistoTrueMesonPtAlpha[fiCut]->Fill(mesonCand->Pt(),mesonCand->GetAlpha(),fWeightJetJetMC);
3429  }
3430  }
3431  }
3432  Bool_t isPrimary = ((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsConversionPrimaryESD( fMCEvent, gamma0MotherLabel, mcProdVtxX, mcProdVtxY, mcProdVtxZ);
3433  if(isPrimary){ // Only primary pi0 for efficiency calculation
3434  Float_t weighted= 1;
3435  if(((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsParticleFromBGEvent(gamma1MotherLabel, fMCEvent, fInputEvent)){
3436  if (((TParticle*)fMCEvent->Particle(gamma1MotherLabel))->Pt()>0.005){
3437  weighted= ((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetWeightForMeson(gamma1MotherLabel, fMCEvent, fInputEvent);
3438  // cout << "rec \t " <<gamma1MotherLabel << "\t" << weighted << endl;
3439  }
3440  }
3441  if (fHistoTruePrimaryMesonInvMassPt[fiCut]) fHistoTruePrimaryMesonInvMassPt[fiCut]->Fill(mesonCand->M(),mesonCand->Pt(),weighted*fWeightJetJetMC);
3443  if (fProfileTruePrimaryMesonWeightsInvMassPt[fiCut]) fProfileTruePrimaryMesonWeightsInvMassPt[fiCut]->Fill(mesonCand->M(),mesonCand->Pt(),weighted*fWeightJetJetMC);
3444 
3445  if (fDoMesonQA > 0 && fIsMC < 2){
3446  if (fHistoTruePrimaryMesonMCPtResolPt[fiCut]) fHistoTruePrimaryMesonMCPtResolPt[fiCut]->Fill(((TParticle*)fMCEvent->Particle(gamma1MotherLabel))->Pt(),(mesonCand->Pt()-((TParticle*)fMCEvent->Particle(gamma1MotherLabel))->Pt())/((TParticle*)fMCEvent->Particle(gamma1MotherLabel))->Pt(),weighted);
3447  }
3448  }
3449  } else if(!isTrueMeson ){ // Background
3450  if (fDoMesonQA > 1 && fIsMC < 2){
3451  if(gamma0MotherLabel>-1 && gamma1MotherLabel>-1){ // Both Tracks are Photons and have a mother but not Pi0 or Eta
3452  if (fHistoTrueBckGGInvMassPt[fiCut])fHistoTrueBckGGInvMassPt[fiCut]->Fill(mesonCand->M(),mesonCand->Pt());
3453  } else { // No photon or without mother
3454  if (fHistoTrueBckContInvMassPt[fiCut]) fHistoTrueBckContInvMassPt[fiCut]->Fill(mesonCand->M(),mesonCand->Pt());
3455  }
3456  }
3457  if (!isTrueMesonDalitz && (gamma0DalitzCand || gamma1DalitzCand)){
3458  if (fDoMesonQA > 1 && fIsMC < 2 && fHistoTrueBckContInvMassPt[fiCut]) fHistoTrueBckContInvMassPt[fiCut]->Fill(mesonCand->M(),mesonCand->Pt());
3459  }
3460  }
3461  }
3462  }
3463 }
3464 
3465 //______________________________________________________________________
3467  AliAODConversionPhoton *TrueGammaCandidate0,
3468  AliAODConversionPhoton *TrueGammaCandidate1)
3469 {
3470  const AliVVertex* primVtxMC = fMCEvent->GetPrimaryVertex();
3471  Double_t mcProdVtxX = primVtxMC->GetX();
3472  Double_t mcProdVtxY = primVtxMC->GetY();
3473  Double_t mcProdVtxZ = primVtxMC->GetZ();
3474 
3475  // Process True Mesons
3476  TClonesArray *AODMCTrackArray = dynamic_cast<TClonesArray*>(fInputEvent->FindListObject(AliAODMCParticle::StdBranchName()));
3477  Bool_t isTrueMeson = kFALSE;
3478  Bool_t isTrueMesonDalitz = kFALSE;
3479  Bool_t gamma0DalitzCand = kFALSE;
3480  Bool_t gamma1DalitzCand = kFALSE;
3481 
3482  if (AODMCTrackArray!=NULL && TrueGammaCandidate0 != NULL){
3483  AliAODMCParticle *positiveMC = static_cast<AliAODMCParticle*>(AODMCTrackArray->At(TrueGammaCandidate0->GetMCLabelPositive()));
3484  AliAODMCParticle *negativeMC = static_cast<AliAODMCParticle*>(AODMCTrackArray->At(TrueGammaCandidate0->GetMCLabelNegative()));
3485 
3486  Int_t gamma0MCLabel = -1;
3487  Int_t gamma0MotherLabel = -1;
3488  if(!positiveMC||!negativeMC)
3489  return;
3490 
3491  if(positiveMC->GetMother()>-1&&(negativeMC->GetMother() == positiveMC->GetMother())){
3492  gamma0MCLabel = positiveMC->GetMother();
3493  }
3494 
3495  if(gamma0MCLabel != -1){ // Gamma is Combinatorial; MC Particles don't belong to the same Mother
3496  // Daughters Gamma 0
3497  AliAODMCParticle * gammaMC0 = static_cast<AliAODMCParticle*>(AODMCTrackArray->At(gamma0MCLabel));
3498  if(TMath::Abs(negativeMC->GetPdgCode())==11 && TMath::Abs(positiveMC->GetPdgCode())==11){ // Electrons ...
3499  if(((positiveMC->GetMCProcessCode())) == 5 && ((negativeMC->GetMCProcessCode())) == 5){ // ... From Conversion ...
3500  if(gammaMC0->GetPdgCode() == 22){ // ... with Gamma Mother
3501  gamma0MotherLabel=gammaMC0->GetMother();
3502  }
3503  }
3504  if(gammaMC0->GetPdgCode() ==fMesonPDG){ // Dalitz candidate
3505  gamma0DalitzCand = kTRUE;
3506  gamma0MotherLabel=-fMesonPDG;
3507  }
3508  }
3509  }
3510  positiveMC = static_cast<AliAODMCParticle*>(AODMCTrackArray->At(TrueGammaCandidate1->GetMCLabelPositive()));
3511  negativeMC = static_cast<AliAODMCParticle*>(AODMCTrackArray->At(TrueGammaCandidate1->GetMCLabelNegative()));
3512 
3513  Int_t gamma1MCLabel = -1;
3514  Int_t gamma1MotherLabel = -1;
3515  if(!positiveMC||!negativeMC)
3516  return;
3517 
3518  if(positiveMC->GetMother()>-1&&(negativeMC->GetMother() == positiveMC->GetMother())){
3519  gamma1MCLabel = positiveMC->GetMother();
3520  }
3521  if(gamma1MCLabel != -1){ // Gamma is Combinatorial; MC Particles don't belong to the same Mother
3522  // Daughters Gamma 1
3523  AliAODMCParticle * gammaMC1 = static_cast<AliAODMCParticle*>(AODMCTrackArray->At(gamma1MCLabel));
3524  if(TMath::Abs(negativeMC->GetPdgCode())==11 && TMath::Abs(positiveMC->GetPdgCode())==11){ // Electrons ...
3525  if(((positiveMC->GetMCProcessCode())) == 5 && ((negativeMC->GetMCProcessCode())) == 5){ // ... From Conversion ...
3526  if(gammaMC1->GetPdgCode() == 22){ // ... with Gamma Mother
3527  gamma1MotherLabel=gammaMC1->GetMother();
3528  }
3529  }
3530  if(gammaMC1->GetPdgCode() ==fMesonPDG ){ // Dalitz candidate
3531  gamma1DalitzCand = kTRUE;
3532  gamma1MotherLabel=-fMesonPDG;
3533  }
3534  }
3535  }
3536  if(gamma0MotherLabel>=0 && gamma0MotherLabel==gamma1MotherLabel){
3537  if(static_cast<AliAODMCParticle*>(AODMCTrackArray->At(gamma1MotherLabel))->GetPdgCode() == fMesonPDG){
3538  isTrueMeson=kTRUE;
3542  }
3543  }
3544  }
3545 
3546  //Identify Dalitz candidate
3547  if (gamma1DalitzCand || gamma0DalitzCand){
3548  if (gamma0DalitzCand && gamma0MCLabel >=0 && gamma0MCLabel==gamma1MotherLabel){
3549  if (gamma0MotherLabel == -fMesonPDG) isTrueMesonDalitz = kTRUE;
3550  }
3551  if (gamma1DalitzCand && gamma1MCLabel >=0 && gamma1MCLabel==gamma0MotherLabel){
3552  if (gamma1MotherLabel == -fMesonPDG) isTrueMesonDalitz = kTRUE;
3553  }
3554  }
3555 
3556  if(isTrueMeson ){// True Meson
3557  if (fHistoTrueMesonInvMassPt[fiCut])fHistoTrueMesonInvMassPt[fiCut]->Fill(mesonCand->M(),mesonCand->Pt(),fWeightJetJetMC);
3558  if (fDoMesonQA > 0){
3559  if ( mesonCand->M() > fMesonInvMassWindow[0] && mesonCand->M() < fMesonInvMassWindow[1]){
3560  if(fIsMC < 2){
3561  if (fHistoTrueMesonPtY[fiCut]) fHistoTrueMesonPtY[fiCut]->Fill(mesonCand->Pt(),mesonCand->Rapidity()-((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetEtaShift());
3562  if (fHistoTrueMesonPtOpenAngle[fiCut]) fHistoTrueMesonPtOpenAngle[fiCut]->Fill(mesonCand->Pt(),mesonCand->GetOpeningAngle());
3563  }
3564  if (fHistoTrueMesonPtAlpha[fiCut]) fHistoTrueMesonPtAlpha[fiCut]->Fill(mesonCand->Pt(),mesonCand->GetAlpha(),fWeightJetJetMC);
3565  }
3566 
3567  }
3568  Bool_t isPrimary = ((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsConversionPrimaryAOD(fInputEvent, static_cast<AliAODMCParticle*>(AODMCTrackArray->At(gamma0MotherLabel)), mcProdVtxX, mcProdVtxY, mcProdVtxZ);
3569 
3570  if(isPrimary){ // Only primary pi0 for efficiency calculation
3571  Float_t weighted= 1;
3572  if(((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsParticleFromBGEvent(gamma1MotherLabel, 0x0, fInputEvent)){
3573  if (static_cast<AliAODMCParticle*>(AODMCTrackArray->At(gamma1MotherLabel))->Pt()>0.005){
3574  weighted= ((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetWeightForMeson(gamma1MotherLabel, 0x0, fInputEvent);
3575  // cout << "rec \t " <<gamma1MotherLabel << "\t" << weighted << endl;
3576  }
3577  }
3578  if (fHistoTruePrimaryMesonInvMassPt[fiCut]) fHistoTruePrimaryMesonInvMassPt[fiCut]->Fill(mesonCand->M(),mesonCand->Pt(),weighted*fWeightJetJetMC);
3580  if (fProfileTruePrimaryMesonWeightsInvMassPt[fiCut]) fProfileTruePrimaryMesonWeightsInvMassPt[fiCut]->Fill(mesonCand->M(),mesonCand->Pt(),weighted*fWeightJetJetMC);
3581 
3582  if (fDoMesonQA > 0 && fIsMC < 2){
3583  if (fHistoTruePrimaryMesonMCPtResolPt[fiCut])fHistoTruePrimaryMesonMCPtResolPt[fiCut]->Fill(static_cast<AliAODMCParticle*>(AODMCTrackArray->At(gamma1MotherLabel))->Pt(),
3584  (mesonCand->Pt()-static_cast<AliAODMCParticle*>(AODMCTrackArray->At(gamma1MotherLabel))->Pt())/static_cast<AliAODMCParticle*>(AODMCTrackArray->At(gamma1MotherLabel))->Pt(),weighted);
3585  }
3586  }
3587  } else if(!isTrueMeson ) { // Background
3588  if (fDoMesonQA > 1 && fIsMC < 2){
3589  if(!(gamma0MotherLabel>-1 && gamma1MotherLabel>-1)){ // Both Tracks are Photons and have a mother but not Pi0 or Eta
3590  if (fHistoTrueBckGGInvMassPt[fiCut])fHistoTrueBckGGInvMassPt[fiCut]->Fill(mesonCand->M(),mesonCand->Pt());
3591  } else { // No photon or without mother
3592  if (fHistoTrueBckContInvMassPt[fiCut]) fHistoTrueBckContInvMassPt[fiCut]->Fill(mesonCand->M(),mesonCand->Pt());
3593  }
3594  }
3595  if (!isTrueMesonDalitz && (gamma0DalitzCand || gamma1DalitzCand)){
3596  if (fDoMesonQA > 1 && fIsMC < 2)if (fHistoTrueBckContInvMassPt[fiCut]) fHistoTrueBckContInvMassPt[fiCut]->Fill(mesonCand->M(),mesonCand->Pt());
3597  }
3598  }
3599  }
3600  return;
3601 }
3602 
3603 
3604 //________________________________________________________________________
3606 
3607  // set current BG handler
3608  AliGammaConversionAODBGHandler* currBGHandler = NULL;
3609  TList* currPhotonList = NULL;
3610  if (fMesonRecoMode == 0){
3611  currBGHandler = fBGHandler[fiCut];
3612  currPhotonList = fGammaCandidates;
3613  } else if (fMesonRecoMode == 1){
3614  currBGHandler = fBGClusHandler[fiCut];
3615  currPhotonList = fGammaCandidates;
3616  } else if (fMesonRecoMode == 2){
3617  currBGHandler = fBGClusHandler[fiCut];
3618  currPhotonList = fClusterCandidates;
3619  }
3620 
3621  Int_t zbin = currBGHandler->GetZBinIndex(fInputEvent->GetPrimaryVertex()->GetZ());
3622  Int_t mbin = 0;
3623  if(((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->UseTrackMultiplicity()){
3624  mbin = currBGHandler->GetMultiplicityBinIndex(fV0Reader->GetNumberOfPrimaryTracks());
3625  } else {
3626  mbin = currBGHandler->GetMultiplicityBinIndex(currPhotonList->GetEntries());
3627  }
3628 
3630  if(((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->UseTrackMultiplicity()){
3631  for(Int_t nEventsInBG=0;nEventsInBG<currBGHandler->GetNBGEvents();nEventsInBG++){
3632  AliGammaConversionAODVector *previousEventV0s = currBGHandler->GetBGGoodV0s(zbin,mbin,nEventsInBG);
3633  if (fMesonRecoMode < 2){
3634  if(fMoveParticleAccordingToVertex == kTRUE || ((AliConversionPhotonCuts*)fCutArray->At(fiCut))->GetInPlaneOutOfPlaneCut() != 0 ){
3635  bgEventVertex = currBGHandler->GetBGEventVertex(zbin,mbin,nEventsInBG);
3636  }
3637  }
3638  for(Int_t iCurrent=0;iCurrent<currPhotonList->GetEntries();iCurrent++){
3639  AliAODConversionPhoton currentEventGoodV0 = *(AliAODConversionPhoton*)(currPhotonList->At(iCurrent));
3640  if (!currentEventGoodV0.GetUseForMesonPair() )continue;
3641  for(UInt_t iPrevious=0;iPrevious<previousEventV0s->size();iPrevious++){
3642  AliAODConversionPhoton previousGoodV0 = (AliAODConversionPhoton)(*(previousEventV0s->at(iPrevious)));
3643  if(fMoveParticleAccordingToVertex == kTRUE){
3644  if (bgEventVertex){
3645  MoveParticleAccordingToVertex(&previousGoodV0,bgEventVertex);
3646  }
3647  }
3648  if (fMesonRecoMode < 2){
3649  if(((AliConversionPhotonCuts*)fCutArray->At(fiCut))->GetInPlaneOutOfPlaneCut() != 0){
3650  if (bgEventVertex){
3651  RotateParticleAccordingToEP(&previousGoodV0,bgEventVertex->fEP,fEventPlaneAngle);
3652  }
3653  }
3654  }
3655  AliAODConversionMother *backgroundCandidate = new AliAODConversionMother(&currentEventGoodV0,&previousGoodV0);
3656  backgroundCandidate->CalculateDistanceOfClossetApproachToPrimVtx(fInputEvent->GetPrimaryVertex());
3657  if((((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->MesonIsSelected(backgroundCandidate, kFALSE,
3658  ((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetEtaShift()), currentEventGoodV0.GetLeadingCellID(), previousGoodV0.GetLeadingCellID(), currentEventGoodV0.GetIsCaloPhoton(), previousGoodV0.GetIsCaloPhoton())){
3659  if (fHistoMotherBackInvMassPt[fiCut]) fHistoMotherBackInvMassPt[fiCut]->Fill(backgroundCandidate->M(),backgroundCandidate->Pt(),fWeightJetJetMC);
3660  if(fDoTHnSparse){
3661  Double_t sparesFill[4] = {backgroundCandidate->M(),backgroundCandidate->Pt(),(Double_t)zbin,(Double_t)mbin};
3662  if (fSparseMotherBackInvMassPtZM[fiCut]) fSparseMotherBackInvMassPtZM[fiCut]->Fill(sparesFill,1);
3663  }
3664  }
3665  }
3666  }
3667  }
3668  }else {
3669  for(Int_t nEventsInBG=0;nEventsInBG <currBGHandler->GetNBGEvents();nEventsInBG++){
3670  AliGammaConversionAODVector *previousEventV0s = currBGHandler->GetBGGoodV0s(zbin,mbin,nEventsInBG);
3671  if(previousEventV0s){
3672  if (fMesonRecoMode < 2){
3673  if(fMoveParticleAccordingToVertex == kTRUE || ((AliConversionPhotonCuts*)fCutArray->At(fiCut))->GetInPlaneOutOfPlaneCut() != 0 ){
3674  bgEventVertex = currBGHandler->GetBGEventVertex(zbin,mbin,nEventsInBG);
3675  }
3676  }
3677  for(Int_t iCurrent=0;iCurrent<currPhotonList->GetEntries();iCurrent++){
3678  AliAODConversionPhoton currentEventGoodV0 = *(AliAODConversionPhoton*)(currPhotonList->At(iCurrent));
3679  if (!currentEventGoodV0.GetUseForMesonPair() )continue;
3680  for(UInt_t iPrevious=0;iPrevious<previousEventV0s->size();iPrevious++){
3681  AliAODConversionPhoton previousGoodV0 = (AliAODConversionPhoton)(*(previousEventV0s->at(iPrevious)));
3682  if(fMoveParticleAccordingToVertex == kTRUE){
3683  if (bgEventVertex){
3684  MoveParticleAccordingToVertex(&previousGoodV0,bgEventVertex);
3685  }
3686  }
3687  if (fMesonRecoMode < 2){
3688  if(((AliConversionPhotonCuts*)fCutArray->At(fiCut))->GetInPlaneOutOfPlaneCut() != 0){
3689  if (bgEventVertex){
3690  RotateParticleAccordingToEP(&previousGoodV0,bgEventVertex->fEP,fEventPlaneAngle);
3691  }
3692  }
3693  }
3694  AliAODConversionMother *backgroundCandidate = new AliAODConversionMother(&currentEventGoodV0,&previousGoodV0);
3695  backgroundCandidate->CalculateDistanceOfClossetApproachToPrimVtx(fInputEvent->GetPrimaryVertex());
3696  if((((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->MesonIsSelected(backgroundCandidate,kFALSE,
3697  ((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetEtaShift()), currentEventGoodV0.GetLeadingCellID(), previousGoodV0.GetLeadingCellID(), currentEventGoodV0.GetIsCaloPhoton(), previousGoodV0.GetIsCaloPhoton())){
3698  if (fHistoMotherBackInvMassPt[fiCut]) fHistoMotherBackInvMassPt[fiCut]->Fill(backgroundCandidate->M(),backgroundCandidate->Pt(),fWeightJetJetMC);
3699  if(fDoTHnSparse){
3700  Double_t sparesFill[4] = {backgroundCandidate->M(),backgroundCandidate->Pt(),(Double_t)zbin,(Double_t)mbin};
3701  if (fSparseMotherBackInvMassPtZM[fiCut]) fSparseMotherBackInvMassPtZM[fiCut]->Fill(sparesFill,1);
3702  }
3703  }
3704  delete backgroundCandidate;
3705  backgroundCandidate = 0x0;
3706  }
3707  }
3708  }
3709  }
3710  }
3711 }
3712 
3713 //________________________________________________________________________
3715 
3716  // set current BG handler
3717  AliConversionAODBGHandlerRP* currBGHandler = NULL;
3718  TList* currPhotonList = NULL;
3719  if (fMesonRecoMode == 0){
3720  currBGHandler = fBGHandlerRP[fiCut];
3721  currPhotonList = fGammaCandidates;
3722  } else if (fMesonRecoMode == 1){
3723  currBGHandler = fBGClusHandlerRP[fiCut];
3724  currPhotonList = fGammaCandidates;
3725  } else if (fMesonRecoMode == 2){
3726  currBGHandler = fBGClusHandlerRP[fiCut];
3727  currPhotonList = fClusterCandidates;
3728  }
3729 
3730  Int_t zbin = currBGHandler->GetZBinIndex(fInputEvent->GetPrimaryVertex()->GetZ());
3731  Int_t mbin = 0;
3732  if(((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->UseTrackMultiplicity()){
3733  mbin = currBGHandler->GetMultiplicityBinIndex(fV0Reader->GetNumberOfPrimaryTracks());
3734  } else {
3735  mbin = currBGHandler->GetMultiplicityBinIndex(currPhotonList->GetEntries());
3736  }
3737 
3738  //Rotation Method
3739  if(((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->UseRotationMethod()){
3740  // Correct for the number of rotations
3741  // BG is for rotation the same, except for factor NRotations
3742  Double_t weight=1./Double_t(((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->GetNumberOfBGEvents());
3743  for(Int_t firstGammaIndex=0;firstGammaIndex<currPhotonList->GetEntries();firstGammaIndex++){
3744  AliAODConversionPhoton *gamma0=dynamic_cast<AliAODConversionPhoton*>(currPhotonList->At(firstGammaIndex));
3745  if (gamma0==NULL) continue;
3746  if (!gamma0->GetUseForMesonPair() ) continue;
3747  for(Int_t secondGammaIndex=firstGammaIndex+1;secondGammaIndex<currPhotonList->GetEntries();secondGammaIndex++){
3748  AliAODConversionPhoton *gamma1=dynamic_cast<AliAODConversionPhoton*>(currPhotonList->At(secondGammaIndex));
3749  if (gamma1 == NULL) continue;
3750  if (!gamma1->GetUseForMesonPair() ) continue;
3751  if(!((AliConversionPhotonCuts*)fCutArray->At(fiCut))->PhotonIsSelected(gamma1,fInputEvent))continue;
3752  for(Int_t nRandom=0;nRandom<((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->GetNumberOfBGEvents();nRandom++){
3753  RotateParticle(gamma1);
3754  AliAODConversionMother backgroundCandidate(gamma0,gamma1);
3755  backgroundCandidate.CalculateDistanceOfClossetApproachToPrimVtx(fInputEvent->GetPrimaryVertex());
3756  if(((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->MesonIsSelected(&backgroundCandidate, kFALSE, ((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetEtaShift(), gamma0->GetLeadingCellID(), gamma1->GetLeadingCellID(), gamma0->GetIsCaloPhoton(), gamma1->GetIsCaloPhoton() )){
3757 
3758  if (fHistoMotherBackInvMassPt[fiCut]) fHistoMotherBackInvMassPt[fiCut]->Fill(backgroundCandidate.M(),backgroundCandidate.Pt(),fWeightJetJetMC);
3759  if(fDoTHnSparse){
3760  Double_t sparesFill[4] = {backgroundCandidate.M(),backgroundCandidate.Pt(),(Double_t)zbin,(Double_t)mbin};
3761  if (fSparseMotherBackInvMassPtZM[fiCut]) fSparseMotherBackInvMassPtZM[fiCut]->Fill(sparesFill,weight);
3762  }
3763  }
3764  }
3765  }
3766  }
3767 
3768  }else {
3769  // Do Event Mixing
3770  for(Int_t nEventsInBG=0;nEventsInBG <currBGHandler->GetNBGEvents(currPhotonList,fInputEvent);nEventsInBG++){
3771  AliGammaConversionPhotonVector *previousEventGammas = currBGHandler->GetBGGoodGammas(currPhotonList,fInputEvent,nEventsInBG);
3772  if(previousEventGammas){
3773  // test weighted background
3774  Double_t weight=1.0;
3775  // Correct for the number of eventmixing:
3776  // N gammas -> (N-1) + (N-2) +(N-3) ...+ (N-(N-1)) using sum formula sum(i)=N*(N-1)/2 -> N*(N-1)/2
3777  // real combinations (since you cannot combine a photon with its own)
3778  // but BG leads to N_{a}*N_{b}combinations
3779  weight*=0.5*(Double_t(currPhotonList->GetEntries()-1))/Double_t(previousEventGammas->size());
3780 
3781  for(Int_t iCurrent=0;iCurrent<currPhotonList->GetEntries();iCurrent++){
3782  AliAODConversionPhoton *gamma0 = (AliAODConversionPhoton*)(currPhotonList->At(iCurrent));
3783  if (!gamma0->GetUseForMesonPair() ) continue;
3784  for(UInt_t iPrevious=0;iPrevious<previousEventGammas->size();iPrevious++){
3785  AliAODConversionPhoton *gamma1 = (AliAODConversionPhoton*)(previousEventGammas->at(iPrevious));
3786  AliAODConversionMother backgroundCandidate(gamma0,gamma1);
3787  backgroundCandidate.CalculateDistanceOfClossetApproachToPrimVtx(fInputEvent->GetPrimaryVertex());
3788  if(((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))
3789  ->MesonIsSelected(&backgroundCandidate,kFALSE,((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetEtaShift(), gamma0->GetLeadingCellID(), gamma1->GetLeadingCellID(), gamma0->GetIsCaloPhoton(), gamma1->GetIsCaloPhoton())){
3790  if (fHistoMotherBackInvMassPt[fiCut]) fHistoMotherBackInvMassPt[fiCut]->Fill(backgroundCandidate.M(),backgroundCandidate.Pt(),fWeightJetJetMC);
3791  if(fDoTHnSparse){
3792  Double_t sparesFill[4] = {backgroundCandidate.M(),backgroundCandidate.Pt(),(Double_t)zbin,(Double_t)mbin};
3793  if (fSparseMotherBackInvMassPtZM[fiCut]) fSparseMotherBackInvMassPtZM[fiCut]->Fill(sparesFill,weight);
3794  }
3795  }
3796  }
3797  }
3798  }
3799  }
3800  }
3801 }
3802 
3803 //________________________________________________________________________
3804 void AliAnalysisTaskHeavyNeutralMesonToGG::RotateParticle(AliAODConversionPhoton *gamma){
3805  Int_t fNDegreesPMBackground= ((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->NDegreesRotation();
3806  Double_t nRadiansPM = fNDegreesPMBackground*TMath::Pi()/180;
3807  Double_t rotationValue = fRandom.Rndm()*2*nRadiansPM + TMath::Pi()-nRadiansPM;
3808  gamma->RotateZ(rotationValue);
3809 }
3810 
3811 //________________________________________________________________________
3812 void AliAnalysisTaskHeavyNeutralMesonToGG::RotateParticleAccordingToEP(AliAODConversionPhoton *gamma, Double_t previousEventEP, Double_t thisEventEP){
3813  previousEventEP=previousEventEP+TMath::Pi();
3814  thisEventEP=thisEventEP+TMath::Pi();
3815  Double_t rotationValue= thisEventEP-previousEventEP;
3816  gamma->RotateZ(rotationValue);
3817 }
3818 
3819 //________________________________________________________________________
3821  //see header file for documentation
3822 
3823  Double_t dx = vertex->fX - fInputEvent->GetPrimaryVertex()->GetX();
3824  Double_t dy = vertex->fY - fInputEvent->GetPrimaryVertex()->GetY();
3825  Double_t dz = vertex->fZ - fInputEvent->GetPrimaryVertex()->GetZ();
3826 
3827  Double_t movedPlace[3] = {particle->GetConversionX() - dx,particle->GetConversionY() - dy,particle->GetConversionZ() - dz};
3828  particle->SetConversionPoint(movedPlace);
3829 }
3830 //________________________________________________________________________
3832  //see header file for documentation
3833  if(fGammaCandidates->GetEntries() >1 && fMesonRecoMode == 0 ){
3834  if(((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->UseTrackMultiplicity()){
3835  fBGHandler[fiCut]->AddEvent(fGammaCandidates,fInputEvent->GetPrimaryVertex()->GetX(), fInputEvent->GetPrimaryVertex()->GetY(), fInputEvent->GetPrimaryVertex()->GetZ(),
3836  fV0Reader->GetNumberOfPrimaryTracks(), fEventPlaneAngle);
3837  }else { // means we use #V0s for multiplicity
3838  fBGHandler[fiCut]->AddEvent(fGammaCandidates, fInputEvent->GetPrimaryVertex()->GetX(), fInputEvent->GetPrimaryVertex()->GetY(), fInputEvent->GetPrimaryVertex()->GetZ(),
3839  fGammaCandidates->GetEntries(), fEventPlaneAngle);
3840  }
3841  } else if((fGammaCandidates->GetEntries() >0 || fClusterCandidates->GetEntries() > 0 ) && fMesonRecoMode == 1 ){
3842  if(((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->UseTrackMultiplicity()){
3843  fBGHandler[fiCut]->AddEvent(fGammaCandidates,fInputEvent->GetPrimaryVertex()->GetX(), fInputEvent->GetPrimaryVertex()->GetY(), fInputEvent->GetPrimaryVertex()->GetZ(),
3844  fV0Reader->GetNumberOfPrimaryTracks(), fEventPlaneAngle);
3845  fBGClusHandler[fiCut]->AddEvent(fClusterCandidates,fInputEvent->GetPrimaryVertex()->GetX(), fInputEvent->GetPrimaryVertex()->GetY(), fInputEvent->GetPrimaryVertex()->GetZ(),
3846  fV0Reader->GetNumberOfPrimaryTracks(), fEventPlaneAngle);
3847  }else { // means we use #V0s for multiplicity
3848  fBGHandler[fiCut]->AddEvent(fGammaCandidates, fInputEvent->GetPrimaryVertex()->GetX(), fInputEvent->GetPrimaryVertex()->GetY(), fInputEvent->GetPrimaryVertex()->GetZ(),
3849  fGammaCandidates->GetEntries(), fEventPlaneAngle);
3850  fBGClusHandler[fiCut]->AddEvent(fClusterCandidates, fInputEvent->GetPrimaryVertex()->GetX(), fInputEvent->GetPrimaryVertex()->GetY(), fInputEvent->GetPrimaryVertex()->GetZ(),
3851  fGammaCandidates->GetEntries(), fEventPlaneAngle);
3852  }
3853  } else if(fClusterCandidates->GetEntries() >1 && fMesonRecoMode == 2 ){
3854  if(((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->UseTrackMultiplicity()){
3855  fBGClusHandler[fiCut]->AddEvent(fClusterCandidates,fInputEvent->GetPrimaryVertex()->GetX(), fInputEvent->GetPrimaryVertex()->GetY(), fInputEvent->GetPrimaryVertex()->GetZ(),
3856  fV0Reader->GetNumberOfPrimaryTracks(), fEventPlaneAngle);
3857  }else { // means we use #V0s for multiplicity
3858  fBGClusHandler[fiCut]->AddEvent(fClusterCandidates, fInputEvent->GetPrimaryVertex()->GetX(), fInputEvent->GetPrimaryVertex()->GetY(), fInputEvent->GetPrimaryVertex()->GetZ(),
3859  fClusterCandidates->GetEntries(), fEventPlaneAngle);
3860  }
3861  }
3862 }
3863 
3864 //________________________________________________________________________
3866 
3867  // Relabeling For AOD Event
3868  // ESDiD -> AODiD
3869  // MCLabel -> AODMCLabel
3870 
3871  if(mode){
3872  fMCEventPos = new Int_t[fReaderGammas->GetEntries()];
3873  fMCEventNeg = new Int_t[fReaderGammas->GetEntries()];
3874  fESDArrayPos = new Int_t[fReaderGammas->GetEntries()];
3875  fESDArrayNeg = new Int_t[fReaderGammas->GetEntries()];
3876  }
3877 
3878  for(Int_t iGamma = 0;iGamma<fReaderGammas->GetEntries();iGamma++){
3879  AliAODConversionPhoton* PhotonCandidate = (AliAODConversionPhoton*) fReaderGammas->At(iGamma);
3880  if(!PhotonCandidate) continue;
3881  if(!mode){// Back to ESD Labels
3882  PhotonCandidate->SetMCLabelPositive(fMCEventPos[iGamma]);
3883  PhotonCandidate->SetMCLabelNegative(fMCEventNeg[iGamma]);
3884  PhotonCandidate->SetLabelPositive(fESDArrayPos[iGamma]);
3885  PhotonCandidate->SetLabelNegative(fESDArrayNeg[iGamma]);
3886  continue;
3887  }
3888  fMCEventPos[iGamma] = PhotonCandidate->GetMCLabelPositive();
3889  fMCEventNeg[iGamma] = PhotonCandidate->GetMCLabelNegative();
3890  fESDArrayPos[iGamma] = PhotonCandidate->GetTrackLabelPositive();
3891  fESDArrayNeg[iGamma] = PhotonCandidate->GetTrackLabelNegative();
3892 
3893  Bool_t AODLabelPos = kFALSE;
3894  Bool_t AODLabelNeg = kFALSE;
3895 
3896  for(Int_t i = 0; i<fInputEvent->GetNumberOfTracks();i++){
3897  AliAODTrack *tempDaughter = static_cast<AliAODTrack*>(fInputEvent->GetTrack(i));
3898  if(!AODLabelPos){
3899  if( tempDaughter->GetID() == PhotonCandidate->GetTrackLabelPositive() ){
3900  PhotonCandidate->SetMCLabelPositive(TMath::Abs(tempDaughter->GetLabel()));
3901  PhotonCandidate->SetLabelPositive(i);
3902  AODLabelPos = kTRUE;
3903  }
3904  }
3905  if(!AODLabelNeg){
3906  if( tempDaughter->GetID() == PhotonCandidate->GetTrackLabelNegative()){
3907  PhotonCandidate->SetMCLabelNegative(TMath::Abs(tempDaughter->GetLabel()));
3908  PhotonCandidate->SetLabelNegative(i);
3909  AODLabelNeg = kTRUE;
3910  }
3911  }
3912  if(AODLabelNeg && AODLabelPos){
3913  break;
3914  }
3915  }
3916  if(!AODLabelPos || !AODLabelNeg){
3917  cout<<"WARNING!!! AOD TRACKS NOT FOUND FOR"<<endl;
3918  if(!AODLabelNeg){
3919  PhotonCandidate->SetMCLabelNegative(-999999);
3920  PhotonCandidate->SetLabelNegative(-999999);
3921  }
3922  if(!AODLabelPos){
3923  PhotonCandidate->SetMCLabelPositive(-999999);
3924  PhotonCandidate->SetLabelPositive(-999999);
3925  }
3926  }
3927  }
3928 
3929  if(!mode){
3930  delete[] fMCEventPos;
3931  delete[] fMCEventNeg;
3932  delete[] fESDArrayPos;
3933  delete[] fESDArrayNeg;
3934  }
3935 }
3936 
3937 //________________________________________________________________________
3939  TAxis *axisafter = histoRebin->GetXaxis();
3940  Int_t bins = axisafter->GetNbins();
3941  Double_t from = axisafter->GetXmin();
3942  Double_t to = axisafter->GetXmax();
3943  Double_t *newbins = new Double_t[bins+1];
3944  newbins[0] = from;
3945  Double_t factor = TMath::Power(to/from, 1./bins);
3946  for(Int_t i=1; i<=bins; ++i) newbins[i] = factor * newbins[i-1];
3947  axisafter->Set(bins, newbins);
3948  delete [] newbins;
3949 }
3950 
3951 //________________________________________________________________________
3953 {
3954  //fOutputContainer->Print(); // Will crash on GRID
3955 }
3956 
3957 //_________________________________________________________________________________
3959  if(tobechecked > -1){
3960  vector<Int_t>::iterator it;
3961  it = find (vec.begin(), vec.end(), tobechecked);
3962  if (it != vec.end()) return true;
3963  else return false;
3964  }
3965  return false;
3966 }
3967 
3968 //_________________________________________________________________________________
3970  if(tobechecked > -1){
3971  vector<Int_t>::iterator it;
3972  it = find (vec.begin(), vec.end(), tobechecked);
3973  if (it != vec.end()) return true;
3974  else{
3975  vec.push_back(tobechecked);
3976  return false;
3977  }
3978  }
3979  return false;
3980 }
3981 
3982 //_________________________________________________________________________________
3984  if( ma.find(tobechecked) != ma.end() ) ma[tobechecked] += 1;
3985  else ma[tobechecked] = 2;
3986  return;
3987 }
3988 
3989 //_________________________________________________________________________________
3991  map<Int_t, Int_t>::iterator it;
3992  for (it = ma.begin(); it != ma.end(); it++){
3993  hist->Fill(it->second, fWeightJetJetMC);
3994  }
3995  ma.clear();
3996  return;
3997 }
TH1F ** fHistoClusOverlapHeadersGammaPt
array of histos with cluster, E
TH1F ** fHistoNGammaConvCandidates
array of histos with vertex y distribution for selected events
TH2F ** fHistoSPDClusterTrackletBackground
array of histos with double counted cluster photons
TList * fClusterCandidates
current list of photon candidates
TH1F ** fHistoVertexX
array of histos with vertex z distribution for selected events
TTree * fTreeBrokenFiles
array of THnSparseF with BG for same event photon pairs, inv Mass, pt
TH2F ** fHistoTrueMesonCaloMergedClusterPartConvInvMassPt
array of histos with validated mothers, merged cluster invMass, pt
double Double_t
Definition: External.C:58
Bool_t CheckVectorForDoubleCount(vector< Int_t > &vec, Int_t tobechecked)
GammaConversionVertex * GetBGEventVertex(Int_t zbin, Int_t mbin, Int_t event)
TH1F ** fHistoMCMesonWOWeightInAccPt
array of histos with weighted meson in acceptance, pT
Definition: External.C:236
TList * fEventCutArray
current list of cluster candidates
TH1F ** fHistoClusGammaE
array of histos with cluster, pt
TH2F ** fHistoTrueMotherMesonConvPhotonEtaPhi
array of histos with validated weighted primary meson, MCpt, resol pt
TH2F ** fHistoTrueBckFullMesonContainedInOneClusterInvMassPt
array of histos with pure gamma gamma combinatorial BG, invMass, pt
TH2F ** fHistoTrueMesonPtOpenAngle
array of histos with validated meson, pt, alpha
TH1F ** fHistoTrueClusGammaPt
array of histos with validated primary
TH1F ** fHistoMultipleCountTrueConvGamma
array of histos how often TrueMesons are counted
TH1F ** fHistoMCMesonPt
array of histos with invariant mass pairs which were rejected
TH2F ** fHistoTrueMesonCaloConvertedPhotonMatchedInvMassPt
array of histos with validated meson converted photon leading and photon, invMass, pt
TH2F ** fHistoMotherMesonPtY
array of histogram with BG for mixed event photon pairs, inv Mass, pt
TH2F ** fHistoMCMesonPtY
array of histos with validated meson, pt, openAngle