]> git.uio.no Git - u/mrichter/AliRoot.git/blame_incremental - PWGGA/GammaConv/AliAnalysisTaskNeutralMesonToPiPlPiMiPiZero.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame_incremental
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (non-incremental) | history | HEAD
Added LHC11h_pass2 calibration
[u/mrichter/AliRoot.git] / PWGGA / GammaConv / AliAnalysisTaskNeutralMesonToPiPlPiMiPiZero.cxx
... / ...
CommitLineData
1/**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: Friederike Bock *
5 * Version 1 *
6 * *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
15
16//
17
18#include <vector>
19#include "TParticle.h"
20#include "TPDGCode.h"
21#include "TMCProcess.h"
22#include "TDatabasePDG.h"
23#include "TList.h"
24#include "TChain.h"
25#include "TDirectory.h"
26#include "TTree.h"
27#include "TH1.h"
28#include "TH1F.h"
29#include "THnSparse.h"
30#include "TH2F.h"
31#include "AliStack.h"
32#include "AliAnalysisManager.h"
33#include "AliESDInputHandler.h"
34#include "AliESDtrack.h"
35#include "AliMCEvent.h"
36#include "AliStack.h"
37#include "AliMCEventHandler.h"
38#include "AliPID.h"
39#include "AliLog.h"
40#include "AliESDtrackCuts.h"
41#include "AliESDpidCuts.h"
42#include "AliMCEvent.h"
43#include "AliESDv0.h"
44#include "AliESDEvent.h"
45#include "AliESDpid.h"
46#include "AliKFParticle.h"
47#include "AliMCEventHandler.h"
48#include "AliKFVertex.h"
49#include "AliTriggerAnalysis.h"
50#include "AliCentrality.h"
51#include "AliMultiplicity.h"
52#include "AliAnalysisTaskNeutralMesonToPiPlPiMiPiZero.h"
53
54
55ClassImp( AliAnalysisTaskNeutralMesonToPiPlPiMiPiZero )
56
57//-----------------------------------------------------------------------------------------------
58AliAnalysisTaskNeutralMesonToPiPlPiMiPiZero::AliAnalysisTaskNeutralMesonToPiPlPiMiPiZero():
59 fV0Reader(NULL),
60 fPionSelector(NULL),
61 fBGHandler(NULL),
62 fESDEvent(NULL),
63 fMCEvent(NULL),
64 fMCStack(NULL),
65 fCutFolder(NULL),
66 fESDList(NULL),
67 fBackList(NULL),
68 fMotherList(NULL),
69 fTrueList(NULL),
70 fMCList(NULL),
71 fOutputContainer(0),
72 fReaderGammas(NULL),
73 fSelectorNegPionIndex(0),
74 fSelectorPosPionIndex(0),
75 fGoodConvGammas(NULL),
76 fClusterCandidates(NULL),
77 fNeutralPionCandidates(NULL),
78 fGoodVirtualParticles(NULL),
79 fEventCutArray(NULL),
80 fGammaCutArray(NULL),
81 fClusterCutArray(NULL),
82 fPionCutArray(NULL),
83 fNeutralPionMesonCutArray(NULL),
84 fMesonCutArray(NULL),
85 fEventCuts(NULL),
86 fConversionCuts(NULL),
87 fClusterCuts(NULL),
88 fHistoConvGammaPt(NULL),
89 fHistoConvGammaEta(NULL),
90 fHistoClusterGammaPt(NULL),
91 fHistoClusterGammaEta(NULL),
92 fHistoNegPionPt(NULL),
93 fHistoPosPionPt(NULL),
94 fHistoNegPionPhi(NULL),
95 fHistoPosPionPhi(NULL),
96 fHistoNegPionEta(NULL),
97 fHistoPosPionEta(NULL),
98 fHistoNegPionClsTPC(NULL),
99 fHistoPosPionClsTPC(NULL),
100 fHistoPionDCAxy(NULL),
101 fHistoPionDCAz(NULL),
102 fHistoPionTPCdEdxNSigma(NULL),
103 fHistoPionTPCdEdx(NULL),
104 fHistoPionPionInvMassPt(NULL),
105 fHistoGammaGammaInvMassPt(NULL),
106 fHistoMotherInvMassPt(NULL),
107 fTHnSparseMotherInvMassPtZM(NULL),
108 fHistoMotherBackInvMassPt(NULL),
109 fTHnSparseMotherBackInvMassPtZM(NULL),
110 fHistoMCAllGammaPt(NULL),
111 fHistoMCConvGammaPt(NULL),
112 fHistoMCAllPosPionsPt(NULL),
113 fHistoMCAllNegPionsPt(NULL),
114 fHistoMCGammaFromNeutralMesonPt(NULL),
115 fHistoMCPosPionsFromNeutralMesonPt(NULL),
116 fHistoMCNegPionsFromNeutralMesonPt(NULL),
117 fHistoMCEtaPiPlPiMiPiZeroPt(NULL),
118 fHistoMCEtaPiPlPiMiPiZeroInAccPt(NULL),
119 fHistoMCOmegaPiPlPiMiPiZeroPt(NULL),
120 fHistoMCOmegaPiPlPiMiPiZeroInAccPt(NULL),
121 fHistoTrueMotherPiPlPiMiPiZeroInvMassPt(NULL),
122 fHistoTrueMotherGammaGammaInvMassPt(NULL),
123 fHistoTrueMotherGammaGammaFromEtaInvMassPt(NULL),
124 fHistoTrueMotherGammaGammaFromOmegaInvMassPt(NULL),
125 fHistoTrueConvGammaPt(NULL),
126 fHistoTrueConvGammaFromNeutralMesonPt(NULL),
127 fHistoTrueClusterGammaPt(NULL),
128 fHistoTrueClusterGammaFromNeutralMesonPt(NULL),
129 fHistoTruePosPionPt(NULL),
130 fHistoTruePosPionFromNeutralMesonPt(NULL),
131 fHistoTrueNegPionPt(NULL),
132 fHistoTrueNegPionFromNeutralMesonPt(NULL),
133 fHistoTruePionPionInvMassPt(NULL),
134 fHistoTruePionPionFromSameMotherInvMassPt(NULL),
135 fHistoTruePionPionFromEtaInvMassPt(NULL),
136 fHistoTruePionPionFromOmegaInvMassPt(NULL),
137 fHistoNEvents(NULL),
138 fHistoNGoodESDTracks(NULL),
139 fProfileEtaShift(NULL),
140 fRandom(0),
141 fnCuts(0),
142 fiCut(0),
143 fNumberOfESDTracks(0),
144 fMoveParticleAccordingToVertex(kFALSE),
145 fIsHeavyIon(0),
146 fDoMesonAnalysis(kTRUE),
147 fDoMesonQA(kFALSE),
148 fIsFromMBHeader(kTRUE),
149 fIsMC(kFALSE),
150 fNeutralPionMode(0)
151{
152
153}
154
155//-----------------------------------------------------------------------------------------------
156AliAnalysisTaskNeutralMesonToPiPlPiMiPiZero::AliAnalysisTaskNeutralMesonToPiPlPiMiPiZero( const char* name ):
157 AliAnalysisTaskSE(name),
158 fV0Reader(NULL),
159 fPionSelector(NULL),
160 fBGHandler(NULL),
161 fESDEvent(NULL),
162 fMCEvent(NULL),
163 fMCStack(NULL),
164 fCutFolder(NULL),
165 fESDList(NULL),
166 fBackList(NULL),
167 fMotherList(NULL),
168 fTrueList(NULL),
169 fMCList(NULL),
170 fOutputContainer(0),
171 fReaderGammas(NULL),
172 fSelectorNegPionIndex(0),
173 fSelectorPosPionIndex(0),
174 fGoodConvGammas(NULL),
175 fClusterCandidates(NULL),
176 fNeutralPionCandidates(NULL),
177 fGoodVirtualParticles(NULL),
178 fEventCutArray(NULL),
179 fGammaCutArray(NULL),
180 fClusterCutArray(NULL),
181 fPionCutArray(NULL),
182 fNeutralPionMesonCutArray(NULL),
183 fMesonCutArray(NULL),
184 fEventCuts(NULL),
185 fConversionCuts(NULL),
186 fClusterCuts(NULL),
187 fHistoConvGammaPt(NULL),
188 fHistoConvGammaEta(NULL),
189 fHistoClusterGammaPt(NULL),
190 fHistoClusterGammaEta(NULL),
191 fHistoNegPionPt(NULL),
192 fHistoPosPionPt(NULL),
193 fHistoNegPionPhi(NULL),
194 fHistoPosPionPhi(NULL),
195 fHistoNegPionEta(NULL),
196 fHistoPosPionEta(NULL),
197 fHistoNegPionClsTPC(NULL),
198 fHistoPosPionClsTPC(NULL),
199 fHistoPionDCAxy(NULL),
200 fHistoPionDCAz(NULL),
201 fHistoPionTPCdEdxNSigma(NULL),
202 fHistoPionTPCdEdx(NULL),
203 fHistoPionPionInvMassPt(NULL),
204 fHistoGammaGammaInvMassPt(NULL),
205 fHistoMotherInvMassPt(NULL),
206 fTHnSparseMotherInvMassPtZM(NULL),
207 fHistoMotherBackInvMassPt(NULL),
208 fTHnSparseMotherBackInvMassPtZM(NULL),
209 fHistoMCAllGammaPt(NULL),
210 fHistoMCConvGammaPt(NULL),
211 fHistoMCAllPosPionsPt(NULL),
212 fHistoMCAllNegPionsPt(NULL),
213 fHistoMCGammaFromNeutralMesonPt(NULL),
214 fHistoMCPosPionsFromNeutralMesonPt(NULL),
215 fHistoMCNegPionsFromNeutralMesonPt(NULL),
216 fHistoMCEtaPiPlPiMiPiZeroPt(NULL),
217 fHistoMCEtaPiPlPiMiPiZeroInAccPt(NULL),
218 fHistoMCOmegaPiPlPiMiPiZeroPt(NULL),
219 fHistoMCOmegaPiPlPiMiPiZeroInAccPt(NULL),
220 fHistoTrueMotherPiPlPiMiPiZeroInvMassPt(NULL),
221 fHistoTrueMotherGammaGammaInvMassPt(NULL),
222 fHistoTrueMotherGammaGammaFromEtaInvMassPt(NULL),
223 fHistoTrueMotherGammaGammaFromOmegaInvMassPt(NULL),
224 fHistoTrueConvGammaPt(NULL),
225 fHistoTrueConvGammaFromNeutralMesonPt(NULL),
226 fHistoTrueClusterGammaPt(NULL),
227 fHistoTrueClusterGammaFromNeutralMesonPt(NULL),
228 fHistoTruePosPionPt(NULL),
229 fHistoTruePosPionFromNeutralMesonPt(NULL),
230 fHistoTrueNegPionPt(NULL),
231 fHistoTrueNegPionFromNeutralMesonPt(NULL),
232 fHistoTruePionPionInvMassPt(NULL),
233 fHistoTruePionPionFromSameMotherInvMassPt(NULL),
234 fHistoTruePionPionFromEtaInvMassPt(NULL),
235 fHistoTruePionPionFromOmegaInvMassPt(NULL),
236 fHistoNEvents(NULL),
237 fHistoNGoodESDTracks(NULL),
238 fProfileEtaShift(NULL),
239 fRandom(0),
240 fnCuts(0),
241 fiCut(0),
242 fNumberOfESDTracks(0),
243 fMoveParticleAccordingToVertex(kFALSE),
244 fIsHeavyIon(0),
245 fDoMesonAnalysis(kTRUE),
246 fDoMesonQA(kFALSE),
247 fIsFromMBHeader(kTRUE),
248 fIsMC(kFALSE),
249 fNeutralPionMode(0)
250{
251 DefineOutput(1, TList::Class());
252}
253
254//-----------------------------------------------------------------------------------------------
255AliAnalysisTaskNeutralMesonToPiPlPiMiPiZero::~AliAnalysisTaskNeutralMesonToPiPlPiMiPiZero()
256{
257 //
258 // virtual destructor
259 //
260 cout<<"Destructor"<<endl;
261 if(fGoodConvGammas){
262 delete fGoodConvGammas;
263 fGoodConvGammas = 0x0;
264 }
265 if(fClusterCandidates){
266 delete fClusterCandidates;
267 fClusterCandidates = 0x0;
268 }
269
270 if(fNeutralPionCandidates){
271 delete fNeutralPionCandidates;
272 fNeutralPionCandidates = 0x0;
273 }
274
275 if(fGoodVirtualParticles){
276 delete fGoodVirtualParticles;
277 fGoodVirtualParticles = 0x0;
278 }
279 if(fBGHandler){
280 delete[] fBGHandler;
281 fBGHandler = 0x0;
282 }
283}
284//___________________________________________________________
285void AliAnalysisTaskNeutralMesonToPiPlPiMiPiZero::InitBack(){
286
287 const Int_t nDim = 4;
288 Int_t nBins[nDim] = {500,250,7,4};
289 Double_t xMin[nDim] = {0.4,0, 0,0};
290 Double_t xMax[nDim] = {0.9,25,7,4};
291
292 fTHnSparseMotherInvMassPtZM = new THnSparseF*[fnCuts];
293 fTHnSparseMotherBackInvMassPtZM = new THnSparseF*[fnCuts];
294
295 fBGHandler = new AliGammaConversionAODBGHandler*[fnCuts];
296
297 for(Int_t iCut = 0; iCut<fnCuts;iCut++){
298
299 TString cutstringEvent = ((AliConvEventCuts*)fEventCutArray->At(iCut))->GetCutNumber();
300 TString cutstringPion = ((AliPrimaryPionCuts*)fPionCutArray->At(iCut))->GetCutNumber();
301 TString cutstringConvGamma = "";
302 if (fNeutralPionMode < 2) cutstringConvGamma = ((AliConversionPhotonCuts*)fGammaCutArray->At(iCut))->GetCutNumber();
303 TString cutstringCaloGamma = "";
304 if (fNeutralPionMode > 0) cutstringCaloGamma = ((AliCaloPhotonCuts*)fClusterCutArray->At(iCut))->GetCutNumber();
305 TString cutstringNeutralPion= ((AliConversionMesonCuts*)fNeutralPionMesonCutArray->At(iCut))->GetCutNumber();
306 TString cutstringMeson = ((AliConversionMesonCuts*)fMesonCutArray->At(iCut))->GetCutNumber();
307
308 TString fullCutString = "";
309 if (fNeutralPionMode == 0) fullCutString = Form("%i_%s_%s_%s_%s_%s",fNeutralPionMode,cutstringEvent.Data(),cutstringConvGamma.Data(),cutstringNeutralPion.Data(), cutstringPion.Data(),cutstringMeson.Data());
310 else if (fNeutralPionMode == 1) fullCutString = Form("%i_%s_%s_%s_%s_%s_%s",fNeutralPionMode,cutstringEvent.Data(),cutstringConvGamma.Data(),cutstringCaloGamma.Data(),cutstringNeutralPion.Data(), cutstringPion.Data(),cutstringMeson.Data());
311 else if (fNeutralPionMode == 2) fullCutString = Form("%i_%s_%s_%s_%s_%s",fNeutralPionMode,cutstringEvent.Data(),cutstringCaloGamma.Data(),cutstringNeutralPion.Data(), cutstringPion.Data(),cutstringMeson.Data());
312
313 TString nameBackList = Form("%s Back histograms",fullCutString.Data());
314 TString nameMotherList = Form("%s Mother histograms",fullCutString.Data());
315
316 Int_t collisionSystem = atoi((TString)(((AliConvEventCuts*)fEventCutArray->At(iCut))->GetCutNumber())(0,1));
317 Int_t centMin = atoi((TString)(((AliConvEventCuts*)fEventCutArray->At(iCut))->GetCutNumber())(1,1));
318 Int_t centMax = atoi((TString)(((AliConvEventCuts*)fEventCutArray->At(iCut))->GetCutNumber())(2,1));
319
320 if(collisionSystem == 1 || collisionSystem == 2 ||
321 collisionSystem == 5 || collisionSystem == 8 ||
322 collisionSystem == 9){
323 centMin = centMin*10;
324 centMax = centMax*10;
325 }
326 else if(collisionSystem == 3 || collisionSystem == 6){
327 centMin = centMin*5;
328 centMax = centMax*5;
329 }
330 else if(collisionSystem == 4 || collisionSystem == 7){
331 centMin = ((centMin*5)+45);
332 centMax = ((centMax*5)+45);
333 }
334
335
336 fBackList[iCut] = new TList();
337
338 fBackList[iCut]->SetName(nameBackList.Data());
339 fBackList[iCut]->SetOwner(kTRUE);
340 fCutFolder[iCut]->Add(fBackList[iCut]);
341
342 fTHnSparseMotherBackInvMassPtZM[iCut] = new THnSparseF("Back_Back_InvMass_Pt_z_m","Back_Back_InvMass_Pt_z_m",nDim,nBins,xMin,xMax);
343 fBackList[iCut]->Add(fTHnSparseMotherBackInvMassPtZM[iCut]);
344
345 fMotherList[iCut] = new TList();
346 fMotherList[iCut]->SetName(nameMotherList.Data());
347 fMotherList[iCut]->SetOwner(kTRUE);
348 fCutFolder[iCut]->Add(fMotherList[iCut]);
349
350 fTHnSparseMotherInvMassPtZM[iCut] = new THnSparseF("Back_Mother_InvMass_Pt_z_m","Back_Mother_InvMass_Pt_z_m",nDim,nBins,xMin,xMax);
351 fMotherList[iCut]->Add(fTHnSparseMotherInvMassPtZM[iCut]);
352
353
354 fBGHandler[iCut] = new AliGammaConversionAODBGHandler( collisionSystem,centMin,centMax,
355 ((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->GetNumberOfBGEvents(),
356 ((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->UseTrackMultiplicity());
357 }
358}
359
360//______________________________________________________________________
361void AliAnalysisTaskNeutralMesonToPiPlPiMiPiZero::UserCreateOutputObjects()
362{
363 //
364 // Create ouput objects
365 //
366
367 // Create the output container
368 if(fOutputContainer != NULL){
369 delete fOutputContainer;
370 fOutputContainer = NULL;
371 }
372 if(fOutputContainer == NULL){
373 fOutputContainer = new TList();
374 fOutputContainer->SetOwner(kTRUE);
375 }
376
377 fGoodConvGammas = new TList();
378 fClusterCandidates = new TList();
379 fNeutralPionCandidates = new TList();
380 fGoodVirtualParticles = new TList();
381
382 fCutFolder = new TList*[fnCuts];
383 fESDList = new TList*[fnCuts];
384 fBackList = new TList*[fnCuts];
385 fMotherList = new TList*[fnCuts];
386 fHistoNEvents = new TH1I*[fnCuts];
387 fHistoNGoodESDTracks = new TH1I*[fnCuts];
388 fProfileEtaShift = new TProfile*[fnCuts];
389 if (fNeutralPionMode < 2){
390 fHistoConvGammaPt = new TH1F*[fnCuts];
391 fHistoConvGammaEta = new TH1F*[fnCuts];
392 }
393 if (fNeutralPionMode > 0){
394 fHistoClusterGammaPt = new TH1F*[fnCuts];
395 fHistoClusterGammaEta = new TH1F*[fnCuts];
396 }
397 fHistoNegPionPt = new TH1F*[fnCuts];
398 fHistoPosPionPt = new TH1F*[fnCuts];
399 fHistoNegPionPhi = new TH1F*[fnCuts];
400 fHistoPosPionPhi = new TH1F*[fnCuts];
401 fHistoPionPionInvMassPt = new TH2F*[fnCuts];
402
403 if( fDoMesonQA ) {
404 fHistoNegPionEta = new TH1F*[fnCuts];
405 fHistoPosPionEta = new TH1F*[fnCuts];
406 fHistoNegPionClsTPC = new TH2F*[fnCuts];
407 fHistoPosPionClsTPC = new TH2F*[fnCuts];
408 fHistoPionDCAxy = new TH2F*[fnCuts];
409 fHistoPionDCAz = new TH2F*[fnCuts];
410 fHistoPionTPCdEdxNSigma = new TH2F*[fnCuts];
411 fHistoPionTPCdEdx = new TH2F*[fnCuts];
412 }
413 fHistoGammaGammaInvMassPt = new TH2F*[fnCuts];
414 fHistoMotherInvMassPt = new TH2F*[fnCuts];
415 fHistoMotherBackInvMassPt = new TH2F*[fnCuts];
416
417 for(Int_t iCut = 0; iCut<fnCuts;iCut++){
418 TString cutstringEvent = ((AliConvEventCuts*)fEventCutArray->At(iCut))->GetCutNumber();
419 TString cutstringPion = ((AliPrimaryPionCuts*)fPionCutArray->At(iCut))->GetCutNumber();
420 TString cutstringConvGamma = "";
421 if (fNeutralPionMode < 2) cutstringConvGamma = ((AliConversionPhotonCuts*)fGammaCutArray->At(iCut))->GetCutNumber();
422 TString cutstringCaloGamma = "";
423 if (fNeutralPionMode > 0) cutstringCaloGamma = ((AliCaloPhotonCuts*)fClusterCutArray->At(iCut))->GetCutNumber();
424 TString cutstringNeutralPion= ((AliConversionMesonCuts*)fNeutralPionMesonCutArray->At(iCut))->GetCutNumber();
425 TString cutstringMeson = ((AliConversionMesonCuts*)fMesonCutArray->At(iCut))->GetCutNumber();
426
427 TString fullCutString = "";
428 if (fNeutralPionMode == 0) fullCutString = Form("%i_%s_%s_%s_%s_%s",fNeutralPionMode,cutstringEvent.Data(),cutstringConvGamma.Data(),cutstringNeutralPion.Data(), cutstringPion.Data(),cutstringMeson.Data());
429 else if (fNeutralPionMode == 1) fullCutString = Form("%i_%s_%s_%s_%s_%s_%s",fNeutralPionMode,cutstringEvent.Data(),cutstringConvGamma.Data(),cutstringCaloGamma.Data(),cutstringNeutralPion.Data(), cutstringPion.Data(),cutstringMeson.Data());
430 else if (fNeutralPionMode == 2) fullCutString = Form("%i_%s_%s_%s_%s_%s",fNeutralPionMode,cutstringEvent.Data(),cutstringCaloGamma.Data(),cutstringNeutralPion.Data(), cutstringPion.Data(),cutstringMeson.Data());
431 TString nameCutFolder = Form("Cut Number %s", fullCutString.Data());
432 TString nameESDList = Form("%s ESD histograms", fullCutString.Data());
433
434
435 fCutFolder[iCut] = new TList();
436 fCutFolder[iCut]->SetName(nameCutFolder.Data());
437 fCutFolder[iCut]->SetOwner(kTRUE);
438 fOutputContainer->Add(fCutFolder[iCut]);
439
440 fESDList[iCut] = new TList();
441 fESDList[iCut]->SetName(nameESDList.Data());
442 fESDList[iCut]->SetOwner(kTRUE);
443
444 fHistoNEvents[iCut] = new TH1I("NEvents","NEvents",9,-0.5,8.5);
445 fHistoNEvents[iCut]->GetXaxis()->SetBinLabel(1,"Accepted");
446 fHistoNEvents[iCut]->GetXaxis()->SetBinLabel(2,"Centrality");
447 fHistoNEvents[iCut]->GetXaxis()->SetBinLabel(3,"Missing MC");
448 fHistoNEvents[iCut]->GetXaxis()->SetBinLabel(4,"Trigger");
449 fHistoNEvents[iCut]->GetXaxis()->SetBinLabel(5,"Vertex Z");
450 fHistoNEvents[iCut]->GetXaxis()->SetBinLabel(6,"Cont. Vertex");
451 fHistoNEvents[iCut]->GetXaxis()->SetBinLabel(7,"Pile-Up");
452 fHistoNEvents[iCut]->GetXaxis()->SetBinLabel(8,"no SDD");
453 fHistoNEvents[iCut]->GetXaxis()->SetBinLabel(9,"no V0AND");
454 fESDList[iCut]->Add(fHistoNEvents[iCut]);
455
456 if(fIsHeavyIon>0) fHistoNGoodESDTracks[iCut] = new TH1I("GoodESDTracks","GoodESDTracks",3000,0,3000);
457 else fHistoNGoodESDTracks[iCut] = new TH1I("GoodESDTracks","GoodESDTracks",200,0,200);
458 fESDList[iCut]->Add(fHistoNGoodESDTracks[iCut]);
459
460 fProfileEtaShift[iCut] = new TProfile("Eta Shift","Eta Shift",1, -0.5,0.5);
461 fESDList[iCut]->Add(fProfileEtaShift[iCut]);
462 if (fNeutralPionMode < 2){
463 fHistoConvGammaPt[iCut] = new TH1F("ESD_ConvGamma_Pt","ESD_ConvGamma_Pt",250,0,25);
464 fESDList[iCut]->Add(fHistoConvGammaPt[iCut]);
465 fHistoConvGammaEta[iCut] = new TH1F("ESD_ConvGamma_Eta","ESD_ConvGamma_Eta",600,-1.5,1.5);
466 fESDList[iCut]->Add(fHistoConvGammaEta[iCut]);
467 }
468 if (fNeutralPionMode > 0){
469 fHistoClusterGammaPt[iCut] = new TH1F("ESD_ClusterGamma_Pt","ESD_ClusterGamma_Pt",250,0,25);
470 fESDList[iCut]->Add(fHistoClusterGammaPt[iCut]);
471 fHistoClusterGammaEta[iCut] = new TH1F("ESD_ClusterGamma_Eta","ESD_ClusterGamma_Eta",600,-1.5,1.5);
472 fESDList[iCut]->Add(fHistoClusterGammaEta[iCut]);
473 }
474 fHistoNegPionPt[iCut] = new TH1F("ESD_PrimaryNegPions_Pt","ESD_PrimaryNegPions_Pt",1000,0,25);
475 fESDList[iCut]->Add(fHistoNegPionPt[iCut]);
476 fHistoPosPionPt[iCut] = new TH1F("ESD_PrimaryPosPions_Pt","ESD_PrimaryPosPions_Pt",1000,0,25);
477 fESDList[iCut]->Add(fHistoPosPionPt[iCut]);
478 fHistoNegPionPhi[iCut] = new TH1F("ESD_PrimaryNegPions_Phi","ESD_PrimaryNegPions_Phi",360,0,2*TMath::Pi());
479 fESDList[iCut]->Add(fHistoNegPionPhi[iCut]);
480 fHistoPosPionPhi[iCut] = new TH1F("ESD_PrimaryPosPions_Phi","ESD_PrimaryPosPions_Phi",360,0,2*TMath::Pi());
481 fESDList[iCut]->Add(fHistoPosPionPhi[iCut]);
482 fHistoPionPionInvMassPt[iCut] = new TH2F("ESD_PiPlusPiNeg_InvMassPt","ESD_PiPlusPiNeg_InvMassPt",2000,0.,2.,200,0.,20.);
483 fESDList[iCut]->Add(fHistoPionPionInvMassPt[iCut]);
484
485 if ( fDoMesonQA ) {
486 fHistoNegPionEta[iCut] = new TH1F("ESD_PrimaryNegPions_Eta","ESD_PrimaryNegPions_Eta",600,-1.5,1.5);
487 fESDList[iCut]->Add(fHistoNegPionEta[iCut]);
488 fHistoPosPionEta[iCut] = new TH1F("ESD_PrimaryPosPions_Eta","ESD_PrimaryPosPions_Eta",600,-1.5,1.5);
489 fESDList[iCut]->Add(fHistoPosPionEta[iCut]);
490 fHistoNegPionClsTPC[iCut] = new TH2F("ESD_PrimaryNegPions_ClsTPC","ESD_PrimaryNegPions_ClsTPC",100,0,1,400,0.,10.);
491 fESDList[iCut]->Add(fHistoNegPionClsTPC[iCut]);
492 fHistoPosPionClsTPC[iCut] = new TH2F("ESD_PrimaryPosPions_ClsTPC","ESD_PrimaryPosPions_ClsTPC",100,0,1,400,0.,10.);
493 fESDList[iCut]->Add(fHistoPosPionClsTPC[iCut]);
494 fHistoPionDCAxy[iCut] = new TH2F("ESD_PrimaryPions_DCAxy","ESD_PrimaryPions_DCAxy",800,-4.0,4.0,400,0.,10.);
495 fESDList[iCut]->Add(fHistoPionDCAxy[iCut]);
496 fHistoPionDCAz[iCut] = new TH2F("ESD_PrimaryPions_DCAz","ESD_PrimaryPions_DCAz",800,-4.0,4.0,400,0.,10.);
497 fESDList[iCut]->Add(fHistoPionDCAz[iCut]);
498 fHistoPionTPCdEdxNSigma[iCut] = new TH2F("ESD_PrimaryPions_TPCdEdx","ESD_PrimaryPions_TPCdEdx",150,0.05,20,400,-10,10);
499 fESDList[iCut]->Add(fHistoPionTPCdEdxNSigma[iCut]);
500 fHistoPionTPCdEdx[iCut] =new TH2F("ESD_PrimaryPions_TPCdEdxSignal","ESD_PrimaryPions_TPCdEdxSignal" ,150,0.05,20.0,800,0.0,200);
501 fESDList[iCut]->Add(fHistoPionTPCdEdx[iCut]);
502 }
503 fHistoGammaGammaInvMassPt[iCut] = new TH2F("ESD_GammaGamma_InvMass_Pt","ESD_GammaGamma_InvMass_Pt",450,0.,0.45,250,0,25);
504 fESDList[iCut]->Add(fHistoGammaGammaInvMassPt[iCut]);
505 fHistoMotherInvMassPt[iCut] = new TH2F("ESD_Mother_InvMass_Pt","ESD_Mother_InvMass_Pt",500,0.4,0.9,250,0,25);
506 fESDList[iCut]->Add(fHistoMotherInvMassPt[iCut]);
507 fHistoMotherBackInvMassPt[iCut] = new TH2F("ESD_Background_InvMass_Pt","ESD_Background_InvMass_Pt",500,0.4,0.9,250,0,25);
508 fESDList[iCut]->Add(fHistoMotherBackInvMassPt[iCut]);
509
510 if ( fDoMesonQA ) {
511 TAxis *AxisAfter = fHistoPionTPCdEdxNSigma[iCut]->GetXaxis();
512 Int_t bins = AxisAfter->GetNbins();
513 Double_t from = AxisAfter->GetXmin();
514 Double_t to = AxisAfter->GetXmax();
515 Double_t *newBins = new Double_t[bins+1];
516 newBins[0] = from;
517 Double_t factor = TMath::Power(to/from, 1./bins);
518 for(Int_t i=1; i<=bins; ++i) newBins[i] = factor * newBins[i-1];
519
520 AxisAfter->Set(bins, newBins);
521 AxisAfter = fHistoPionTPCdEdx[iCut]->GetXaxis();
522 AxisAfter->Set(bins, newBins);
523
524 delete [] newBins;
525 }
526
527 fCutFolder[iCut]->Add(fESDList[iCut]);
528
529 }
530
531 if( fIsMC ){
532 // MC Histogramms
533 fMCList = new TList*[fnCuts];
534 // True Histogramms
535 fTrueList = new TList*[fnCuts];
536 if (fNeutralPionMode < 2){
537 fHistoTrueConvGammaPt = new TH1F*[fnCuts];
538 fHistoTrueConvGammaFromNeutralMesonPt = new TH1F*[fnCuts];
539 }
540 if (fNeutralPionMode > 0){
541 fHistoTrueClusterGammaPt = new TH1F*[fnCuts];
542 fHistoTrueClusterGammaFromNeutralMesonPt = new TH1F*[fnCuts];
543 }
544 fHistoTruePosPionPt = new TH1F*[fnCuts];
545 fHistoTrueNegPionPt = new TH1F*[fnCuts];
546 fHistoTruePosPionFromNeutralMesonPt = new TH1F*[fnCuts];
547 fHistoTrueNegPionFromNeutralMesonPt = new TH1F*[fnCuts];
548
549
550 fHistoMCAllGammaPt = new TH1F*[fnCuts];
551 if (fNeutralPionMode < 2){
552 fHistoMCConvGammaPt = new TH1F*[fnCuts];
553 }
554 fHistoMCAllPosPionsPt = new TH1F*[fnCuts];
555 fHistoMCAllNegPionsPt = new TH1F*[fnCuts];
556 fHistoMCGammaFromNeutralMesonPt = new TH1F*[fnCuts];
557 fHistoMCPosPionsFromNeutralMesonPt = new TH1F*[fnCuts];
558 fHistoMCNegPionsFromNeutralMesonPt = new TH1F*[fnCuts];
559
560// hMCPi0DalitzGammaPt = new TH1F*[fnCuts];
561// hMCPi0DalitzElectronPt = new TH1F*[fnCuts];
562// hMCPi0DalitzPositronPt = new TH1F*[fnCuts];
563
564 fHistoMCEtaPiPlPiMiPiZeroPt = new TH1F*[fnCuts];
565 fHistoMCEtaPiPlPiMiPiZeroInAccPt = new TH1F*[fnCuts];
566 fHistoMCOmegaPiPlPiMiPiZeroPt = new TH1F*[fnCuts];
567 fHistoMCOmegaPiPlPiMiPiZeroInAccPt = new TH1F*[fnCuts];
568
569 fHistoTrueMotherPiPlPiMiPiZeroInvMassPt = new TH2F*[fnCuts];
570 fHistoTrueMotherGammaGammaInvMassPt = new TH2F*[fnCuts];
571 fHistoTrueMotherGammaGammaFromEtaInvMassPt = new TH2F*[fnCuts];
572 fHistoTrueMotherGammaGammaFromOmegaInvMassPt = new TH2F*[fnCuts];
573// if (fDoMesonQA){
574 fHistoTruePionPionInvMassPt = new TH2F*[fnCuts];
575 fHistoTruePionPionFromSameMotherInvMassPt = new TH2F*[fnCuts];
576 fHistoTruePionPionFromEtaInvMassPt = new TH2F*[fnCuts];
577 fHistoTruePionPionFromOmegaInvMassPt = new TH2F*[fnCuts];
578// }
579
580 for(Int_t iCut = 0; iCut<fnCuts;iCut++){
581 TString cutstringEvent = ((AliConvEventCuts*)fEventCutArray->At(iCut))->GetCutNumber();
582 TString cutstringPion = ((AliPrimaryPionCuts*)fPionCutArray->At(iCut))->GetCutNumber();
583 TString cutstringConvGamma = "";
584 if (fNeutralPionMode < 2) cutstringConvGamma = ((AliConversionPhotonCuts*)fGammaCutArray->At(iCut))->GetCutNumber();
585 TString cutstringCaloGamma = "";
586 if (fNeutralPionMode > 0) cutstringCaloGamma = ((AliCaloPhotonCuts*)fClusterCutArray->At(iCut))->GetCutNumber();
587 TString cutstringNeutralPion= ((AliConversionMesonCuts*)fNeutralPionMesonCutArray->At(iCut))->GetCutNumber();
588 TString cutstringMeson = ((AliConversionMesonCuts*)fMesonCutArray->At(iCut))->GetCutNumber();
589
590 TString fullCutString = "";
591 if (fNeutralPionMode == 0) fullCutString = Form("%i_%s_%s_%s_%s_%s",fNeutralPionMode,cutstringEvent.Data(),cutstringConvGamma.Data(),cutstringNeutralPion.Data(), cutstringPion.Data(),cutstringMeson.Data());
592 else if (fNeutralPionMode == 1) fullCutString = Form("%i_%s_%s_%s_%s_%s_%s",fNeutralPionMode,cutstringEvent.Data(),cutstringConvGamma.Data(),cutstringCaloGamma.Data(),cutstringNeutralPion.Data(), cutstringPion.Data(),cutstringMeson.Data());
593 else if (fNeutralPionMode == 2) fullCutString = Form("%i_%s_%s_%s_%s_%s",fNeutralPionMode,cutstringEvent.Data(),cutstringCaloGamma.Data(),cutstringNeutralPion.Data(), cutstringPion.Data(),cutstringMeson.Data());
594 TString nameMCList = Form("%s MC histograms", fullCutString.Data());
595 TString nameTrueRecList = Form("%s True histograms", fullCutString.Data());
596
597 fMCList[iCut] = new TList();
598 fMCList[iCut]->SetName(nameMCList.Data());
599 fMCList[iCut]->SetOwner(kTRUE);
600 fCutFolder[iCut]->Add(fMCList[iCut]);
601
602 fHistoMCAllGammaPt[iCut] = new TH1F("MC_AllGamma_Pt","MC_AllGamma_Pt",250,0,25);
603 fMCList[iCut]->Add(fHistoMCAllGammaPt[iCut]);
604 if (fNeutralPionMode < 2){
605 fHistoMCConvGammaPt[iCut] = new TH1F("MC_ConvGamma_Pt","MC_ConvGamma_Pt",250,0,25);
606 fMCList[iCut]->Add(fHistoMCConvGammaPt[iCut]);
607 }
608
609 fHistoMCAllPosPionsPt[iCut] = new TH1F("MC_AllPosPions_Pt","MC_AllPosPions_Pt",1000,0,25);
610 fMCList[iCut]->Add(fHistoMCAllPosPionsPt[iCut]);
611 fHistoMCAllNegPionsPt[iCut] = new TH1F("MC_AllNegPions_Pt","MC_AllNegPions_Pt",1000,0,25);
612 fMCList[iCut]->Add(fHistoMCAllNegPionsPt[iCut]);
613 fHistoMCGammaFromNeutralMesonPt[iCut] = new TH1F("MC_GammaFromNeutralMeson_Pt","MC_GammaFromNeutralMeson_Pt",250,0,25);
614 fMCList[iCut]->Add(fHistoMCGammaFromNeutralMesonPt[iCut]);
615 fHistoMCPosPionsFromNeutralMesonPt[iCut] = new TH1F("MC_PosPionsFromNeutralMeson_Pt","MC_PosPionsFromNeutralMeson_Pt",1000,0,25);
616 fMCList[iCut]->Add(fHistoMCPosPionsFromNeutralMesonPt[iCut]);
617 fHistoMCNegPionsFromNeutralMesonPt[iCut] = new TH1F("MC_NegPionsFromNeutralMeson_Pt","MC_NegPionsFromNeutralMeson_Pt",1000,0,25);
618 fMCList[iCut]->Add(fHistoMCNegPionsFromNeutralMesonPt[iCut]);
619
620 fHistoMCEtaPiPlPiMiPiZeroPt[iCut] = new TH1F("MC_Eta_Pt","MC_Eta_Pt",250,0,25);
621 fHistoMCEtaPiPlPiMiPiZeroPt[iCut]->Sumw2();
622 fMCList[iCut]->Add(fHistoMCEtaPiPlPiMiPiZeroPt[iCut]);
623
624 fHistoMCEtaPiPlPiMiPiZeroInAccPt[iCut] = new TH1F("MC_EtaInAcc_Pt","MC_EtaInAcc_Pt",250,0,25);
625 fHistoMCEtaPiPlPiMiPiZeroInAccPt[iCut]->Sumw2();
626 fMCList[iCut]->Add(fHistoMCEtaPiPlPiMiPiZeroInAccPt[iCut]);
627
628 fHistoMCOmegaPiPlPiMiPiZeroPt[iCut] = new TH1F("MC_Omega_Pt","MC_Omega_Pt",250,0,25);
629 fHistoMCOmegaPiPlPiMiPiZeroPt[iCut]->Sumw2();
630 fMCList[iCut]->Add(fHistoMCOmegaPiPlPiMiPiZeroPt[iCut]);
631
632 fHistoMCOmegaPiPlPiMiPiZeroInAccPt[iCut] = new TH1F("MC_OmegaInAcc_Pt","MC_OmegaInAcc_Pt",250,0,25);
633 fHistoMCOmegaPiPlPiMiPiZeroInAccPt[iCut]->Sumw2();
634 fMCList[iCut]->Add(fHistoMCOmegaPiPlPiMiPiZeroInAccPt[iCut]);
635
636 fTrueList[iCut] = new TList();
637 fTrueList[iCut]->SetName(nameTrueRecList.Data());
638 fTrueList[iCut]->SetOwner(kTRUE);
639 fCutFolder[iCut]->Add(fTrueList[iCut]);
640
641 if (fNeutralPionMode < 2){
642 fHistoTrueConvGammaPt[iCut] = new TH1F("ESD_TrueConvGamma_Pt","ESD_TrueConvGamma_Pt",250,0,25);
643 fTrueList[iCut]->Add(fHistoTrueConvGammaPt[iCut]);
644 fHistoTrueConvGammaFromNeutralMesonPt[iCut] = new TH1F("ESD_TrueConvGammaFromNeutralMeson_Pt","ESD_TrueConvGammaFromNeutralMeson_Pt",250,0,25);
645 fTrueList[iCut]->Add(fHistoTrueConvGammaFromNeutralMesonPt[iCut]);
646 }
647 if (fNeutralPionMode > 0){
648 fHistoTrueClusterGammaPt[iCut] = new TH1F("ESD_TrueClusterGamma_Pt","ESD_TrueClusterGamma_Pt",250,0,25);
649 fTrueList[iCut]->Add(fHistoTrueClusterGammaPt[iCut]);
650 fHistoTrueClusterGammaFromNeutralMesonPt[iCut] = new TH1F("ESD_TrueClusterGammaFromNeutralMeson_Pt","ESD_TrueClusterGammaFromNeutralMeson_Pt",250,0,25);
651 fTrueList[iCut]->Add(fHistoTrueClusterGammaFromNeutralMesonPt[iCut]);
652 }
653 fHistoTruePosPionPt[iCut] = new TH1F("ESD_TruePosPion_Pt","ESD_TruePosPion_Pt",1000,0,25);
654 fTrueList[iCut]->Add(fHistoTruePosPionPt[iCut]);
655 fHistoTrueNegPionPt[iCut] = new TH1F("ESD_TrueNegPion_Pt","ESD_TrueNegPion_Pt",1000,0,25);
656 fTrueList[iCut]->Add(fHistoTrueNegPionPt[iCut]);
657
658 fHistoTrueNegPionFromNeutralMesonPt[iCut] = new TH1F("ESD_TrueNegPionFromNeutralMeson_Pt","ESD_TrueNegPionFromNeutralMeson_Pt",1000,0,25);
659 fTrueList[iCut]->Add(fHistoTrueNegPionFromNeutralMesonPt[iCut]);
660 fHistoTruePosPionFromNeutralMesonPt[iCut] = new TH1F("ESD_TruePosPionFromNeutralMeson_Pt","ESD_TruePosPionFromNeutralMeson_Pt",1000,0,25);
661 fTrueList[iCut]->Add(fHistoTruePosPionFromNeutralMesonPt[iCut]);
662
663 fHistoTrueMotherPiPlPiMiPiZeroInvMassPt[iCut] = new TH2F("ESD_TrueMotherPiPlPiMiPiZero_InvMass_Pt","ESD_TrueMotherPiPlPiMiPiZero_InvMass_Pt",500,0.4,0.9,250,0,25);
664 fHistoTrueMotherPiPlPiMiPiZeroInvMassPt[iCut]->Sumw2();
665 fTrueList[iCut]->Add(fHistoTrueMotherPiPlPiMiPiZeroInvMassPt[iCut]);
666
667 fHistoTrueMotherGammaGammaInvMassPt[iCut] = new TH2F("ESD_TrueMotherGG_InvMass_Pt","ESD_TrueMotherGG_InvMass_Pt",450,0.,0.45,250,0,25);
668 fTrueList[iCut]->Add(fHistoTrueMotherGammaGammaInvMassPt[iCut]);
669 fHistoTrueMotherGammaGammaFromEtaInvMassPt[iCut] = new TH2F("ESD_TrueMotherGGFromEta_InvMass_Pt","ESD_TrueMotherGGFromEta_InvMass_Pt",450,0.,0.45,250,0,25);
670 fTrueList[iCut]->Add(fHistoTrueMotherGammaGammaFromEtaInvMassPt[iCut]);
671 fHistoTrueMotherGammaGammaFromOmegaInvMassPt[iCut] = new TH2F("ESD_TrueMotherGGFromOmega_InvMass_Pt","ESD_TrueMotherGGFromOmega_InvMass_Pt",450,0.,0.45,250,0,25);
672 fTrueList[iCut]->Add(fHistoTrueMotherGammaGammaFromOmegaInvMassPt[iCut]);
673
674
675// if (fDoMesonQA){
676 fHistoTruePionPionInvMassPt[iCut] = new TH2F("ESD_TruePiPlusPiNeg_InvMassPt","ESD_TruePiPlusPiNeg_InvMassPt",2000,0.,2.,200,0.,20.);
677 fTrueList[iCut]->Add(fHistoTruePionPionInvMassPt[iCut]);
678 fHistoTruePionPionFromSameMotherInvMassPt[iCut] = new TH2F("ESD_TruePiPlusPiNegFromSameMother_InvMassPt","ESD_TruePiPlusPiNegFromSameMother_InvMassPt",2000,0.,2.,200,0.,20.);
679 fTrueList[iCut]->Add(fHistoTruePionPionFromSameMotherInvMassPt[iCut]);
680 fHistoTruePionPionFromEtaInvMassPt[iCut] = new TH2F("ESD_TruePiPlusPiNegFromEta_InvMassPt","ESD_TruePiPlusPiNegFromEta_InvMassPt",2000,0.,2.,200,0.,20.);
681 fTrueList[iCut]->Add(fHistoTruePionPionFromEtaInvMassPt[iCut]);
682 fHistoTruePionPionFromOmegaInvMassPt[iCut] = new TH2F("ESD_TruePiPlusPiNegFromOmega_InvMassPt","ESD_TruePiPlusPiNegFromOmega_InvMassPt",2000,0.,2.,200,0.,20.);
683 fTrueList[iCut]->Add(fHistoTruePionPionFromOmegaInvMassPt[iCut]);
684
685 // }
686 }
687 }
688
689
690
691 InitBack(); // Init Background Handler
692
693 fV0Reader=(AliV0ReaderV1*)AliAnalysisManager::GetAnalysisManager()->GetTask("V0ReaderV1");
694 if(!fV0Reader){printf("Error: No V0 Reader");return;} // GetV0Reader
695
696 if(fV0Reader)
697 if((AliConversionPhotonCuts*)fV0Reader->GetConversionCuts())
698 if(((AliConversionPhotonCuts*)fV0Reader->GetConversionCuts())->GetCutHistograms())
699 fOutputContainer->Add(((AliConversionPhotonCuts*)fV0Reader->GetConversionCuts())->GetCutHistograms());
700
701
702
703 fPionSelector=(AliPrimaryPionSelector*)AliAnalysisManager::GetAnalysisManager()->GetTask("PionSelector");
704 if(!fPionSelector){printf("Error: No PionSelector");return;} // GetV0Reader
705
706 if( fPionSelector ){
707 if ( ((AliPrimaryPionCuts*)fPionSelector->GetPrimaryPionCuts())->GetCutHistograms() ){
708 fOutputContainer->Add( ((AliPrimaryPionCuts*)fPionSelector->GetPrimaryPionCuts())->GetCutHistograms() );
709 }
710 }
711
712 for(Int_t iCut = 0; iCut<fnCuts;iCut++){
713 if( fPionCutArray ){
714 if( ((AliPrimaryPionCuts*)fPionCutArray->At(iCut))->GetCutHistograms() ) {
715 fCutFolder[iCut]->Add( ((AliPrimaryPionCuts*)fPionCutArray->At(iCut))->GetCutHistograms() );
716 }
717 }
718 if (fNeutralPionMode < 2){
719 if( fGammaCutArray ) {
720 if( ((AliConversionPhotonCuts*)fGammaCutArray->At(iCut))->GetCutHistograms() ) {
721 fCutFolder[iCut]->Add( ((AliConversionPhotonCuts*)fGammaCutArray->At(iCut))->GetCutHistograms() );
722 }
723 }
724 }
725 if (fNeutralPionMode > 0){
726 if( fClusterCutArray ) {
727 if( ((AliCaloPhotonCuts*)fClusterCutArray->At(iCut))->GetCutHistograms() ) {
728 fCutFolder[iCut]->Add( ((AliCaloPhotonCuts*)fClusterCutArray->At(iCut))->GetCutHistograms() );
729 }
730 }
731 }
732 if( fNeutralPionMesonCutArray ) {
733 if( ((AliConversionMesonCuts*)fNeutralPionMesonCutArray->At(iCut))->GetCutHistograms() ) {
734 fCutFolder[iCut]->Add( ((AliConversionMesonCuts*)fNeutralPionMesonCutArray->At(iCut))->GetCutHistograms());
735 }
736 }
737 if( fMesonCutArray ) {
738 if( ((AliConversionMesonCuts*)fMesonCutArray->At(iCut))->GetCutHistograms() ) {
739 fCutFolder[iCut]->Add( ((AliConversionMesonCuts*)fMesonCutArray->At(iCut))->GetCutHistograms());
740 }
741 }
742 }
743
744 PostData(1, fOutputContainer);
745
746}
747
748//______________________________________________________________________
749void AliAnalysisTaskNeutralMesonToPiPlPiMiPiZero::UserExec(Option_t *){
750
751 //
752 // Execute analysis for current event
753 //
754
755 fV0Reader=(AliV0ReaderV1*)AliAnalysisManager::GetAnalysisManager()->GetTask("V0ReaderV1");
756 if(!fV0Reader){printf("Error: No V0 Reader");return;} // GetV0Reader
757
758 Int_t eventQuality = ((AliConvEventCuts*)fV0Reader->GetEventCuts())->GetEventQuality();
759
760 if(eventQuality == 2 || eventQuality == 3){// Event Not Accepted due to MC event missing or wrong trigger for V0ReaderV1
761 for(Int_t iCut = 0; iCut<fnCuts; iCut++){
762 fHistoNEvents[iCut]->Fill(eventQuality);
763 }
764 return;
765 }
766
767 fPionSelector=(AliPrimaryPionSelector*)AliAnalysisManager::GetAnalysisManager()->GetTask("PionSelector");
768 if(!fPionSelector){printf("Error: No PionSelector");return;} // GetV0Reader
769
770 if(fIsMC) fMCEvent = MCEvent();
771 fESDEvent = (AliESDEvent*)InputEvent();
772 fReaderGammas = fV0Reader->GetReconstructedGammas(); // Gammas from default Cut
773 fSelectorNegPionIndex = fPionSelector->GetReconstructedNegPionIndex(); // Electrons from default Cut
774 fSelectorPosPionIndex = fPionSelector->GetReconstructedPosPionIndex(); // Positrons from default Cut
775
776 fNumberOfESDTracks = fV0Reader->GetNumberOfPrimaryTracks();
777 //AddTaskContainers(); //Add conatiner
778
779 for(Int_t iCut = 0; iCut<fnCuts; iCut++){
780 fiCut = iCut;
781 Int_t eventNotAccepted = ((AliConvEventCuts*)fEventCutArray->At(iCut))->IsEventAcceptedByCut(fV0Reader->GetEventCuts(),fInputEvent,fMCEvent,fIsHeavyIon);
782
783 if(eventNotAccepted){
784 // cout << "event rejected due to wrong trigger: " <<eventNotAccepted << endl;
785 fHistoNEvents[iCut]->Fill(eventNotAccepted); // Check Centrality, PileUp, SDD and V0AND --> Not Accepted => eventQuality = 1
786 continue;
787 }
788
789 if(eventQuality != 0){// Event Not Accepted
790 // cout << "event rejected due to: " <<eventQuality << endl;
791 fHistoNEvents[iCut]->Fill(eventQuality);
792 continue;
793 }
794
795 fHistoNEvents[iCut]->Fill(eventQuality);
796 fHistoNGoodESDTracks[iCut]->Fill(fNumberOfESDTracks);
797
798 if(fMCEvent){ // Process MC Particle
799 fMCStack = fMCEvent->Stack();
800 if(((AliConvEventCuts*)fEventCutArray->At(iCut))->GetSignalRejection() != 0){
801 ((AliConvEventCuts*)fEventCutArray->At(iCut))->GetNotRejectedParticles(((AliConvEventCuts*)fEventCutArray->At(iCut))->GetSignalRejection(),
802 ((AliConvEventCuts*)fEventCutArray->At(iCut))->GetAcceptedHeader(),
803 fMCEvent);
804 }
805 ProcessMCParticles();
806 }
807
808 if (fNeutralPionMode < 2){
809 ProcessConversionPhotonCandidates(); // Process this cuts conversion gammas
810 }
811 if (fNeutralPionMode > 0){
812 ProcessCaloPhotonCandidates(); // Process this cuts calo gammas
813 }
814
815 if (fNeutralPionMode == 0 ){
816 ProcessNeutralPionCandidatesPureConversions(); // Process neutral pion candidates purely from conversions
817 }
818 if (fNeutralPionMode == 1){
819 ProcessNeutralPionCandidatesMixedConvCalo(); // Process neutral pion candidates mixed conv and calo
820 }
821 if (fNeutralPionMode == 2){
822 ProcessNeutralPionCandidatesPureCalo(); // Process neutral pion candidates purely from calo
823 }
824
825 ProcessPionCandidates(); // Process this cuts gammas
826
827 CalculateMesonCandidates();
828 CalculateBackground();
829 UpdateEventByEventData();
830
831
832 if (fGoodConvGammas->GetEntries()>0) fGoodConvGammas->Clear();
833 if (fClusterCandidates->GetEntries()>0) fClusterCandidates->Clear();
834 if (fNeutralPionCandidates->GetEntries()>0) fNeutralPionCandidates->Clear();
835 fGoodVirtualParticles->Clear(); // delete this cuts good gammas
836 }
837
838 fSelectorNegPionIndex.clear();
839 fSelectorPosPionIndex.clear();
840
841 PostData( 1, fOutputContainer );
842}
843
844Bool_t AliAnalysisTaskNeutralMesonToPiPlPiMiPiZero::Notify(){
845 for(Int_t iCut = 0; iCut<fnCuts;iCut++){
846 if( !((AliConvEventCuts*)fEventCutArray->At(iCut))->GetDoEtaShift() ){
847 fProfileEtaShift[iCut]->Fill(0.,0.);
848 continue; // No Eta Shift requested, continue
849 }
850 if( ((AliConvEventCuts*)fEventCutArray->At(iCut))->GetEtaShift() == 0.0){ // Eta Shift requested but not set, get shift automatically
851 ((AliConvEventCuts*)fEventCutArray->At(iCut))->GetCorrectEtaShiftFromPeriod(fV0Reader->GetPeriodName());
852 ((AliConvEventCuts*)fEventCutArray->At(iCut))->DoEtaShift(kFALSE); // Eta Shift Set, make sure that it is called only once
853 ((AliPrimaryPionCuts*)fPionCutArray->At(iCut))->SetEtaShift( ((AliConvEventCuts*)fEventCutArray->At(iCut))->GetEtaShift() );
854 fProfileEtaShift[iCut]->Fill(0.,(((AliConvEventCuts*)fEventCutArray->At(iCut))->GetEtaShift()));
855 continue;
856 } else {
857 printf(" Eta t PiPlusPiMinus Gamma Task %s :: Eta Shift Manually Set to %f \n\n",
858 (((AliConvEventCuts*)fEventCutArray->At(iCut))->GetCutNumber()).Data(),((AliConvEventCuts*)fEventCutArray->At(iCut))->GetEtaShift());
859 ((AliConvEventCuts*)fEventCutArray->At(iCut))->DoEtaShift(kFALSE); // Eta Shift Set, make sure that it is called only once
860 ((AliPrimaryPionCuts*)fPionCutArray->At(iCut))->SetEtaShift( ((AliConvEventCuts*)fEventCutArray->At(iCut))->GetEtaShift() );
861 fProfileEtaShift[iCut]->Fill(0.,(((AliConvEventCuts*)fEventCutArray->At(iCut))->GetEtaShift()));
862 }
863 }
864 return kTRUE;
865}
866
867
868void AliAnalysisTaskNeutralMesonToPiPlPiMiPiZero::Terminate(const Option_t *){
869///Grid
870}
871
872
873//________________________________________________________________________
874void AliAnalysisTaskNeutralMesonToPiPlPiMiPiZero::ProcessCaloPhotonCandidates()
875{
876
877 Int_t nclus = 0;
878 nclus = fInputEvent->GetNumberOfCaloClusters();
879
880// cout << nclus << endl;
881
882 if(nclus == 0) return;
883
884 // vertex
885 Double_t vertex[3] = {0};
886 InputEvent()->GetPrimaryVertex()->GetXYZ(vertex);
887
888 // Loop over EMCal clusters
889 for(Long_t i = 0; i < nclus; i++){
890
891 AliVCluster* clus = NULL;
892 clus = fInputEvent->GetCaloCluster(i);
893 if (!clus) continue;
894 if(!((AliCaloPhotonCuts*)fClusterCutArray->At(fiCut))->ClusterIsSelected(clus,fInputEvent,fIsMC)) continue;
895 // TLorentzvector with cluster
896 TLorentzVector clusterVector;
897 clus->GetMomentum(clusterVector,vertex);
898
899 TLorentzVector* tmpvec = new TLorentzVector();
900 tmpvec->SetPxPyPzE(clusterVector.Px(),clusterVector.Py(),clusterVector.Pz(),clusterVector.E());
901
902 // convert to AODConversionPhoton
903 AliAODConversionPhoton *PhotonCandidate=new AliAODConversionPhoton(tmpvec);
904 if(!PhotonCandidate) continue;
905
906 // Flag Photon as CaloPhoton
907 PhotonCandidate->SetIsCaloPhoton();
908 PhotonCandidate->SetCaloClusterRef(i);
909 // get MC label
910 if(fIsMC){
911 Int_t* mclabelsCluster = clus->GetLabels();
912 PhotonCandidate->SetNCaloPhotonMCLabels(clus->GetNLabels());
913// cout << clus->GetNLabels() << endl;
914 if (clus->GetNLabels()>0){
915 for (Int_t k =0; k< (Int_t)clus->GetNLabels(); k++){
916 if (k< 20)PhotonCandidate->SetCaloPhotonMCLabel(k,mclabelsCluster[k]);
917// Int_t pdgCode = fMCStack->Particle(mclabelsCluster[k])->GetPdgCode();
918// cout << "label " << k << "\t" << mclabelsCluster[k] << " pdg code: " << pdgCode << endl;
919 }
920 }
921 }
922
923 fIsFromMBHeader = kTRUE;
924 // test whether largest contribution to cluster orginates in added signals
925 if (fIsMC && ((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsParticleFromBGEvent(PhotonCandidate->GetCaloPhotonMCLabel(0), fMCStack, fInputEvent) == 0) fIsFromMBHeader = kFALSE;
926
927 if (fIsFromMBHeader){
928 fHistoClusterGammaPt[fiCut]->Fill(PhotonCandidate->Pt());
929 fHistoClusterGammaEta[fiCut]->Fill(PhotonCandidate->Eta());
930 }
931 fClusterCandidates->Add(PhotonCandidate); // if no second loop is required add to events good gammas
932
933 if(fIsMC){
934// if(fInputEvent->IsA()==AliESDEvent::Class()){
935 ProcessTrueCaloPhotonCandidates(PhotonCandidate);
936// } else {
937// ProcessTrueClusterCandidatesAOD(PhotonCandidate);
938// }
939 }
940
941 delete tmpvec;
942 }
943
944}
945
946//________________________________________________________________________
947void AliAnalysisTaskNeutralMesonToPiPlPiMiPiZero::ProcessTrueCaloPhotonCandidates(AliAODConversionPhoton *TruePhotonCandidate)
948{
949 TParticle *Photon = NULL;
950 if (!TruePhotonCandidate->GetIsCaloPhoton()) AliFatal("CaloPhotonFlag has not been set task will abort");
951// fHistoTrueNLabelsInClus[fiCut]->Fill(TruePhotonCandidate->GetNCaloPhotonMCLabels());
952
953 if (TruePhotonCandidate->GetNCaloPhotonMCLabels()>0)Photon = fMCStack->Particle(TruePhotonCandidate->GetCaloPhotonMCLabel(0));
954 else return;
955
956 if(Photon == NULL){
957 // cout << "no photon" << endl;
958 return;
959 }
960
961// Int_t pdgCodeParticle = Photon->GetPdgCode();
962 TruePhotonCandidate->SetCaloPhotonMCFlags(fMCStack);
963
964 // True Photon
965 if(fIsFromMBHeader){
966 if(TruePhotonCandidate->GetCaloPhotonMCLabel(0) <= fMCStack->GetNprimary()){
967 if (TruePhotonCandidate->IsLargestComponentPhoton()){
968 fHistoTrueClusterGammaPt[fiCut]->Fill(TruePhotonCandidate->Pt());
969 if (GammaIsNeutralMesonPiPlPiMiPiZeroDaughter(TruePhotonCandidate->GetCaloPhotonMCLabel(0))){
970 fHistoTrueClusterGammaFromNeutralMesonPt[fiCut]->Fill(TruePhotonCandidate->Pt());
971 }
972 }
973 if (TruePhotonCandidate->IsLargestComponentElectron() && TruePhotonCandidate->IsConversion()){
974 fHistoTrueClusterGammaPt[fiCut]->Fill(TruePhotonCandidate->Pt());
975 if (GammaIsNeutralMesonPiPlPiMiPiZeroDaughter(TruePhotonCandidate->GetCaloPhotonMCLabel(0))){
976 fHistoTrueClusterGammaFromNeutralMesonPt[fiCut]->Fill(TruePhotonCandidate->Pt());
977 }
978 }
979 }
980 }
981 return;
982}
983
984
985
986//________________________________________________________________________
987void AliAnalysisTaskNeutralMesonToPiPlPiMiPiZero::ProcessConversionPhotonCandidates(){
988 Int_t nV0 = 0;
989 TList *GoodGammasStepOne = new TList();
990 TList *GoodGammasStepTwo = new TList();
991 // Loop over Photon Candidates allocated by ReaderV1
992
993 for(Int_t i = 0; i < fReaderGammas->GetEntriesFast(); i++){
994 AliAODConversionPhoton* PhotonCandidate = (AliAODConversionPhoton*) fReaderGammas->At(i);
995 if(!PhotonCandidate) continue;
996
997 fIsFromMBHeader = kTRUE;
998
999 if( fMCEvent && ((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetSignalRejection() != 0 ){
1000 Int_t isPosFromMBHeader
1001 = ((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsParticleFromBGEvent(PhotonCandidate->GetMCLabelPositive(), fMCStack, fInputEvent);
1002 if(isPosFromMBHeader == 0 && ((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetSignalRejection() != 3) continue;
1003 Int_t isNegFromMBHeader
1004 = ((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsParticleFromBGEvent(PhotonCandidate->GetMCLabelNegative(), fMCStack,fInputEvent);
1005 if(isNegFromMBHeader == 0 && ((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetSignalRejection() != 3) continue;
1006 if( (isNegFromMBHeader+isPosFromMBHeader) != 4) fIsFromMBHeader = kFALSE;
1007 }
1008
1009 if(!((AliConversionPhotonCuts*)fGammaCutArray->At(fiCut))->PhotonIsSelected(PhotonCandidate,fESDEvent)) continue;
1010
1011 if(!((AliConversionPhotonCuts*)fGammaCutArray->At(fiCut))->UseElecSharingCut() &&
1012 !((AliConversionPhotonCuts*)fGammaCutArray->At(fiCut))->UseToCloseV0sCut()){ // if no post reader loop is required add to events good gammas
1013
1014 fGoodConvGammas->Add(PhotonCandidate);
1015
1016 if(fIsFromMBHeader){
1017 fHistoConvGammaPt[fiCut]->Fill(PhotonCandidate->Pt());
1018 fHistoConvGammaEta[fiCut]->Fill(PhotonCandidate->Eta());
1019 }
1020
1021 if(fMCEvent){
1022 ProcessTrueConversionPhotonCandidates(PhotonCandidate);
1023 }
1024 } else if(((AliConversionPhotonCuts*)fGammaCutArray->At(fiCut))->UseElecSharingCut()){ // if Shared Electron cut is enabled, Fill array, add to step one
1025 ((AliConversionPhotonCuts*)fGammaCutArray->At(fiCut))->FillElectonLabelArray(PhotonCandidate,nV0);
1026 nV0++;
1027 GoodGammasStepOne->Add(PhotonCandidate);
1028 } else if(!((AliConversionPhotonCuts*)fGammaCutArray->At(fiCut))->UseElecSharingCut() &&
1029 ((AliConversionPhotonCuts*)fGammaCutArray->At(fiCut))->UseToCloseV0sCut()){ // shared electron is disabled, step one not needed -> step two
1030 GoodGammasStepTwo->Add(PhotonCandidate);
1031 }
1032 }
1033
1034
1035 if(((AliConversionPhotonCuts*)fGammaCutArray->At(fiCut))->UseElecSharingCut()){
1036 for(Int_t i = 0;i<GoodGammasStepOne->GetEntries();i++){
1037 AliAODConversionPhoton *PhotonCandidate= (AliAODConversionPhoton*) GoodGammasStepOne->At(i);
1038 if(!PhotonCandidate) continue;
1039 fIsFromMBHeader = kTRUE;
1040 if(fMCEvent && ((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetSignalRejection() != 0){
1041 Int_t isPosFromMBHeader
1042 = ((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsParticleFromBGEvent(PhotonCandidate->GetMCLabelPositive(), fMCStack,fInputEvent);
1043 Int_t isNegFromMBHeader
1044 = ((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsParticleFromBGEvent(PhotonCandidate->GetMCLabelNegative(), fMCStack,fInputEvent);
1045 if( (isNegFromMBHeader+isPosFromMBHeader) != 4) fIsFromMBHeader = kFALSE;
1046 }
1047 if(!((AliConversionPhotonCuts*)fGammaCutArray->At(fiCut))->RejectSharedElectronV0s(PhotonCandidate,i,GoodGammasStepOne->GetEntries())) continue;
1048 if(!((AliConversionPhotonCuts*)fGammaCutArray->At(fiCut))->UseToCloseV0sCut()){ // To Colse v0s cut diabled, step two not needed
1049 fGoodConvGammas->Add(PhotonCandidate);
1050 if(fIsFromMBHeader){
1051 fHistoConvGammaPt[fiCut]->Fill(PhotonCandidate->Pt());
1052 fHistoConvGammaEta[fiCut]->Fill(PhotonCandidate->Eta());
1053 }
1054 if(fMCEvent){
1055 ProcessTrueConversionPhotonCandidates(PhotonCandidate);
1056 }
1057 }
1058 else GoodGammasStepTwo->Add(PhotonCandidate); // Close v0s cut enabled -> add to list two
1059 }
1060 }
1061 if(((AliConversionPhotonCuts*)fGammaCutArray->At(fiCut))->UseToCloseV0sCut()){
1062 for(Int_t i = 0;i<GoodGammasStepTwo->GetEntries();i++){
1063 AliAODConversionPhoton* PhotonCandidate = (AliAODConversionPhoton*) GoodGammasStepTwo->At(i);
1064 if(!PhotonCandidate) continue;
1065
1066 if(fMCEvent && ((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetSignalRejection() != 0){
1067 Int_t isPosFromMBHeader
1068 = ((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsParticleFromBGEvent(PhotonCandidate->GetMCLabelPositive(), fMCStack,fInputEvent);
1069 Int_t isNegFromMBHeader
1070 = ((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsParticleFromBGEvent(PhotonCandidate->GetMCLabelNegative(), fMCStack,fInputEvent);
1071 if( (isNegFromMBHeader+isPosFromMBHeader) != 4) fIsFromMBHeader = kFALSE;
1072 }
1073
1074 if(!((AliConversionPhotonCuts*)fGammaCutArray->At(fiCut))->RejectToCloseV0s(PhotonCandidate,GoodGammasStepTwo,i)) continue;
1075 fGoodConvGammas->Add(PhotonCandidate); // Add gamma to current cut TList
1076
1077 if(fIsFromMBHeader){
1078 fHistoConvGammaPt[fiCut]->Fill(PhotonCandidate->Pt()); // Differences to old V0Reader in p_t due to conversion KF->TLorentzVector
1079 fHistoConvGammaEta[fiCut]->Fill(PhotonCandidate->Eta());
1080 }
1081
1082 if(fMCEvent){
1083 ProcessTrueConversionPhotonCandidates(PhotonCandidate);
1084 }
1085 }
1086 }
1087
1088 delete GoodGammasStepOne;
1089 GoodGammasStepOne = 0x0;
1090 delete GoodGammasStepTwo;
1091 GoodGammasStepTwo = 0x0;
1092}
1093
1094//________________________________________________________________________
1095void AliAnalysisTaskNeutralMesonToPiPlPiMiPiZero::ProcessTrueConversionPhotonCandidates(AliAODConversionPhoton *TruePhotonCandidate)
1096{
1097 // Process True Photons
1098 AliStack *MCStack = fMCEvent->Stack();
1099 TParticle *posDaughter = TruePhotonCandidate->GetPositiveMCDaughter(MCStack);
1100 TParticle *negDaughter = TruePhotonCandidate->GetNegativeMCDaughter(MCStack);
1101
1102 if(posDaughter == NULL || negDaughter == NULL) return; // One particle does not exist
1103 if(posDaughter->GetMother(0) != negDaughter->GetMother(0)){ // Not Same Mother == Combinatorial Bck
1104 return;
1105 }
1106
1107 else if (posDaughter->GetMother(0) == -1){
1108 return;
1109 }
1110
1111 if(TMath::Abs(posDaughter->GetPdgCode())!=11 || TMath::Abs(negDaughter->GetPdgCode())!=11) return; //One Particle is not electron
1112 if(posDaughter->GetPdgCode()==negDaughter->GetPdgCode()) return; // Same Charge
1113 if(posDaughter->GetUniqueID() != 5 || negDaughter->GetUniqueID() !=5) return;// check if the daughters come from a conversion
1114
1115 TParticle *Photon = TruePhotonCandidate->GetMCParticle(MCStack);
1116 if(Photon->GetPdgCode() != 22) return; // Mother is no Photon
1117
1118 // True Photon
1119
1120 Int_t labelGamma = TruePhotonCandidate->GetMCParticleLabel(MCStack);
1121
1122 if( labelGamma < MCStack->GetNprimary() ){
1123 if( fIsFromMBHeader ){
1124 fHistoTrueConvGammaPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1125 if (GammaIsNeutralMesonPiPlPiMiPiZeroDaughter(labelGamma)){
1126 fHistoTrueConvGammaFromNeutralMesonPt[fiCut]->Fill(TruePhotonCandidate->Pt());
1127 }
1128 }
1129 }
1130}
1131
1132//________________________________________________________________________
1133void AliAnalysisTaskNeutralMesonToPiPlPiMiPiZero::ProcessNeutralPionCandidatesPureConversions(){
1134 // Conversion Gammas
1135 if(fGoodConvGammas->GetEntries()>1){
1136 for(Int_t firstGammaIndex=0;firstGammaIndex<fGoodConvGammas->GetEntries()-1;firstGammaIndex++){
1137 AliAODConversionPhoton *gamma0=dynamic_cast<AliAODConversionPhoton*>(fGoodConvGammas->At(firstGammaIndex));
1138 if (gamma0==NULL) continue;
1139 for(Int_t secondGammaIndex=firstGammaIndex+1;secondGammaIndex<fGoodConvGammas->GetEntries();secondGammaIndex++){
1140 AliAODConversionPhoton *gamma1=dynamic_cast<AliAODConversionPhoton*>(fGoodConvGammas->At(secondGammaIndex));
1141 //Check for same Electron ID
1142 if (gamma1==NULL) continue;
1143 if(gamma0->GetTrackLabelPositive() == gamma1->GetTrackLabelPositive() ||
1144 gamma0->GetTrackLabelNegative() == gamma1->GetTrackLabelNegative() ||
1145 gamma0->GetTrackLabelNegative() == gamma1->GetTrackLabelPositive() ||
1146 gamma0->GetTrackLabelPositive() == gamma1->GetTrackLabelNegative() ) continue;
1147
1148 AliAODConversionMother *pi0cand = new AliAODConversionMother(gamma0,gamma1);
1149 pi0cand->SetLabels(firstGammaIndex,secondGammaIndex);
1150
1151 pi0cand->CalculateDistanceOfClossetApproachToPrimVtx(fInputEvent->GetPrimaryVertex());
1152 if((((AliConversionMesonCuts*)fNeutralPionMesonCutArray->At(fiCut))->MesonIsSelected(pi0cand,kTRUE,((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetEtaShift()))){
1153 fHistoGammaGammaInvMassPt[fiCut]->Fill(pi0cand->M(),pi0cand->Pt());
1154 if (pi0cand->M() > ((AliConversionMesonCuts*)fNeutralPionMesonCutArray->At(fiCut))->GetSelectionLow() && pi0cand->M() < ((AliConversionMesonCuts*)fNeutralPionMesonCutArray->At(fiCut))->GetSelectionHigh()){
1155 fNeutralPionCandidates->Add(pi0cand);
1156// cout << "Pi0 candidate " << pi0cand->M() << "\t" << pi0cand->Pt() << endl;
1157 }
1158
1159 if(fIsMC){
1160 if(fInputEvent->IsA()==AliESDEvent::Class())
1161 ProcessTrueNeutralPionCandidatesPureConversions(pi0cand,gamma0,gamma1);
1162 if(fInputEvent->IsA()==AliAODEvent::Class())
1163 ProcessTrueNeutralPionCandidatesPureConversionsAOD(pi0cand,gamma0,gamma1);
1164 }
1165 }
1166// delete pi0cand;
1167// pi0cand=0x0;
1168 }
1169 }
1170 }
1171}
1172
1173
1174//________________________________________________________________________
1175void AliAnalysisTaskNeutralMesonToPiPlPiMiPiZero::ProcessNeutralPionCandidatesPureCalo(){
1176
1177 // Conversion Gammas
1178 if(fClusterCandidates->GetEntries()>0){
1179
1180 // vertex
1181 Double_t vertex[3] = {0};
1182 InputEvent()->GetPrimaryVertex()->GetXYZ(vertex);
1183
1184 for(Int_t firstGammaIndex=0;firstGammaIndex<fClusterCandidates->GetEntries();firstGammaIndex++){
1185 AliAODConversionPhoton *gamma0=dynamic_cast<AliAODConversionPhoton*>(fClusterCandidates->At(firstGammaIndex));
1186 if (gamma0==NULL) continue;
1187
1188 for(Int_t secondGammaIndex=0;secondGammaIndex<fClusterCandidates->GetEntries();secondGammaIndex++){
1189 if (firstGammaIndex == secondGammaIndex) continue;
1190 AliAODConversionPhoton *gamma1=dynamic_cast<AliAODConversionPhoton*>(fClusterCandidates->At(secondGammaIndex));
1191 if (gamma1==NULL) continue;
1192
1193 AliAODConversionMother *pi0cand = new AliAODConversionMother(gamma0,gamma1);
1194 pi0cand->SetLabels(firstGammaIndex,secondGammaIndex);
1195
1196 if((((AliConversionMesonCuts*)fNeutralPionMesonCutArray->At(fiCut))->MesonIsSelected(pi0cand,kTRUE,((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetEtaShift()))){
1197 fHistoGammaGammaInvMassPt[fiCut]->Fill(pi0cand->M(),pi0cand->Pt());
1198 if (pi0cand->M() > ((AliConversionMesonCuts*)fNeutralPionMesonCutArray->At(fiCut))->GetSelectionLow() && pi0cand->M() < ((AliConversionMesonCuts*)fNeutralPionMesonCutArray->At(fiCut))->GetSelectionHigh()){
1199 fNeutralPionCandidates->Add(pi0cand);
1200// cout << "Pi0 candidate " << pi0cand->M() << "\t" << pi0cand->Pt() << endl;
1201 }
1202
1203 if(fIsMC){
1204// if(fInputEvent->IsA()==AliESDEvent::Class())
1205 ProcessTrueNeutralPionCandidatesPureCalo(pi0cand,gamma0,gamma1);
1206// if(fInputEvent->IsA()==AliAODEvent::Class())
1207// ProcessTrueNeutralPionCandidatesPureConversionsAOD(pi0cand,gamma0,gamma1);
1208 }
1209 }
1210 }
1211 }
1212 }
1213}
1214
1215//______________________________________________________________________
1216void AliAnalysisTaskNeutralMesonToPiPlPiMiPiZero::ProcessTrueNeutralPionCandidatesPureCalo( AliAODConversionMother *Pi0Candidate, AliAODConversionPhoton *TrueGammaCandidate0, AliAODConversionPhoton *TrueGammaCandidate1)
1217{
1218 // Process True Mesons
1219 AliStack *MCStack = fMCEvent->Stack();
1220
1221 Bool_t isTruePi0 = kFALSE;
1222 Int_t gamma0MCLabel = TrueGammaCandidate0->GetCaloPhotonMCLabel(0); // get most probable MC label
1223 Int_t gamma0MotherLabel = -1;
1224 Int_t motherRealLabel = -1;
1225
1226 if(gamma0MCLabel != -1){ // Gamma is Combinatorial; MC Particles don't belong to the same Mother
1227 TParticle * gammaMC0 = (TParticle*)MCStack->Particle(gamma0MCLabel);
1228 if (TrueGammaCandidate0->IsLargestComponentPhoton() || TrueGammaCandidate0->IsLargestComponentElectron()){ // largest component is electro magnetic
1229 // get mother of interest (pi0 or eta)
1230 if (TrueGammaCandidate0->IsLargestComponentPhoton()){ // for photons its the direct mother
1231 gamma0MotherLabel=gammaMC0->GetMother(0);
1232 motherRealLabel=gammaMC0->GetFirstMother();
1233 } else if (TrueGammaCandidate0->IsLargestComponentElectron()){ // for electrons its either the direct mother or for conversions the grandmother
1234 if (TrueGammaCandidate0->IsConversion()){
1235 gamma0MotherLabel=MCStack->Particle(gammaMC0->GetMother(0))->GetMother(0);
1236 motherRealLabel=MCStack->Particle(gammaMC0->GetMother(0))->GetMother(0);
1237 } else {
1238 gamma0MotherLabel=gammaMC0->GetMother(0);
1239 motherRealLabel=gammaMC0->GetMother(0);
1240 }
1241 }
1242 }
1243 }
1244
1245 if (!TrueGammaCandidate1->GetIsCaloPhoton()) AliFatal("CaloPhotonFlag has not been set. Aborting");
1246
1247 Int_t gamma1MCLabel = TrueGammaCandidate1->GetCaloPhotonMCLabel(0); // get most probable MC label
1248 Int_t gamma1MotherLabel = -1;
1249 // check if
1250 if(gamma1MCLabel != -1){ // Gamma is Combinatorial; MC Particles don't belong to the same Mother
1251 // Daughters Gamma 1
1252 TParticle * gammaMC1 = (TParticle*)MCStack->Particle(gamma1MCLabel);
1253 if (TrueGammaCandidate1->IsLargestComponentPhoton() || TrueGammaCandidate1->IsLargestComponentElectron()){ // largest component is electro magnetic
1254 // get mother of interest (pi0 or eta)
1255 if (TrueGammaCandidate1->IsLargestComponentPhoton()){ // for photons its the direct mother
1256 gamma1MotherLabel=gammaMC1->GetMother(0);
1257 } else if (TrueGammaCandidate1->IsLargestComponentElectron()){ // for electrons its either the direct mother or for conversions the grandmother
1258 if (TrueGammaCandidate1->IsConversion()) gamma1MotherLabel=MCStack->Particle(gammaMC1->GetMother(0))->GetMother(0);
1259 else gamma1MotherLabel=gammaMC1->GetMother(0);
1260 }
1261 }
1262 }
1263
1264 if(gamma0MotherLabel>=0 && gamma0MotherLabel==gamma1MotherLabel){
1265 if(((TParticle*)MCStack->Particle(gamma1MotherLabel))->GetPdgCode() == 111){
1266 isTruePi0=kTRUE;
1267 }
1268 }
1269
1270 if(isTruePi0){// True Pion
1271 Pi0Candidate->SetTrueMesonValue(1);
1272 Pi0Candidate->SetMCLabel(motherRealLabel);
1273 fHistoTrueMotherGammaGammaInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt());
1274 if( IsEtaPiPlPiMiPiZeroDaughter(motherRealLabel) ) {
1275 fHistoTrueMotherGammaGammaFromEtaInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt());
1276 }
1277 if( IsOmegaPiPlPiMiPiZeroDaughter(motherRealLabel) ) {
1278 fHistoTrueMotherGammaGammaFromOmegaInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt());
1279 }
1280 }
1281}
1282
1283
1284
1285//______________________________________________________________________
1286void AliAnalysisTaskNeutralMesonToPiPlPiMiPiZero::ProcessTrueNeutralPionCandidatesPureConversions(AliAODConversionMother *Pi0Candidate, AliAODConversionPhoton *TrueGammaCandidate0, AliAODConversionPhoton *TrueGammaCandidate1)
1287{
1288 // Process True Mesons
1289 AliStack *MCStack = fMCEvent->Stack();
1290 if(TrueGammaCandidate0->GetV0Index()<fInputEvent->GetNumberOfV0s()){
1291 Bool_t isTruePi0 = kFALSE;
1292 Bool_t isTruePi0Dalitz = kFALSE;
1293 Bool_t gamma0DalitzCand = kFALSE;
1294 Bool_t gamma1DalitzCand = kFALSE;
1295 Int_t gamma0MCLabel = TrueGammaCandidate0->GetMCParticleLabel(MCStack);
1296 Int_t gamma0MotherLabel = -1;
1297 Int_t motherRealLabel = -1;
1298 if(gamma0MCLabel != -1){ // Gamma is Combinatorial; MC Particles don't belong to the same Mother
1299 // Daughters Gamma 0
1300 TParticle * negativeMC = (TParticle*)TrueGammaCandidate0->GetNegativeMCDaughter(MCStack);
1301 TParticle * positiveMC = (TParticle*)TrueGammaCandidate0->GetPositiveMCDaughter(MCStack);
1302 TParticle * gammaMC0 = (TParticle*)MCStack->Particle(gamma0MCLabel);
1303 if(abs(negativeMC->GetPdgCode())==11 && abs(positiveMC->GetPdgCode())==11){ // Electrons ...
1304 if(negativeMC->GetUniqueID() == 5 && positiveMC->GetUniqueID() ==5){ // ... From Conversion ...
1305 if(gammaMC0->GetPdgCode() == 22){ // ... with Gamma Mother
1306 gamma0MotherLabel=gammaMC0->GetFirstMother();
1307 motherRealLabel=gammaMC0->GetFirstMother();
1308 }
1309 }
1310 if(gammaMC0->GetPdgCode() ==111){ // Dalitz candidate
1311 gamma0DalitzCand = kTRUE;
1312 gamma0MotherLabel=-111;
1313 motherRealLabel=gamma0MCLabel;
1314 }
1315 }
1316 }
1317 if(TrueGammaCandidate1->GetV0Index()<fInputEvent->GetNumberOfV0s()){
1318 Int_t gamma1MCLabel = TrueGammaCandidate1->GetMCParticleLabel(MCStack);
1319 Int_t gamma1MotherLabel = -1;
1320 if(gamma1MCLabel != -1){ // Gamma is Combinatorial; MC Particles don't belong to the same Mother
1321 // Daughters Gamma 1
1322 TParticle * negativeMC = (TParticle*)TrueGammaCandidate1->GetNegativeMCDaughter(MCStack);
1323 TParticle * positiveMC = (TParticle*)TrueGammaCandidate1->GetPositiveMCDaughter(MCStack);
1324 TParticle * gammaMC1 = (TParticle*)MCStack->Particle(gamma1MCLabel);
1325 if(abs(negativeMC->GetPdgCode())==11 && abs(positiveMC->GetPdgCode())==11){ // Electrons ...
1326 if(negativeMC->GetUniqueID() == 5 && positiveMC->GetUniqueID() ==5){ // ... From Conversion ...
1327 if(gammaMC1->GetPdgCode() == 22){ // ... with Gamma Mother
1328 gamma1MotherLabel=gammaMC1->GetFirstMother();
1329 }
1330 }
1331 if(gammaMC1->GetPdgCode() ==111 ){ // Dalitz candidate
1332 gamma1DalitzCand = kTRUE;
1333 gamma1MotherLabel=-111;
1334 }
1335 }
1336 }
1337 if(gamma0MotherLabel>=0 && gamma0MotherLabel==gamma1MotherLabel){
1338 if(((TParticle*)MCStack->Particle(gamma1MotherLabel))->GetPdgCode() == 111){
1339 isTruePi0=kTRUE;
1340 }
1341 }
1342
1343 //Identify Dalitz candidate
1344 if (gamma1DalitzCand || gamma0DalitzCand){
1345 if (gamma0DalitzCand && gamma0MCLabel >=0 && gamma0MCLabel==gamma1MotherLabel){
1346 if (gamma0MotherLabel == -111) isTruePi0Dalitz = kTRUE;
1347 }
1348 if (gamma1DalitzCand && gamma1MCLabel >=0 && gamma1MCLabel==gamma0MotherLabel){
1349 if (gamma1MotherLabel == -111) isTruePi0Dalitz = kTRUE;
1350 }
1351 }
1352
1353
1354 if(isTruePi0 || isTruePi0Dalitz){// True Pion
1355 Pi0Candidate->SetTrueMesonValue(1);
1356 Pi0Candidate->SetMCLabel(motherRealLabel);
1357 fHistoTrueMotherGammaGammaInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt());
1358 if( IsEtaPiPlPiMiPiZeroDaughter(motherRealLabel) ) {
1359 fHistoTrueMotherGammaGammaFromEtaInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt());
1360 }
1361 if( IsOmegaPiPlPiMiPiZeroDaughter(motherRealLabel) ) {
1362 fHistoTrueMotherGammaGammaFromOmegaInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt());
1363 }
1364 }
1365 }
1366 }
1367}
1368
1369//______________________________________________________________________
1370void AliAnalysisTaskNeutralMesonToPiPlPiMiPiZero::ProcessTrueNeutralPionCandidatesPureConversionsAOD(AliAODConversionMother *Pi0Candidate, AliAODConversionPhoton *TrueGammaCandidate0, AliAODConversionPhoton *TrueGammaCandidate1)
1371{
1372
1373 // Process True Mesons
1374 TClonesArray *AODMCTrackArray = dynamic_cast<TClonesArray*>(fInputEvent->FindListObject(AliAODMCParticle::StdBranchName()));
1375 Bool_t isTruePi0 = kFALSE;
1376 Bool_t isTruePi0Dalitz = kFALSE;
1377 Bool_t gamma0DalitzCand = kFALSE;
1378 Bool_t gamma1DalitzCand = kFALSE;
1379 Int_t motherRealLabel = -1;
1380
1381 if (AODMCTrackArray!=NULL && TrueGammaCandidate0 != NULL){
1382 AliAODMCParticle *positiveMC = static_cast<AliAODMCParticle*>(AODMCTrackArray->At(TrueGammaCandidate0->GetMCLabelPositive()));
1383 AliAODMCParticle *negativeMC = static_cast<AliAODMCParticle*>(AODMCTrackArray->At(TrueGammaCandidate0->GetMCLabelNegative()));
1384
1385 Int_t gamma0MCLabel = -1;
1386 Int_t gamma0MotherLabel = -1;
1387 if(!positiveMC||!negativeMC)
1388 return;
1389
1390 if(positiveMC->GetMother()>-1&&(negativeMC->GetMother() == positiveMC->GetMother())){
1391 gamma0MCLabel = positiveMC->GetMother();
1392 }
1393
1394 if(gamma0MCLabel != -1){ // Gamma is Combinatorial; MC Particles don't belong to the same Mother
1395 // Daughters Gamma 0
1396 AliAODMCParticle * gammaMC0 = static_cast<AliAODMCParticle*>(AODMCTrackArray->At(gamma0MCLabel));
1397 if(abs(negativeMC->GetPdgCode())==11 && abs(positiveMC->GetPdgCode())==11){ // Electrons ...
1398 if(((positiveMC->GetMCProcessCode())) == 5 && ((negativeMC->GetMCProcessCode())) == 5){ // ... From Conversion ...
1399 if(gammaMC0->GetPdgCode() == 22){ // ... with Gamma Mother
1400 gamma0MotherLabel=gammaMC0->GetMother();
1401 motherRealLabel=gammaMC0->GetMother();
1402 }
1403 }
1404 if(gammaMC0->GetPdgCode() ==111){ // Dalitz candidate
1405 gamma0DalitzCand = kTRUE;
1406 gamma0MotherLabel=-111;
1407 motherRealLabel=gamma0MCLabel;
1408 }
1409 }
1410 }
1411 positiveMC = static_cast<AliAODMCParticle*>(AODMCTrackArray->At(TrueGammaCandidate1->GetMCLabelPositive()));
1412 negativeMC = static_cast<AliAODMCParticle*>(AODMCTrackArray->At(TrueGammaCandidate1->GetMCLabelNegative()));
1413
1414 Int_t gamma1MCLabel = -1;
1415 Int_t gamma1MotherLabel = -1;
1416 if(!positiveMC||!negativeMC)
1417 return;
1418
1419 if(positiveMC->GetMother()>-1&&(negativeMC->GetMother() == positiveMC->GetMother())){
1420 gamma1MCLabel = positiveMC->GetMother();
1421 }
1422 if(gamma1MCLabel != -1){ // Gamma is Combinatorial; MC Particles don't belong to the same Mother
1423 // Daughters Gamma 1
1424 AliAODMCParticle * gammaMC1 = static_cast<AliAODMCParticle*>(AODMCTrackArray->At(gamma1MCLabel));
1425 if(abs(negativeMC->GetPdgCode())==11 && abs(positiveMC->GetPdgCode())==11){ // Electrons ...
1426 if(((positiveMC->GetMCProcessCode())) == 5 && ((negativeMC->GetMCProcessCode())) == 5){ // ... From Conversion ...
1427 if(gammaMC1->GetPdgCode() == 22){ // ... with Gamma Mother
1428 gamma1MotherLabel=gammaMC1->GetMother();
1429 }
1430 }
1431 if(gammaMC1->GetPdgCode() ==111 ){ // Dalitz candidate
1432 gamma1DalitzCand = kTRUE;
1433 gamma1MotherLabel=-111;
1434 }
1435 }
1436 }
1437 if(gamma0MotherLabel>=0 && gamma0MotherLabel==gamma1MotherLabel){
1438 if(static_cast<AliAODMCParticle*>(AODMCTrackArray->At(gamma1MotherLabel))->GetPdgCode() == 111){
1439 isTruePi0=kTRUE;
1440 }
1441 }
1442
1443 //Identify Dalitz candidate
1444 if (gamma1DalitzCand || gamma0DalitzCand){
1445 if (gamma0DalitzCand && gamma0MCLabel >=0 && gamma0MCLabel==gamma1MotherLabel){
1446 if (gamma0MotherLabel == -111) isTruePi0Dalitz = kTRUE;
1447 }
1448 if (gamma1DalitzCand && gamma1MCLabel >=0 && gamma1MCLabel==gamma0MotherLabel){
1449 if (gamma1MotherLabel == -111) isTruePi0Dalitz = kTRUE;
1450 }
1451 }
1452
1453 if(isTruePi0 || isTruePi0Dalitz){// True Pion
1454 Pi0Candidate->SetTrueMesonValue(1);
1455 Pi0Candidate->SetMCLabel(motherRealLabel);
1456 fHistoTrueMotherGammaGammaInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt());
1457 if( IsEtaPiPlPiMiPiZeroDaughter(motherRealLabel) ) {
1458 fHistoTrueMotherGammaGammaFromEtaInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt());
1459 }
1460 if( IsOmegaPiPlPiMiPiZeroDaughter(motherRealLabel) ) {
1461 fHistoTrueMotherGammaGammaFromOmegaInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt());
1462 }
1463 }
1464 }
1465 return;
1466}
1467
1468
1469//________________________________________________________________________
1470void AliAnalysisTaskNeutralMesonToPiPlPiMiPiZero::ProcessNeutralPionCandidatesMixedConvCalo(){
1471
1472 // Conversion Gammas
1473 if(fGoodConvGammas->GetEntries()>0){
1474 // vertex
1475 Double_t vertex[3] = {0};
1476 InputEvent()->GetPrimaryVertex()->GetXYZ(vertex);
1477
1478 for(Int_t firstGammaIndex=0;firstGammaIndex<fGoodConvGammas->GetEntries();firstGammaIndex++){
1479 AliAODConversionPhoton *gamma0=dynamic_cast<AliAODConversionPhoton*>(fGoodConvGammas->At(firstGammaIndex));
1480 if (gamma0==NULL) continue;
1481
1482 for(Int_t secondGammaIndex=0;secondGammaIndex<fClusterCandidates->GetEntries();secondGammaIndex++){
1483 Bool_t matched = kFALSE;
1484 AliAODConversionPhoton *gamma1=dynamic_cast<AliAODConversionPhoton*>(fClusterCandidates->At(secondGammaIndex));
1485 if (gamma1==NULL) continue;
1486
1487 if (gamma1->GetIsCaloPhoton()){
1488 AliVCluster* cluster = fInputEvent->GetCaloCluster(gamma1->GetCaloClusterRef());
1489 matched = ((AliCaloPhotonCuts*)fClusterCutArray->At(fiCut))->MatchConvPhotonToCluster(gamma0,cluster, fInputEvent );
1490 }
1491
1492 AliAODConversionMother *pi0cand = new AliAODConversionMother(gamma0,gamma1);
1493 pi0cand->SetLabels(firstGammaIndex,secondGammaIndex);
1494
1495 if((((AliConversionMesonCuts*)fNeutralPionMesonCutArray->At(fiCut))->MesonIsSelected(pi0cand,kTRUE,((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetEtaShift()))){
1496 if (!matched){
1497 fHistoGammaGammaInvMassPt[fiCut]->Fill(pi0cand->M(),pi0cand->Pt());
1498 if (pi0cand->M() > ((AliConversionMesonCuts*)fNeutralPionMesonCutArray->At(fiCut))->GetSelectionLow() && pi0cand->M() < ((AliConversionMesonCuts*)fNeutralPionMesonCutArray->At(fiCut))->GetSelectionHigh()){
1499 fNeutralPionCandidates->Add(pi0cand);
1500 // cout << "Pi0 candidate " << pi0cand->M() << "\t" << pi0cand->Pt() << endl;
1501 }
1502
1503 if(fIsMC){
1504// if(fInputEvent->IsA()==AliESDEvent::Class())
1505 ProcessTrueNeutralPionCandidatesMixedConvCalo(pi0cand,gamma0,gamma1);
1506// if(fInputEvent->IsA()==AliAODEvent::Class())
1507// ProcessTrueMesonCandidatesAOD(pi0cand,gamma0,gamma1, matched);
1508 }
1509 }
1510 }
1511 }
1512 }
1513 }
1514}
1515
1516//______________________________________________________________________
1517void AliAnalysisTaskNeutralMesonToPiPlPiMiPiZero::ProcessTrueNeutralPionCandidatesMixedConvCalo( AliAODConversionMother *Pi0Candidate, AliAODConversionPhoton *TrueGammaCandidate0, AliAODConversionPhoton *TrueGammaCandidate1)
1518{
1519 // Process True Mesons
1520 AliStack *MCStack = fMCEvent->Stack();
1521 if(TrueGammaCandidate0->GetV0Index()<fInputEvent->GetNumberOfV0s()){
1522 Bool_t isTruePi0 = kFALSE;
1523 Bool_t isTruePi0Dalitz = kFALSE;
1524 Bool_t gamma0DalitzCand = kFALSE;
1525
1526 Int_t gamma0MCLabel = TrueGammaCandidate0->GetMCParticleLabel(MCStack);
1527 Int_t gamma0MotherLabel = -1;
1528 Int_t motherRealLabel = -1;
1529 if(gamma0MCLabel != -1){ // Gamma is Combinatorial; MC Particles don't belong to the same Mother
1530 // Daughters Gamma 0
1531 TParticle * negativeMC = (TParticle*)TrueGammaCandidate0->GetNegativeMCDaughter(MCStack);
1532 TParticle * positiveMC = (TParticle*)TrueGammaCandidate0->GetPositiveMCDaughter(MCStack);
1533 TParticle * gammaMC0 = (TParticle*)MCStack->Particle(gamma0MCLabel);
1534 if(abs(negativeMC->GetPdgCode())==11 && abs(positiveMC->GetPdgCode())==11){ // Electrons ...
1535 if(negativeMC->GetUniqueID() == 5 && positiveMC->GetUniqueID() ==5){ // ... From Conversion ...
1536 if(gammaMC0->GetPdgCode() == 22){ // ... with Gamma Mother
1537 gamma0MotherLabel=gammaMC0->GetFirstMother();
1538 motherRealLabel=gammaMC0->GetFirstMother();
1539 }
1540 }
1541 if(gammaMC0->GetPdgCode() ==111){ // Dalitz candidate
1542 gamma0DalitzCand = kTRUE;
1543 gamma0MotherLabel=-111;
1544 motherRealLabel=gamma0MCLabel;
1545 }
1546
1547 }
1548 }
1549
1550 if (!TrueGammaCandidate1->GetIsCaloPhoton()) AliFatal("CaloPhotonFlag has not been set. Aborting");
1551
1552 Int_t gamma1MCLabel = TrueGammaCandidate1->GetCaloPhotonMCLabel(0); // get most probable MC label
1553 Int_t gamma1MotherLabel = -1;
1554 // check if
1555
1556 if(gamma1MCLabel != -1){ // Gamma is Combinatorial; MC Particles don't belong to the same Mother
1557 // Daughters Gamma 1
1558 TParticle * gammaMC1 = (TParticle*)MCStack->Particle(gamma1MCLabel);
1559 if (TrueGammaCandidate1->IsLargestComponentPhoton() || TrueGammaCandidate1->IsLargestComponentElectron()){ // largest component is electro magnetic
1560 // get mother of interest (pi0 or eta)
1561 if (TrueGammaCandidate1->IsLargestComponentPhoton()){ // for photons its the direct mother
1562 gamma1MotherLabel=gammaMC1->GetMother(0);
1563 } else if (TrueGammaCandidate1->IsLargestComponentElectron()){ // for electrons its either the direct mother or for conversions the grandmother
1564 if (TrueGammaCandidate1->IsConversion()) gamma1MotherLabel=MCStack->Particle(gammaMC1->GetMother(0))->GetMother(0);
1565 else gamma1MotherLabel=gammaMC1->GetMother(0);
1566 }
1567 }
1568 }
1569
1570 if(gamma0MotherLabel>=0 && gamma0MotherLabel==gamma1MotherLabel){
1571 if(((TParticle*)MCStack->Particle(gamma1MotherLabel))->GetPdgCode() == 111){
1572 isTruePi0=kTRUE;
1573 }
1574 }
1575
1576 if (gamma0DalitzCand ){
1577 if (gamma0DalitzCand && gamma0MCLabel >=0 && gamma0MCLabel==gamma1MotherLabel){
1578 if (gamma0MotherLabel == -111) isTruePi0Dalitz = kTRUE;
1579 }
1580 }
1581
1582 if(isTruePi0 || isTruePi0Dalitz ){
1583 Pi0Candidate->SetTrueMesonValue(1);
1584 Pi0Candidate->SetMCLabel(motherRealLabel);
1585 fHistoTrueMotherGammaGammaInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt());
1586 if( IsEtaPiPlPiMiPiZeroDaughter(motherRealLabel) ) {
1587 fHistoTrueMotherGammaGammaFromEtaInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt());
1588 }
1589 if( IsOmegaPiPlPiMiPiZeroDaughter(motherRealLabel) ) {
1590 fHistoTrueMotherGammaGammaFromOmegaInvMassPt[fiCut]->Fill(Pi0Candidate->M(),Pi0Candidate->Pt());
1591 }
1592 }
1593 }
1594}
1595
1596
1597
1598//________________________________________________________________________
1599void AliAnalysisTaskNeutralMesonToPiPlPiMiPiZero::ProcessPionCandidates(){
1600
1601 Double_t magField = fInputEvent->GetMagneticField();
1602 if( magField < 0.0 ){
1603 magField = 1.0;
1604 } else {
1605 magField = -1.0;
1606 }
1607
1608 vector<Int_t> lGoodNegPionIndexPrev(0);
1609 vector<Int_t> lGoodPosPionIndexPrev(0);
1610
1611 for(UInt_t i = 0; i < fSelectorNegPionIndex.size(); i++){
1612 AliESDtrack* negPionCandidate = fESDEvent->GetTrack(fSelectorNegPionIndex[i]);
1613 if(! ((AliPrimaryPionCuts*)fPionCutArray->At(fiCut))->PionIsSelected(negPionCandidate) ) continue;
1614 lGoodNegPionIndexPrev.push_back( fSelectorNegPionIndex[i] );
1615 fHistoNegPionPt[fiCut]->Fill(negPionCandidate->Pt());
1616 fHistoNegPionPhi[fiCut]->Fill(negPionCandidate->Phi());
1617 if( fMCEvent ) {
1618 Int_t labelNegPion = TMath::Abs( negPionCandidate->GetLabel() );
1619 if( labelNegPion < fMCStack->GetNtrack() ){
1620 TParticle* negPion = fMCStack->Particle(labelNegPion);
1621 if( negPion->GetPdgCode() == -211 ){
1622 if( labelNegPion < fMCStack->GetNprimary() ){
1623 fHistoTrueNegPionPt[fiCut]->Fill(negPionCandidate->Pt()); //primary negPion
1624 }
1625 if( IsEtaPiPlPiMiPiZeroDaughter(labelNegPion) || IsOmegaPiPlPiMiPiZeroDaughter(labelNegPion) ) {
1626 if( labelNegPion < fMCStack->GetNprimary() ) {
1627 fHistoTrueNegPionFromNeutralMesonPt[fiCut]->Fill(negPionCandidate->Pt());
1628 }
1629 }
1630 }
1631 }
1632 }
1633 }
1634
1635 for(UInt_t i = 0; i < fSelectorPosPionIndex.size(); i++){
1636 AliESDtrack* posPionCandidate = fESDEvent->GetTrack( fSelectorPosPionIndex[i] );
1637 if(! ((AliPrimaryPionCuts*)fPionCutArray->At(fiCut))->PionIsSelected(posPionCandidate) ) continue;
1638 lGoodPosPionIndexPrev.push_back( fSelectorPosPionIndex[i] );
1639 fHistoPosPionPt[fiCut]->Fill( posPionCandidate->Pt() );
1640 fHistoPosPionPhi[fiCut]->Fill( posPionCandidate->Phi() );
1641
1642 if( fMCEvent ) {
1643 Int_t labelPosPion = TMath::Abs( posPionCandidate->GetLabel() );
1644 if( labelPosPion < fMCStack->GetNtrack() ) {
1645 TParticle* posPion = fMCStack->Particle(labelPosPion);
1646 if( posPion->GetPdgCode() == 211 ){
1647 if( labelPosPion < fMCStack->GetNprimary() ){
1648 fHistoTruePosPionPt[fiCut]->Fill(posPionCandidate->Pt());
1649 }
1650 if( IsEtaPiPlPiMiPiZeroDaughter(labelPosPion) || IsOmegaPiPlPiMiPiZeroDaughter(labelPosPion) ) {
1651 if( labelPosPion < fMCStack->GetNprimary() ){
1652 fHistoTruePosPionFromNeutralMesonPt[fiCut]->Fill(posPionCandidate->Pt());
1653 }
1654 }
1655 }
1656 }
1657 }
1658 }
1659
1660
1661 for(UInt_t i = 0; i < lGoodNegPionIndexPrev.size(); i++){
1662
1663 AliESDtrack *negPionCandidate = fESDEvent->GetTrack(lGoodNegPionIndexPrev[i]);
1664 AliKFParticle negPionCandidateKF( *negPionCandidate->GetConstrainedParam(), 211 );
1665
1666 for(UInt_t j = 0; j < lGoodPosPionIndexPrev.size(); j++){
1667
1668 AliESDtrack *posPionCandidate = fESDEvent->GetTrack(lGoodPosPionIndexPrev[j]);
1669 AliKFParticle posPionCandidateKF( *posPionCandidate->GetConstrainedParam(), 211 );
1670
1671 AliKFConversionPhoton* virtualPhoton = NULL;
1672 virtualPhoton = new AliKFConversionPhoton(negPionCandidateKF,posPionCandidateKF);
1673 AliKFVertex primaryVertexImproved(*fInputEvent->GetPrimaryVertex());
1674// primaryVertexImproved+=*virtualPhoton;
1675 virtualPhoton->SetProductionVertex(primaryVertexImproved);
1676 virtualPhoton->SetTrackLabels( lGoodPosPionIndexPrev[j], lGoodNegPionIndexPrev[i]);
1677
1678 Int_t labeln=0;
1679 Int_t labelp=0;
1680 Int_t motherlabelp = 0;
1681 Int_t motherlabeln = 0;
1682 TParticle *fNegativeMCParticle =NULL;
1683 TParticle *fPositiveMCParticle =NULL;
1684 if( fMCEvent ) {
1685 labeln=TMath::Abs(negPionCandidate->GetLabel());
1686 labelp=TMath::Abs(posPionCandidate->GetLabel());
1687 fNegativeMCParticle = fMCStack->Particle(labeln);
1688 fPositiveMCParticle = fMCStack->Particle(labelp);
1689 motherlabeln = fNegativeMCParticle->GetMother(0);
1690 motherlabelp = fPositiveMCParticle->GetMother(0);
1691 if( fPositiveMCParticle && fNegativeMCParticle) {
1692 virtualPhoton->SetMCLabelPositive(labelp);
1693 virtualPhoton->SetMCLabelNegative(labeln);
1694 }
1695 }
1696
1697 AliAODConversionPhoton *vParticle = new AliAODConversionPhoton(virtualPhoton); //To apply mass 2 pion mass cut
1698 if (((AliPrimaryPionCuts*)fPionCutArray->At(fiCut))->DoMassCut()){
1699 if (vParticle->GetMass() < ((AliPrimaryPionCuts*)fPionCutArray->At(fiCut))->GetMassCut()){
1700 fGoodVirtualParticles->Add( vParticle );
1701 fHistoPionPionInvMassPt[fiCut]->Fill( vParticle->GetMass(),vParticle->Pt());
1702 }
1703 } else {
1704 fGoodVirtualParticles->Add( vParticle );
1705 fHistoPionPionInvMassPt[fiCut]->Fill( vParticle->GetMass(),vParticle->Pt());
1706 }
1707
1708 Double_t clsToFPos = -1.0;
1709 Double_t clsToFNeg = -1.0;
1710
1711 Float_t dcaToVertexXYPos = -1.0;
1712 Float_t dcaToVertexZPos = -1.0;
1713 Float_t dcaToVertexXYNeg = -1.0;
1714 Float_t dcaToVertexZNeg = -1.0;
1715
1716
1717
1718 if ( fDoMesonQA ) {
1719 clsToFPos = ((AliPrimaryPionCuts*)fPionCutArray->At(fiCut))->GetNFindableClustersTPC(posPionCandidate);
1720 clsToFNeg = ((AliPrimaryPionCuts*)fPionCutArray->At(fiCut))->GetNFindableClustersTPC(negPionCandidate);
1721
1722 Float_t bPos[2];
1723 Float_t bCovPos[3];
1724 posPionCandidate->GetImpactParameters(bPos,bCovPos);
1725 if (bCovPos[0]<=0 || bCovPos[2]<=0) {
1726 AliDebug(1, "Estimated b resolution lower or equal zero!");
1727 bCovPos[0]=0; bCovPos[2]=0;
1728 }
1729
1730 Float_t bNeg[2];
1731 Float_t bCovNeg[3];
1732 posPionCandidate->GetImpactParameters(bNeg,bCovNeg);
1733 if (bCovNeg[0]<=0 || bCovNeg[2]<=0) {
1734 AliDebug(1, "Estimated b resolution lower or equal zero!");
1735 bCovNeg[0]=0; bCovNeg[2]=0;
1736 }
1737
1738 dcaToVertexXYPos = bPos[0];
1739 dcaToVertexZPos = bPos[1];
1740 dcaToVertexXYNeg = bNeg[0];
1741 dcaToVertexZNeg = bNeg[1];
1742
1743 fHistoNegPionEta[fiCut]->Fill( negPionCandidate->Eta() );
1744 fHistoPosPionEta[fiCut]->Fill( posPionCandidate->Eta() );
1745
1746 fHistoNegPionClsTPC[fiCut]->Fill(clsToFNeg,negPionCandidate->Pt());
1747 fHistoPosPionClsTPC[fiCut]->Fill(clsToFPos,posPionCandidate->Pt());
1748
1749 fHistoPionDCAxy[fiCut]->Fill( dcaToVertexXYNeg, negPionCandidate->Pt() );
1750 fHistoPionDCAz[fiCut]->Fill( dcaToVertexZNeg, negPionCandidate->Pt() );
1751 fHistoPionDCAxy[fiCut]->Fill( dcaToVertexXYPos, posPionCandidate->Pt() );
1752 fHistoPionDCAz[fiCut]->Fill( dcaToVertexZPos, posPionCandidate->Pt() );
1753
1754 fHistoPionTPCdEdxNSigma[fiCut]->Fill( posPionCandidate->P(),((AliPrimaryPionCuts*)fPionCutArray->At(fiCut))->GetPIDResponse()->NumberOfSigmasTPC(posPionCandidate, AliPID::kPion) );
1755 fHistoPionTPCdEdxNSigma[fiCut]->Fill( negPionCandidate->P(),((AliPrimaryPionCuts*)fPionCutArray->At(fiCut))->GetPIDResponse()->NumberOfSigmasTPC(negPionCandidate, AliPID::kPion) );
1756
1757 fHistoPionTPCdEdx[fiCut]->Fill( posPionCandidate->P(), TMath::Abs(posPionCandidate->GetTPCsignal()));
1758 fHistoPionTPCdEdx[fiCut]->Fill( negPionCandidate->P(), TMath::Abs(negPionCandidate->GetTPCsignal()));
1759 }
1760
1761 if (fMCEvent){
1762 if (fPositiveMCParticle && fNegativeMCParticle ) {
1763 if (((AliPrimaryPionCuts*)fPionCutArray->At(fiCut))->DoMassCut()){
1764 if (vParticle->GetMass() < ((AliPrimaryPionCuts*)fPionCutArray->At(fiCut))->GetMassCut()){
1765 if(TMath::Abs(fNegativeMCParticle->GetPdgCode())==211 && TMath::Abs(fPositiveMCParticle->GetPdgCode())==211){ // Pions ...
1766 fHistoTruePionPionInvMassPt[fiCut]->Fill(vParticle->GetMass(),vParticle->Pt());
1767 if (motherlabeln == motherlabelp){
1768 fHistoTruePionPionFromSameMotherInvMassPt[fiCut]->Fill(vParticle->GetMass(),vParticle->Pt());
1769 if( IsEtaPiPlPiMiPiZeroDaughter(labeln) ) { //|| IsOmegaPiPlPiMiPiZeroDaughter(labeln)
1770 fHistoTruePionPionFromEtaInvMassPt[fiCut]->Fill(vParticle->GetMass(),vParticle->Pt());
1771 }
1772 if( IsOmegaPiPlPiMiPiZeroDaughter(labeln) ) { //||
1773 fHistoTruePionPionFromOmegaInvMassPt[fiCut]->Fill(vParticle->GetMass(),vParticle->Pt());
1774 }
1775 }
1776 }
1777 }
1778 } else {
1779 if(TMath::Abs(fNegativeMCParticle->GetPdgCode())==211 && TMath::Abs(fPositiveMCParticle->GetPdgCode())==211){ // Pions ...
1780 fHistoTruePionPionInvMassPt[fiCut]->Fill(vParticle->GetMass(),vParticle->Pt());
1781 if (motherlabeln == motherlabelp){
1782 fHistoTruePionPionFromSameMotherInvMassPt[fiCut]->Fill(vParticle->GetMass(),vParticle->Pt());
1783 if( IsEtaPiPlPiMiPiZeroDaughter(labeln) ) { //|| IsOmegaPiPlPiMiPiZeroDaughter(labeln)
1784 fHistoTruePionPionFromEtaInvMassPt[fiCut]->Fill(vParticle->GetMass(),vParticle->Pt());
1785 }
1786 if( IsOmegaPiPlPiMiPiZeroDaughter(labeln) ) { //||
1787 fHistoTruePionPionFromOmegaInvMassPt[fiCut]->Fill(vParticle->GetMass(),vParticle->Pt());
1788 }
1789 }
1790 }
1791 }
1792 }
1793 }
1794
1795 delete virtualPhoton;
1796 virtualPhoton=NULL;
1797
1798 }
1799 }
1800}
1801
1802//_____________________________________________________________________________
1803void AliAnalysisTaskNeutralMesonToPiPlPiMiPiZero::ProcessMCParticles(){
1804
1805 // Loop over all primary MC particle
1806
1807 for(Int_t i = 0; i < fMCStack->GetNprimary(); i++) {
1808
1809 TParticle* particle = (TParticle *)fMCStack->Particle(i);
1810 if (!particle) continue;
1811
1812 Int_t isMCFromMBHeader = -1;
1813 if(((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetSignalRejection() != 0){
1814 isMCFromMBHeader
1815 = ((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsParticleFromBGEvent(i, fMCStack, fInputEvent);
1816 if(isMCFromMBHeader == 0 && ((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetSignalRejection() != 3) continue;
1817 }
1818
1819 if(((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsParticleFromBGEvent(i, fMCStack,fInputEvent)){
1820 // find MC photons
1821 if (fNeutralPionMode < 2){
1822 if(((AliConversionPhotonCuts*)fGammaCutArray->At(fiCut))->PhotonIsSelectedMC(particle,fMCStack,kFALSE)){
1823 fHistoMCAllGammaPt[fiCut]->Fill(particle->Pt()); // All MC Gamma
1824 if(particle->GetMother(0) >-1){
1825 if (fMCStack->Particle(particle->GetMother(0))->GetPdgCode() ==111){
1826 if (fMCStack->Particle(particle->GetMother(0))->GetMother(0) > -1){
1827 if ( fMCStack->Particle((fMCStack->Particle(particle->GetMother(0)))->GetMother(0))->GetPdgCode() == 221 ||
1828 fMCStack->Particle((fMCStack->Particle(particle->GetMother(0)))->GetMother(0))->GetPdgCode() == 223 ){
1829 if ( fMCStack->Particle(particle->GetMother(0))->GetNDaughters()==3 )
1830 fHistoMCGammaFromNeutralMesonPt[fiCut]->Fill(particle->Pt()); // All photons from eta or omega via pi0
1831 }
1832 }
1833 }
1834 }
1835 }
1836 } else if (fNeutralPionMode == 2){
1837 if(((AliCaloPhotonCuts*)fClusterCutArray->At(fiCut))->ClusterIsSelectedMC(particle,fMCStack)){
1838 fHistoMCAllGammaPt[fiCut]->Fill(particle->Pt()); // All MC Gamma
1839 if(particle->GetMother(0) >-1){
1840 if (fMCStack->Particle(particle->GetMother(0))->GetPdgCode() ==111){
1841 if (fMCStack->Particle(particle->GetMother(0))->GetMother(0) > -1){
1842 if ( fMCStack->Particle((fMCStack->Particle(particle->GetMother(0)))->GetMother(0))->GetPdgCode() == 221 ||
1843 fMCStack->Particle((fMCStack->Particle(particle->GetMother(0)))->GetMother(0))->GetPdgCode() == 223 ){
1844 if ( fMCStack->Particle(particle->GetMother(0))->GetNDaughters()==3 )
1845 fHistoMCGammaFromNeutralMesonPt[fiCut]->Fill(particle->Pt()); // All photons from eta or omega via pi0
1846 }
1847 }
1848 }
1849 }
1850 }
1851 }
1852 if (fNeutralPionMode < 2){
1853 if (((AliConversionPhotonCuts*)fGammaCutArray->At(fiCut))->PhotonIsSelectedMC(particle,fMCStack,kTRUE)){
1854 fHistoMCConvGammaPt[fiCut]->Fill(particle->Pt());
1855 } // Converted MC Gamma
1856 }
1857 if(((AliPrimaryPionCuts*)fPionCutArray->At(fiCut))->PionIsSelectedMC(i,fMCStack)){
1858 if( particle->GetPdgCode() == 211){
1859 fHistoMCAllPosPionsPt[fiCut]->Fill(particle->Pt()); // All pos pions
1860 if(particle->GetMother(0) >-1){
1861 if (fMCStack->Particle(particle->GetMother(0))->GetPdgCode() ==221 || fMCStack->Particle(particle->GetMother(0))->GetPdgCode() ==223)
1862 fHistoMCPosPionsFromNeutralMesonPt[fiCut]->Fill(particle->Pt()); // All pos from eta or omega
1863 }
1864 }
1865 if( particle->GetPdgCode() == -211){
1866 fHistoMCAllNegPionsPt[fiCut]->Fill(particle->Pt()); // All neg pions
1867 if(particle->GetMother(0) >-1){
1868 if (fMCStack->Particle(particle->GetMother(0))->GetPdgCode() ==221 || fMCStack->Particle(particle->GetMother(0))->GetPdgCode() ==223 )
1869 fHistoMCNegPionsFromNeutralMesonPt[fiCut]->Fill(particle->Pt()); // All pos from eta or omega
1870 }
1871 }
1872 }
1873
1874
1875 // \eta -> pi+ pi- \gamma
1876 Int_t labelNeutPion = -1;
1877 Int_t labelNegPion = -1;
1878 Int_t labelPosPion = -1;
1879
1880 if( ((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->MesonIsSelectedMCPiPlPiMiPiZero(particle,fMCStack,labelNegPion,labelPosPion,labelNeutPion,((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetEtaShift())){
1881 Float_t weighted= 1;
1882 if( ((AliPrimaryPionCuts*) fPionCutArray->At(fiCut))->DoWeights() ) {
1883 if(((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsParticleFromBGEvent(i, fMCStack,fInputEvent)){
1884 if (particle->Pt()>0.005){
1885 weighted= ((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetWeightForMeson(fV0Reader->GetPeriodName(),i, fMCStack,fInputEvent);
1886 }
1887 }
1888 }
1889 if(particle->GetPdgCode() == 221)fHistoMCEtaPiPlPiMiPiZeroPt[fiCut]->Fill(particle->Pt(), weighted); // All MC Eta in respective decay channel
1890 if(particle->GetPdgCode() == 223)fHistoMCOmegaPiPlPiMiPiZeroPt[fiCut]->Fill(particle->Pt(), weighted); // All MC Omega in respective decay channel
1891
1892 TParticle *neutPion = fMCStack->Particle(labelNeutPion);
1893 TParticle *gamma1 = fMCStack->Particle(neutPion->GetDaughter(0));
1894 TParticle *gamma2 = fMCStack->Particle(neutPion->GetDaughter(1));
1895 if (fNeutralPionMode < 2){
1896 if(
1897 ((AliConversionPhotonCuts*)fGammaCutArray->At(fiCut))->PhotonIsSelectedMC(gamma1,fMCStack,kFALSE) && // test first daugther of pi0
1898 ((AliConversionPhotonCuts*)fGammaCutArray->At(fiCut))->PhotonIsSelectedMC(gamma2,fMCStack,kFALSE) && // test second daughter of pi0
1899 ((AliPrimaryPionCuts*)fPionCutArray->At(fiCut))->PionIsSelectedMC(labelNegPion,fMCStack) && // test negative pion
1900 ((AliPrimaryPionCuts*)fPionCutArray->At(fiCut))->PionIsSelectedMC(labelPosPion,fMCStack) // test positive pion
1901 ) {
1902 if(particle->GetPdgCode() == 221) fHistoMCEtaPiPlPiMiPiZeroInAccPt[fiCut]->Fill(particle->Pt(), weighted ); // MC Eta pi+ pi- pi0 with gamma's and e+e- in acc
1903 if(particle->GetPdgCode() == 223) fHistoMCOmegaPiPlPiMiPiZeroInAccPt[fiCut]->Fill(particle->Pt(), weighted ); // MC Omega pi+ pi- pi0 with gamma's and e+e- in acc
1904 }
1905 } else if (fNeutralPionMode == 2){
1906 if(
1907 ((AliCaloPhotonCuts*)fClusterCutArray->At(fiCut))->ClusterIsSelectedMC(gamma1,fMCStack) && // test first daugther of pi0
1908 ((AliCaloPhotonCuts*)fClusterCutArray->At(fiCut))->ClusterIsSelectedMC(gamma2,fMCStack) && // test second daughter of pi0
1909 ((AliPrimaryPionCuts*)fPionCutArray->At(fiCut))->PionIsSelectedMC(labelNegPion,fMCStack) && // test negative pion
1910 ((AliPrimaryPionCuts*)fPionCutArray->At(fiCut))->PionIsSelectedMC(labelPosPion,fMCStack) // test positive pion
1911 ) {
1912 if(particle->GetPdgCode() == 221) fHistoMCEtaPiPlPiMiPiZeroInAccPt[fiCut]->Fill(particle->Pt(), weighted ); // MC Eta pi+ pi- pi0 with gamma's and e+e- in acc
1913 if(particle->GetPdgCode() == 223) fHistoMCOmegaPiPlPiMiPiZeroInAccPt[fiCut]->Fill(particle->Pt(), weighted ); // MC Omega pi+ pi- pi0 with gamma's and e+e- in acc
1914 }
1915 }
1916
1917 }
1918 }
1919 }
1920}
1921
1922
1923//________________________________________________________________________
1924void AliAnalysisTaskNeutralMesonToPiPlPiMiPiZero::CalculateMesonCandidates(){
1925
1926 // Conversion Gammas
1927 if( fNeutralPionCandidates->GetEntries() > 0 && fGoodVirtualParticles->GetEntries() > 0 ){
1928
1929 for(Int_t mesonIndex=0; mesonIndex<fNeutralPionCandidates->GetEntries(); mesonIndex++){
1930 AliAODConversionMother *neutralPion=dynamic_cast<AliAODConversionMother*>(fNeutralPionCandidates->At(mesonIndex));
1931 if (neutralPion==NULL) continue;
1932 for(Int_t virtualParticleIndex=0;virtualParticleIndex<fGoodVirtualParticles->GetEntries();virtualParticleIndex++){
1933
1934 AliAODConversionPhoton *vParticle=dynamic_cast<AliAODConversionPhoton*>(fGoodVirtualParticles->At(virtualParticleIndex));
1935 if (vParticle==NULL) continue;
1936 //Check for same Electron ID
1937
1938 AliAODConversionMother *mesoncand = new AliAODConversionMother(neutralPion,vParticle);
1939 mesoncand->SetLabels(mesonIndex,virtualParticleIndex);
1940
1941 if( ( ((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->MesonIsSelected(mesoncand,kTRUE,((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetEtaShift())) ){
1942
1943// cout<< "Meson Accepted "<<endl;
1944
1945 Int_t zbin= fBGHandler[fiCut]->GetZBinIndex(fESDEvent->GetPrimaryVertex()->GetZ());
1946 Int_t mbin = 0;
1947 if( ((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->UseTrackMultiplicity()){
1948 mbin = fBGHandler[fiCut]->GetMultiplicityBinIndex(fNumberOfESDTracks);
1949 } else {
1950 if (fNeutralPionMode < 2) mbin = fBGHandler[fiCut]->GetMultiplicityBinIndex(fGoodConvGammas->GetEntries());
1951 else mbin = fBGHandler[fiCut]->GetMultiplicityBinIndex(fClusterCandidates->GetEntries());
1952 }
1953
1954 fHistoMotherInvMassPt[fiCut]->Fill(mesoncand->M(),mesoncand->Pt());
1955 Double_t sparesFill[4] = {mesoncand->M(),mesoncand->Pt(),(Double_t)zbin,(Double_t)mbin};
1956 fTHnSparseMotherInvMassPtZM[fiCut]->Fill(sparesFill,1);
1957
1958 if(fMCEvent){
1959 ProcessTrueMesonCandidates(mesoncand,neutralPion,vParticle);
1960 }
1961 }
1962 delete mesoncand;
1963 mesoncand=0x0;
1964 }
1965 }
1966 }
1967}
1968
1969//________________________________________________________________________
1970void AliAnalysisTaskNeutralMesonToPiPlPiMiPiZero::CalculateBackground(){
1971
1972 Int_t zbin= fBGHandler[fiCut]->GetZBinIndex(fESDEvent->GetPrimaryVertex()->GetZ());
1973 Int_t mbin = 0;
1974
1975
1976 if(((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->UseTrackMultiplicity()){
1977 mbin = fBGHandler[fiCut]->GetMultiplicityBinIndex(fNumberOfESDTracks);
1978 } else {
1979 if (fNeutralPionMode < 2) mbin = fBGHandler[fiCut]->GetMultiplicityBinIndex(fGoodConvGammas->GetEntries());
1980 else mbin = fBGHandler[fiCut]->GetMultiplicityBinIndex(fClusterCandidates->GetEntries());
1981 }
1982
1983 Int_t method = 1;
1984 AliGammaConversionAODBGHandler::GammaConversionVertex *bgEventVertex = NULL;
1985 if( ((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->UseTrackMultiplicity() ) {
1986 for(Int_t nEventsInBG=0;nEventsInBG<fBGHandler[fiCut]->GetNBGEvents();nEventsInBG++){
1987 AliGammaConversionMotherAODVector *previousEventMesons = fBGHandler[fiCut]->GetBGGoodMesons(zbin,mbin,nEventsInBG);
1988 if(fMoveParticleAccordingToVertex == kTRUE && method == 1){
1989 bgEventVertex = fBGHandler[fiCut]->GetBGEventVertex(zbin,mbin,nEventsInBG);
1990 }
1991
1992 for(Int_t iCurrent=0;iCurrent<fGoodVirtualParticles->GetEntries();iCurrent++){
1993 AliAODConversionPhoton currentEventGoodV0 = *(AliAODConversionPhoton*)(fGoodVirtualParticles->At(iCurrent));
1994
1995 for(UInt_t iPrevious=0;iPrevious<previousEventMesons->size();iPrevious++){
1996 AliAODConversionMother previousGoodMeson = (AliAODConversionMother)(*(previousEventMesons->at(iPrevious)));
1997
1998 if(fMoveParticleAccordingToVertex == kTRUE && method == 1 ){
1999 MoveParticleAccordingToVertex(&previousGoodMeson,bgEventVertex);
2000 }
2001
2002 AliAODConversionMother *backgroundCandidate = new AliAODConversionMother(&previousGoodMeson,&currentEventGoodV0);
2003
2004
2005 if( ( ((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->MesonIsSelected(backgroundCandidate,kFALSE, ((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetEtaShift()))){
2006 fHistoMotherBackInvMassPt[fiCut]->Fill(backgroundCandidate->M(),backgroundCandidate->Pt());
2007 Double_t sparesFill[4] = {backgroundCandidate->M(),backgroundCandidate->Pt(),(Double_t)zbin,(Double_t)mbin};
2008 fTHnSparseMotherBackInvMassPtZM[fiCut]->Fill(sparesFill,1);
2009 }
2010 delete backgroundCandidate;
2011 backgroundCandidate = 0x0;
2012 }
2013 }
2014 }
2015 } else {
2016 for(Int_t nEventsInBG=0;nEventsInBG <fBGHandler[fiCut]->GetNBGEvents();nEventsInBG++){
2017 AliGammaConversionMotherAODVector *previousEventMesons = fBGHandler[fiCut]->GetBGGoodMesons(zbin,mbin,nEventsInBG);
2018 if(previousEventMesons){
2019 if(fMoveParticleAccordingToVertex == kTRUE && method == 1){
2020 bgEventVertex = fBGHandler[fiCut]->GetBGEventVertex(zbin,mbin,nEventsInBG);
2021 }
2022 for(Int_t iCurrent=0;iCurrent<fGoodVirtualParticles->GetEntries();iCurrent++){
2023
2024 AliAODConversionPhoton currentEventGoodV0 = *(AliAODConversionPhoton*)(fGoodVirtualParticles->At(iCurrent));
2025
2026 for(UInt_t iPrevious=0;iPrevious<previousEventMesons->size();iPrevious++){
2027
2028 AliAODConversionMother previousGoodMeson = (AliAODConversionMother)(*(previousEventMesons->at(iPrevious)));
2029
2030 if(fMoveParticleAccordingToVertex == kTRUE && method ==1){
2031 MoveParticleAccordingToVertex(&previousGoodMeson,bgEventVertex);
2032 }
2033
2034 AliAODConversionMother *backgroundCandidate = new AliAODConversionMother(&previousGoodMeson,&currentEventGoodV0);
2035
2036 if((((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->MesonIsSelected(backgroundCandidate,kFALSE,((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetEtaShift()))){
2037 fHistoMotherBackInvMassPt[fiCut]->Fill(backgroundCandidate->M(),backgroundCandidate->Pt());
2038 Double_t sparesFill[4] = {backgroundCandidate->M(),backgroundCandidate->Pt(),(Double_t)zbin,(Double_t)mbin};
2039 fTHnSparseMotherBackInvMassPtZM[fiCut]->Fill(sparesFill,1);
2040 }
2041 delete backgroundCandidate;
2042 backgroundCandidate = 0x0;
2043 }
2044 }
2045 }
2046 }
2047 }
2048}
2049
2050//______________________________________________________________________
2051void AliAnalysisTaskNeutralMesonToPiPlPiMiPiZero::ProcessTrueMesonCandidates(AliAODConversionMother *mesoncand, AliAODConversionMother *TrueNeutralPionCandidate, AliAODConversionPhoton *TrueVirtualParticleCandidate){
2052
2053 // Process True Mesons
2054 AliStack *MCStack = fMCEvent->Stack();
2055
2056 Bool_t isTrueEta = kFALSE;
2057 Bool_t isTrueOmega = kFALSE;
2058 Int_t trueMesonFlag = TrueNeutralPionCandidate->GetTrueMesonValue();
2059 Int_t pi0MCLabel= TrueNeutralPionCandidate->GetMCLabel();
2060
2061
2062 if ( !(trueMesonFlag == 1 && pi0MCLabel != -1)) return;
2063// cout << trueMesonFlag << "\t" << pi0MCLabel << endl;
2064
2065
2066 Int_t virtualParticleMCLabel = TrueVirtualParticleCandidate->GetMCParticleLabel(MCStack);
2067 Int_t virtualParticleMotherLabel = -1;
2068 Bool_t isPiPiDecay = kFALSE;
2069
2070// if (fDoMesonQA){
2071 TParticle * negativeMC = (TParticle*)TrueVirtualParticleCandidate->GetNegativeMCDaughter(MCStack);
2072 TParticle * positiveMC = (TParticle*)TrueVirtualParticleCandidate->GetPositiveMCDaughter(MCStack);
2073// }
2074
2075 if(virtualParticleMCLabel != -1){ // if virtualParticleMCLabel==-1 particles don't have same mother
2076// TParticle * negativeMC = (TParticle*)TrueVirtualParticleCandidate->GetNegativeMCDaughter(MCStack);
2077// TParticle * positiveMC = (TParticle*)TrueVirtualParticleCandidate->GetPositiveMCDaughter(MCStack);
2078// TParticle * virtualParticleMotherMC = (TParticle*)MCStack->Particle(virtualParticleMCLabel);
2079// cout << "pdg code same mother - " << virtualParticleMotherMC->GetPdgCode() << endl;
2080
2081 if(TMath::Abs(negativeMC->GetPdgCode())==211 && TMath::Abs(positiveMC->GetPdgCode())==211){ // Pions ...
2082 virtualParticleMotherLabel=virtualParticleMCLabel;
2083 isPiPiDecay=kTRUE;
2084// } else if(TMath::Abs(negativeMC->GetPdgCode())==11 && TMath::Abs(positiveMC->GetPdgCode())==11){ // Electrons ...
2085// if( virtualParticleMotherMC->GetPdgCode() != 22 ){
2086// virtualParticleMotherLabel=virtualParticleMCLabel;
2087// isDalitz = kTRUE;
2088// } else if(negativeMC->GetUniqueID() == 5 && positiveMC->GetUniqueID() ==5){ // ... From Conversion ...
2089// virtualParticleMotherLabel=virtualParticleMotherMC->GetFirstMother();
2090// isRealGamma = kTRUE; //no virtual gamma
2091// }
2092 }
2093 }
2094 if (IsEtaPiPlPiMiPiZeroDaughter(pi0MCLabel) || IsOmegaPiPlPiMiPiZeroDaughter(pi0MCLabel)){
2095 Int_t pi0MotherMCLabel = ((TParticle*)MCStack->Particle(pi0MCLabel))->GetMother(0);
2096 if(virtualParticleMCLabel != -1){
2097// cout << "pi+pi- mother: "<< virtualParticleMCLabel << endl;
2098// cout << "pi0 mother: "<< pi0MotherMCLabel << endl;
2099
2100// TParticle * virtualParticleMotherMC = (TParticle*)MCStack->Particle(virtualParticleMCLabel);
2101// cout << "pdg code same mother - " << virtualParticleMotherMC->GetPdgCode() << endl;
2102 }
2103 if( pi0MotherMCLabel == virtualParticleMotherLabel ){
2104 if(((TParticle*)MCStack->Particle(virtualParticleMotherLabel))->GetPdgCode() == 221){
2105// cout << "found eta" << endl;
2106 isTrueEta=kTRUE;
2107 }
2108 if(((TParticle*)MCStack->Particle(virtualParticleMotherLabel))->GetPdgCode() == 223){
2109// cout << "found omega" << endl;
2110 isTrueOmega=kTRUE;
2111 }
2112 }
2113 }
2114
2115
2116 if( isTrueEta || isTrueOmega ){ // True Eta or Omega
2117 if ( isPiPiDecay) { //real eta -> Pi+ Pi- Pi0
2118 Float_t weighted= 1;
2119// if( ((AliPrimaryPionCuts*) fPionCutArray->At(fiCut))->DoWeights() ) {
2120// if(((AliConvEventCuts*)fEventCutArray->At(fiCut))->IsParticleFromBGEvent(gammaMotherLabel, fMCStack,fInputEvent)){
2121// if (((TParticle*)MCStack->Particle(gammaMotherLabel))->Pt()>0.005){
2122// weighted= ((AliConvEventCuts*)fEventCutArray->At(fiCut))->GetWeightForMeson(fV0Reader->GetPeriodName(),gammaMotherLabel,fMCStack,fInputEvent);
2123// }
2124// }
2125// }
2126 fHistoTrueMotherPiPlPiMiPiZeroInvMassPt[fiCut]->Fill(mesoncand->M(),mesoncand->Pt(),weighted);
2127 }
2128 }
2129
2130}
2131
2132
2133//________________________________________________________________________
2134void AliAnalysisTaskNeutralMesonToPiPlPiMiPiZero::UpdateEventByEventData(){
2135 //see header file for documentation
2136
2137 Int_t method = 1;
2138 if( method == 1 ) {
2139 if(fNeutralPionCandidates->GetEntries() >0 ){
2140 if(((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->UseTrackMultiplicity()){
2141 fBGHandler[fiCut]->AddMesonEvent(fNeutralPionCandidates,fInputEvent->GetPrimaryVertex()->GetX(),fInputEvent->GetPrimaryVertex()->GetY(),fInputEvent->GetPrimaryVertex()->GetZ(),fV0Reader->GetNumberOfPrimaryTracks(),0);
2142 } else { // means we use #V0s for multiplicity
2143 if (fNeutralPionMode < 2) fBGHandler[fiCut]->AddMesonEvent(fNeutralPionCandidates,fInputEvent->GetPrimaryVertex()->GetX(),fInputEvent->GetPrimaryVertex()->GetY(),fInputEvent->GetPrimaryVertex()->GetZ(),fGoodConvGammas->GetEntries(),0);
2144 else fBGHandler[fiCut]->AddMesonEvent(fNeutralPionCandidates,fInputEvent->GetPrimaryVertex()->GetX(),fInputEvent->GetPrimaryVertex()->GetY(),fInputEvent->GetPrimaryVertex()->GetZ(),fClusterCandidates->GetEntries(),0);
2145 }
2146 }
2147 } else if ( method == 2 ){
2148 if(fGoodVirtualParticles->GetEntries() > 0 ){
2149 if(((AliConversionMesonCuts*)fMesonCutArray->At(fiCut))->UseTrackMultiplicity()){
2150 fBGHandler[fiCut]->AddEvent(fGoodVirtualParticles,fInputEvent->GetPrimaryVertex()->GetX(),fInputEvent->GetPrimaryVertex()->GetY(),fInputEvent->GetPrimaryVertex()->GetZ(),fV0Reader->GetNumberOfPrimaryTracks(),0);
2151 } else{ // means we use #V0s for multiplicity
2152 fBGHandler[fiCut]->AddEvent(fGoodVirtualParticles,fInputEvent->GetPrimaryVertex()->GetX(),fInputEvent->GetPrimaryVertex()->GetY(),fInputEvent->GetPrimaryVertex()->GetZ(),fGoodVirtualParticles->GetEntries(),0);
2153 }
2154 }
2155 }
2156}
2157
2158//________________________________________________________________________
2159void AliAnalysisTaskNeutralMesonToPiPlPiMiPiZero::MoveParticleAccordingToVertex(AliAODConversionMother* particle,const AliGammaConversionAODBGHandler::GammaConversionVertex *vertex){
2160 //see header file for documentation
2161
2162 Double_t dx = vertex->fX - fInputEvent->GetPrimaryVertex()->GetX();
2163 Double_t dy = vertex->fY - fInputEvent->GetPrimaryVertex()->GetY();
2164 Double_t dz = vertex->fZ - fInputEvent->GetPrimaryVertex()->GetZ();
2165
2166 Double_t movedPlace[3] = {particle->GetProductionX() - dx,particle->GetProductionY() - dy,particle->GetProductionZ() - dz};
2167 particle->SetProductionPoint(movedPlace);
2168}
2169
2170//_____________________________________________________________________________________
2171Bool_t AliAnalysisTaskNeutralMesonToPiPlPiMiPiZero::IsEtaPiPlPiMiPiZeroDaughter( Int_t label ) const {
2172//
2173// Returns true if the particle comes from eta -> pi+ pi- gamma
2174//
2175 Int_t motherLabel = fMCStack->Particle( label )->GetMother(0);
2176 if( motherLabel < 0 || motherLabel >= fMCStack->GetNtrack() ) return kFALSE;
2177
2178 TParticle* mother = fMCStack->Particle( motherLabel );
2179// cout << "found eta? " << endl;
2180 if( mother->GetPdgCode() != 221 ) return kFALSE;
2181// else cout << "YES" << endl;
2182 if( IsPiPlPiMiPiZeroDecay( mother ) ) return kTRUE;
2183 return kFALSE;
2184}
2185
2186//_____________________________________________________________________________________
2187Bool_t AliAnalysisTaskNeutralMesonToPiPlPiMiPiZero::IsOmegaPiPlPiMiPiZeroDaughter( Int_t label ) const {
2188//
2189// Returns true if the particle comes from eta -> pi+ pi- gamma
2190//
2191 Int_t motherLabel = fMCStack->Particle( label )->GetMother(0);
2192 if( motherLabel < 0 || motherLabel >= fMCStack->GetNtrack() ) return kFALSE;
2193
2194 TParticle* mother = fMCStack->Particle( motherLabel );
2195// cout << "found omega? " << endl;
2196 if( mother->GetPdgCode() != 223 ) return kFALSE;
2197// else cout << "YES" << endl;
2198 if( IsPiPlPiMiPiZeroDecay( mother ) ) return kTRUE;
2199 return kFALSE;
2200}
2201
2202
2203//_____________________________________________________________________________
2204Bool_t AliAnalysisTaskNeutralMesonToPiPlPiMiPiZero::IsPiPlPiMiPiZeroDecay(TParticle *fMCMother) const
2205{
2206// cout << fMCMother->GetNDaughters() << endl;
2207 if( fMCMother->GetNDaughters() != 3 ) return kFALSE;
2208// cout << fMCMother->GetPdgCode() << endl;
2209 if( !(fMCMother->GetPdgCode() == 221 || fMCMother->GetPdgCode() == 223) ) return kFALSE;
2210// cout << "made it til here" << endl;
2211
2212 TParticle *posPion = 0x0;
2213 TParticle *negPion = 0x0;
2214 TParticle *neutPion = 0x0;
2215
2216 for(Int_t index= fMCMother->GetFirstDaughter();index<= fMCMother->GetLastDaughter();index++){
2217 TParticle* temp = (TParticle*)fMCStack->Particle( index );
2218
2219 switch( temp->GetPdgCode() ) {
2220 case 211:
2221 posPion = temp;
2222 break;
2223 case -211:
2224 negPion = temp;
2225 break;
2226 case 111:
2227 neutPion = temp;
2228 break;
2229 }
2230 }
2231 if( posPion && negPion && neutPion) return kTRUE;
2232
2233 return kFALSE;
2234}
2235
2236//_____________________________________________________________________________________
2237Bool_t AliAnalysisTaskNeutralMesonToPiPlPiMiPiZero::GammaIsNeutralMesonPiPlPiMiPiZeroDaughter( Int_t label ) const {
2238//
2239// Returns true if the particle comes from eta -> pi+ pi- gamma
2240//
2241 Int_t motherLabel = fMCStack->Particle( label )->GetMother(0);
2242 if( motherLabel < 0 || motherLabel >= fMCStack->GetNtrack() ) return kFALSE;
2243
2244 TParticle* mother = fMCStack->Particle( motherLabel );
2245// cout << "found omega? " << endl;
2246 if( mother->GetPdgCode() != 111 ) return kFALSE;
2247// else cout << "YES" << endl;
2248 Int_t grandMotherLabel = mother->GetMother(0);
2249 if( grandMotherLabel < 0 || grandMotherLabel >= fMCStack->GetNtrack() ) return kFALSE;
2250 TParticle* grandmother = fMCStack->Particle( grandMotherLabel );
2251
2252 if( IsPiPlPiMiPiZeroDecay( grandmother ) ) return kTRUE;
2253 return kFALSE;
2254}