]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGLF/SPECTRA/ChargedHadrons/dNdPt/AlidNdPtAnalysisPbPbAOD.cxx
improved parameter for deuteron PID
[u/mrichter/AliRoot.git] / PWGLF / SPECTRA / ChargedHadrons / dNdPt / AlidNdPtAnalysisPbPbAOD.cxx
CommitLineData
d25bcbe6 1/**************************************************************************
2* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3* *
4* Author: The ALICE Off-line Project. *
5* Contributors are mentioned in the code where appropriate. *
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// AlidNdPtAnalysisPbPbAOD class.
17//
18// Author: P. Luettig, 15.05.2013
19//------------------------------------------------------------------------------
20
21
22#include "AlidNdPtAnalysisPbPbAOD.h"
23
24
25using namespace std;
26
27ClassImp(AlidNdPtAnalysisPbPbAOD)
28
29// dummy constructor
30AlidNdPtAnalysisPbPbAOD::AlidNdPtAnalysisPbPbAOD() : AliAnalysisTaskSE(),
31fOutputList(0),
32// Histograms
33hPt(0),
34hMCPt(0),
35hnZvPtEtaCent(0),
36hnMCRecPrimZvPtEtaCent(0),
37hnMCGenZvPtEtaCent(0),
5747f2c6 38hnMCRecSecZvPtEtaCent(0),
d25bcbe6 39hEventStatistics(0),
40hEventStatisticsCentrality(0),
41hAllEventStatisticsCentrality(0),
fad9b70b 42hEventStatisticsCentralityTrigger(0),
d25bcbe6 43hnZvMultCent(0),
44hTriggerStatistics(0),
45hMCTrackPdgCode(0),
46hMCTrackStatusCode(0),
47hCharge(0),
48hMCCharge(0),
49hMCPdgPt(0),
50hMCHijingPrim(0),
51hAccNclsTPC(0),
95f26ffa 52hAllCrossedRowsTPC(0),
53hFilterCrossedRowsTPC(0),
d25bcbe6 54hAccCrossedRowsTPC(0),
5747f2c6 55hDCAPtAll(0),
56hDCAPtAccepted(0),
57hMCDCAPtSecondary(0),
58hMCDCAPtPrimary(0),
d25bcbe6 59//global
60bIsMonteCarlo(0),
61// event cut variables
62dCutMaxZVertex(10.),
63// track kinematic cut variables
64dCutPtMin(0.15),
65dCutPtMax(1000.),
66dCutEtaMin(-0.8),
67dCutEtaMax(0.8),
68// track quality cut variables
69bCutRequireTPCRefit(kTRUE),
70dCutMinNumberOfCrossedRows(120.),
71dCutMinRatioCrossedRowsOverFindableClustersTPC(0.8),
72dCutMaxChi2PerClusterTPC(4.),
73dCutMaxFractionSharedTPCClusters(0.4),
74dCutMaxDCAToVertexZ(3.0),
75dCutMaxDCAToVertexXY(3.0),
76bCutRequireITSRefit(kTRUE),
77dCutMaxChi2PerClusterITS(36.),
78dCutDCAToVertex2D(kFALSE),
79dCutRequireSigmaToVertex(kFALSE),
80dCutMaxDCAToVertexXYPtDepPar0(0.0182),
81dCutMaxDCAToVertexXYPtDepPar1(0.0350),
82dCutMaxDCAToVertexXYPtDepPar2(1.01),
83bCutAcceptKinkDaughters(kFALSE),
84dCutMaxChi2TPCConstrainedGlobal(36.),
85// binning for THnSparse
86fMultNbins(0),
87fPtNbins(0),
88fPtCorrNbins(0),
89fEtaNbins(0),
90fZvNbins(0),
91fCentralityNbins(0),
92fBinsMult(0),
93fBinsPt(0),
94fBinsPtCorr(0),
95fBinsEta(0),
96fBinsZv(0),
97fBinsCentrality(0)
98{
99
100 fMultNbins = 0;
101 fPtNbins = 0;
102 fPtCorrNbins = 0;
103 fEtaNbins = 0;
104 fZvNbins = 0;
105 fCentralityNbins = 0;
106 fBinsMult = 0;
107 fBinsPt = 0;
108 fBinsPtCorr = 0;
109 fBinsEta = 0;
110 fBinsEta = 0;
111 fBinsZv = 0;
112 fBinsCentrality = 0;
113
114}
115
116AlidNdPtAnalysisPbPbAOD::AlidNdPtAnalysisPbPbAOD(const char *name) : AliAnalysisTaskSE(name),
117fOutputList(0),
118// Histograms
119hPt(0),
120hMCPt(0),
121hnZvPtEtaCent(0),
122hnMCRecPrimZvPtEtaCent(0),
123hnMCGenZvPtEtaCent(0),
5747f2c6 124hnMCRecSecZvPtEtaCent(0),
d25bcbe6 125hEventStatistics(0),
126hEventStatisticsCentrality(0),
127hAllEventStatisticsCentrality(0),
fad9b70b 128hEventStatisticsCentralityTrigger(0),
d25bcbe6 129hnZvMultCent(0),
130hTriggerStatistics(0),
131hMCTrackPdgCode(0),
132hMCTrackStatusCode(0),
133hCharge(0),
134hMCCharge(0),
135hMCPdgPt(0),
136hMCHijingPrim(0),
137hAccNclsTPC(0),
95f26ffa 138hAllCrossedRowsTPC(0),
139hFilterCrossedRowsTPC(0),
d25bcbe6 140hAccCrossedRowsTPC(0),
5747f2c6 141hDCAPtAll(0),
142hDCAPtAccepted(0),
143hMCDCAPtSecondary(0),
144hMCDCAPtPrimary(0),
d25bcbe6 145//global
146bIsMonteCarlo(0),
147// event cut variables
148dCutMaxZVertex(10.),
149// track kinematic cut variables
150dCutPtMin(0.15),
151dCutPtMax(200.),
152dCutEtaMin(-0.8),
153dCutEtaMax(0.8),
154// track quality cut variables
155bCutRequireTPCRefit(kTRUE),
156dCutMinNumberOfCrossedRows(120.),
157dCutMinRatioCrossedRowsOverFindableClustersTPC(0.8),
158dCutMaxChi2PerClusterTPC(4.),
159dCutMaxFractionSharedTPCClusters(0.4),
160dCutMaxDCAToVertexZ(3.0),
161dCutMaxDCAToVertexXY(3.0),
162bCutRequireITSRefit(kTRUE),
163dCutMaxChi2PerClusterITS(36.),
164dCutDCAToVertex2D(kFALSE),
165dCutRequireSigmaToVertex(kFALSE),
166dCutMaxDCAToVertexXYPtDepPar0(0.0182),
167dCutMaxDCAToVertexXYPtDepPar1(0.0350),
168dCutMaxDCAToVertexXYPtDepPar2(1.01),
169bCutAcceptKinkDaughters(kFALSE),
170dCutMaxChi2TPCConstrainedGlobal(36.),
171// binning for THnSparse
172fMultNbins(0),
173fPtNbins(0),
174fPtCorrNbins(0),
175fEtaNbins(0),
176fZvNbins(0),
177fCentralityNbins(0),
178fBinsMult(0),
179fBinsPt(0),
180fBinsPtCorr(0),
181fBinsEta(0),
182fBinsZv(0),
183fBinsCentrality(0)
184{
185 fMultNbins = 0;
186 fPtNbins = 0;
187 fPtCorrNbins = 0;
188 fEtaNbins = 0;
189 fZvNbins = 0;
190 fCentralityNbins = 0;
191 fBinsMult = 0;
192 fBinsPt = 0;
193 fBinsPtCorr = 0;
194 fBinsEta = 0;
195 fBinsEta = 0;
196 fBinsZv = 0;
197 fBinsCentrality = 0;
198
199 DefineOutput(1, TList::Class());
200}
201
202// destructor
203AlidNdPtAnalysisPbPbAOD::~AlidNdPtAnalysisPbPbAOD()
204{
205 if(hnZvPtEtaCent) delete hnZvPtEtaCent; hnZvPtEtaCent = 0;
206 if(hPt) delete hPt; hPt = 0;
5747f2c6 207 if(hnMCRecPrimZvPtEtaCent) delete hnMCRecPrimZvPtEtaCent; hnMCRecPrimZvPtEtaCent = 0;
208 if(hnMCGenZvPtEtaCent) delete hnMCGenZvPtEtaCent; hnMCGenZvPtEtaCent = 0;
209 if(hnMCRecSecZvPtEtaCent) delete hnMCRecSecZvPtEtaCent; hnMCRecSecZvPtEtaCent = 0;
210 if(hMCPt) delete hMCPt; hMCPt = 0;
211 if(hEventStatistics) delete hEventStatistics; hEventStatistics = 0;
212 if(hEventStatisticsCentrality) delete hEventStatisticsCentrality; hEventStatisticsCentrality = 0;
213 if(hAllEventStatisticsCentrality) delete hAllEventStatisticsCentrality; hAllEventStatisticsCentrality = 0;
214 if(hnZvMultCent) delete hnZvMultCent; hnZvMultCent = 0;
215 if(hTriggerStatistics) delete hTriggerStatistics; hTriggerStatistics = 0;
216 if(hMCTrackPdgCode) delete hMCTrackPdgCode; hMCTrackPdgCode = 0;
217 if(hMCTrackStatusCode) delete hMCTrackStatusCode; hMCTrackStatusCode = 0;
218 if(hCharge) delete hCharge; hCharge = 0;
219 if(hMCCharge) delete hMCCharge; hMCCharge = 0;
220 if(hMCPdgPt) delete hMCPdgPt; hMCPdgPt = 0;
221 if(hMCHijingPrim) delete hMCHijingPrim; hMCHijingPrim = 0;
222 if(hAccNclsTPC) delete hAccNclsTPC; hAccNclsTPC = 0;
223 if(hAccCrossedRowsTPC) delete hAccCrossedRowsTPC; hAccCrossedRowsTPC = 0;
224 if(hDCAPtAll) delete hDCAPtAll; hDCAPtAll = 0;
225 if(hDCAPtAccepted) delete hDCAPtAccepted; hDCAPtAccepted = 0;
226 if(hMCDCAPtSecondary) delete hMCDCAPtSecondary; hMCDCAPtSecondary = 0;
227 if(hMCDCAPtPrimary) delete hMCDCAPtPrimary; hMCDCAPtPrimary = 0;
d25bcbe6 228}
229
230void AlidNdPtAnalysisPbPbAOD::UserCreateOutputObjects()
231{
232 // create all output histograms here
233 OpenFile(1, "RECREATE");
234
235 fOutputList = new TList();
236 fOutputList->SetOwner();
237
238 //define default binning
239 Double_t binsMultDefault[48] = {-0.5, 0.5 , 1.5 , 2.5 , 3.5 , 4.5 , 5.5 , 6.5 , 7.5 , 8.5,9.5, 10.5, 11.5, 12.5, 13.5, 14.5, 15.5, 16.5, 17.5, 18.5,19.5, 20.5, 30.5, 40.5 , 50.5 , 60.5 , 70.5 , 80.5 , 90.5 , 100.5,200.5, 300.5, 400.5, 500.5, 600.5, 700.5, 800.5, 900.5, 1000.5, 2000.5, 3000.5, 4000.5, 5000.5, 6000.5, 7000.5, 8000.5, 9000.5, 10000.5 };
240 Double_t binsPtDefault[82] = {0.0, 0.05, 0.1, 0.15, 0.2, 0.25, 0.3, 0.35, 0.4, 0.45, 0.5, 0.55, 0.6, 0.65, 0.7, 0.75, 0.8, 0.85, 0.9, 0.95, 1.0, 1.1, 1.2, 1.3, 1.4, 1.5, 1.6, 1.7, 1.8, 1.9, 2.0, 2.2, 2.4, 2.6, 2.8, 3.0, 3.2, 3.4, 3.6, 3.8, 4.0, 4.5, 5.0, 5.5, 6.0, 6.5, 7.0, 8.0, 9.0, 10.0, 11.0, 12.0, 13.0, 14.0, 15.0, 16.0, 18.0, 20.0, 22.0, 24.0, 26.0, 28.0, 30.0, 32.0, 34.0, 36.0, 40.0, 45.0, 50.0, 60.0, 70.0, 80.0, 90.0, 100.0, 110.0, 120.0, 130.0, 140.0, 150.0, 160.0, 180.0, 200.0};
241 Double_t binsPtCorrDefault[37] = {0.,0.05,0.1,0.15,0.2,0.25,0.3,0.35,0.4,0.45,0.5,0.55,0.6,0.65,0.7,0.75,0.8,0.85,0.9,0.95,1.0,1.1,1.2,1.3,1.4,1.5,1.6,1.7,1.8,1.9,2.0,2.2,2.4,2.6,3.0,4.0,200.0};
242 Double_t binsEtaDefault[31] = {-1.5,-1.4,-1.3,-1.2,-1.1,-1.0,-0.9,-0.8,-0.7,-0.6,-0.5,-0.4,-0.3,-0.2,-0.1,0.,0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,0.9,1.0,1.1,1.2,1.3,1.4,1.5};
243 Double_t binsZvDefault[13] = {-30.,-25.,-20.,-15.,-10.,-5.,0.,5.,10.,15.,20.,25.,30.};
244 Double_t binsCentralityDefault[12] = {0., 5., 10., 20., 30., 40., 50., 60., 70., 80., 90., 100.};
245
246 // if no binning is set, use the default
247 if (!fBinsMult) { SetBinsMult(48,binsMultDefault); }
248 if (!fBinsPt) { SetBinsPt(82,binsPtDefault); }
249 if (!fBinsPtCorr) { SetBinsPtCorr(37,binsPtCorrDefault); }
250 if (!fBinsEta) { SetBinsEta(31,binsEtaDefault); }
251 if (!fBinsZv) { SetBinsZv(13,binsZvDefault); }
252 if (!fBinsCentrality) { SetBinsCentrality(12,binsCentralityDefault); }
253
254 Int_t binsZvPtEtaCent[4]={fZvNbins-1,fPtNbins-1,fEtaNbins-1,fCentralityNbins-1};
255 Int_t binsZvMultCent[3]={fZvNbins-1,fMultNbins-1,fCentralityNbins-1};
256
257 // define Histograms
258 hnZvPtEtaCent = new THnSparseF("hnZvPtEtaCent","Zv:Pt:Eta:Centrality",4,binsZvPtEtaCent);
259 hnZvPtEtaCent->SetBinEdges(0,fBinsZv);
260 hnZvPtEtaCent->SetBinEdges(1,fBinsPt);
261 hnZvPtEtaCent->SetBinEdges(2,fBinsEta);
262 hnZvPtEtaCent->SetBinEdges(3,fBinsCentrality);
263 hnZvPtEtaCent->GetAxis(0)->SetTitle("Zv (cm)");
264 hnZvPtEtaCent->GetAxis(1)->SetTitle("Pt (GeV/c)");
265 hnZvPtEtaCent->GetAxis(2)->SetTitle("Eta");
266 hnZvPtEtaCent->GetAxis(3)->SetTitle("Centrality");
267 hnZvPtEtaCent->Sumw2();
268
269 hnMCRecPrimZvPtEtaCent = new THnSparseF("hnMCRecPrimZvPtEtaCent","mcZv:mcPt:mcEta:Centrality",4,binsZvPtEtaCent);
270 hnMCRecPrimZvPtEtaCent->SetBinEdges(0,fBinsZv);
271 hnMCRecPrimZvPtEtaCent->SetBinEdges(1,fBinsPt);
272 hnMCRecPrimZvPtEtaCent->SetBinEdges(2,fBinsEta);
273 hnMCRecPrimZvPtEtaCent->SetBinEdges(3,fBinsCentrality);
274 hnMCRecPrimZvPtEtaCent->GetAxis(0)->SetTitle("MC Zv (cm)");
275 hnMCRecPrimZvPtEtaCent->GetAxis(1)->SetTitle("MC Pt (GeV/c)");
276 hnMCRecPrimZvPtEtaCent->GetAxis(2)->SetTitle("MC Eta");
277 hnMCRecPrimZvPtEtaCent->GetAxis(3)->SetTitle("Centrality");
278 hnMCRecPrimZvPtEtaCent->Sumw2();
279
280 hnMCGenZvPtEtaCent = new THnSparseF("hnMCGenZvPtEtaCent","mcZv:mcPt:mcEta:Centrality",4,binsZvPtEtaCent);
281 hnMCGenZvPtEtaCent->SetBinEdges(0,fBinsZv);
282 hnMCGenZvPtEtaCent->SetBinEdges(1,fBinsPt);
283 hnMCGenZvPtEtaCent->SetBinEdges(2,fBinsEta);
284 hnMCGenZvPtEtaCent->SetBinEdges(3,fBinsCentrality);
285 hnMCGenZvPtEtaCent->GetAxis(0)->SetTitle("MC Zv (cm)");
286 hnMCGenZvPtEtaCent->GetAxis(1)->SetTitle("MC Pt (GeV/c)");
287 hnMCGenZvPtEtaCent->GetAxis(2)->SetTitle("MC Eta");
288 hnMCGenZvPtEtaCent->GetAxis(3)->SetTitle("Centrality");
289 hnMCGenZvPtEtaCent->Sumw2();
290
5747f2c6 291 hnMCRecSecZvPtEtaCent = new THnSparseF("hnMCRecSecZvPtEtaCent","mcZv:mcPt:mcEta:Centrality",4,binsZvPtEtaCent);
292 hnMCRecSecZvPtEtaCent->SetBinEdges(0,fBinsZv);
293 hnMCRecSecZvPtEtaCent->SetBinEdges(1,fBinsPt);
294 hnMCRecSecZvPtEtaCent->SetBinEdges(2,fBinsEta);
295 hnMCRecSecZvPtEtaCent->SetBinEdges(3,fBinsCentrality);
296 hnMCRecSecZvPtEtaCent->GetAxis(0)->SetTitle("MC Sec Zv (cm)");
297 hnMCRecSecZvPtEtaCent->GetAxis(1)->SetTitle("MC Sec Pt (GeV/c)");
298 hnMCRecSecZvPtEtaCent->GetAxis(2)->SetTitle("MC Sec Eta");
299 hnMCRecSecZvPtEtaCent->GetAxis(3)->SetTitle("Centrality");
300 hnMCRecSecZvPtEtaCent->Sumw2();
d25bcbe6 301
302 hPt = new TH1F("hPt","hPt",2000,0,200);
303 hPt->GetXaxis()->SetTitle("p_{T} (GeV/c)");
304 hPt->GetYaxis()->SetTitle("dN/dp_{T}");
305 hPt->Sumw2();
306
307 hMCPt = new TH1F("hMCPt","hMCPt",2000,0,200);
308 hMCPt->GetXaxis()->SetTitle("MC p_{T} (GeV/c)");
309 hMCPt->GetYaxis()->SetTitle("dN/dp_{T}");
310 hMCPt->Sumw2();
311
312 hEventStatistics = new TH1F("hEventStatistics","hEventStatistics",10,0,10);
313 hEventStatistics->GetYaxis()->SetTitle("number of events");
314 hEventStatistics->SetBit(TH1::kCanRebin);
315
316 hEventStatisticsCentrality = new TH1F("hEventStatisticsCentrality","hEventStatisticsCentrality",fCentralityNbins-1, fBinsCentrality);
317 hEventStatisticsCentrality->GetYaxis()->SetTitle("number of events");
318
319 hAllEventStatisticsCentrality = new TH1F("hAllEventStatisticsCentrality","hAllEventStatisticsCentrality",fCentralityNbins-1, fBinsCentrality);
320 hAllEventStatisticsCentrality->GetYaxis()->SetTitle("number of events");
321
fad9b70b 322 hEventStatisticsCentralityTrigger = new TH2F("hEventStatisticsCentralityTrigger","hEventStatisticsCentralityTrigger;centrality;trigger",100,0,100,3,0,3);
323 hEventStatisticsCentralityTrigger->Sumw2();
324
d25bcbe6 325 hnZvMultCent = new THnSparseF("hnZvMultCent","Zv:mult:Centrality",3,binsZvMultCent);
326 hnZvMultCent->SetBinEdges(0,fBinsZv);
327 hnZvMultCent->SetBinEdges(1,fBinsMult);
328 hnZvMultCent->SetBinEdges(2,fBinsCentrality);
329 hnZvMultCent->GetAxis(0)->SetTitle("Zv (cm)");
330 hnZvMultCent->GetAxis(1)->SetTitle("N_{acc}");
331 hnZvMultCent->GetAxis(2)->SetTitle("Centrality");
332 hnZvMultCent->Sumw2();
333
334 hTriggerStatistics = new TH1F("hTriggerStatistics","hTriggerStatistics",10,0,10);
335 hTriggerStatistics->GetYaxis()->SetTitle("number of events");
336
337 hMCTrackPdgCode = new TH1F("hMCTrackPdgCode","hMCTrackPdgCode",100,0,10);
338 hMCTrackPdgCode->GetYaxis()->SetTitle("number of tracks");
339 hMCTrackPdgCode->SetBit(TH1::kCanRebin);
340
341 hMCTrackStatusCode = new TH1F("hMCTrackStatusCode","hMCTrackStatusCode",100,0,10);
342 hMCTrackStatusCode->GetYaxis()->SetTitle("number of tracks");
343 hMCTrackStatusCode->SetBit(TH1::kCanRebin);
344
345 hCharge = new TH1F("hCharge","hCharge",30, -5, 5);
346 hCharge->GetXaxis()->SetTitle("Charge");
347 hCharge->GetYaxis()->SetTitle("number of tracks");
348
349 hMCCharge = new TH1F("hMCCharge","hMCCharge",30, -5, 5);
350 hMCCharge->GetXaxis()->SetTitle("MC Charge");
351 hMCCharge->GetYaxis()->SetTitle("number of tracks");
352
353 hMCPdgPt = new TH2F("hMCPdgPt","hMCPdgPt",fPtNbins-1, fBinsPt, 100,0,100);
354 hMCPdgPt->GetYaxis()->SetTitle("particle");
355 hMCPdgPt->GetXaxis()->SetTitle("Pt (GeV/c)");
356
357 hMCHijingPrim = new TH1F("hMCHijingPrim","hMCHijingPrim",2,0,2);
358 hMCPdgPt->GetYaxis()->SetTitle("number of particles");
359
360 hAccNclsTPC = new TH1F("hAccNclsTPC","hAccNclsTPC",160,0,159);
361 hAccNclsTPC->GetXaxis()->SetTitle("number of clusters per track after cut");
362
95f26ffa 363 hAllCrossedRowsTPC = new TH1F("hAllCrossedRowsTPC","hAllCrossedRowsTPC",160,0,159);
364 hAllCrossedRowsTPC->GetXaxis()->SetTitle("number of crossed rows per track for all tracks");
365
366 hFilterCrossedRowsTPC = new TH1F("hFilterCrossedRowsTPC","hFilterCrossedRowsTPC",160,0,159);
367 hFilterCrossedRowsTPC->GetXaxis()->SetTitle("number of crossed rows per track after filter bit");
368
d25bcbe6 369 hAccCrossedRowsTPC = new TH1F("hAccCrossedRowsTPC","hAccCrossedRowsTPC",160,0,159);
370 hAccCrossedRowsTPC->GetXaxis()->SetTitle("number of crossed rows per track after cut");
371
d0483ba3 372 Int_t binsDCAxyDCAzPt[3] = { 200,200, fPtNbins-1};
fad9b70b 373 Double_t minDCAxyDCAzPt[3] = { -5, -5, 0};
374 Double_t maxDCAxyDCAzPt[3] = { 5., 5., 100};
d0483ba3 375
fad9b70b 376 hDCAPtAll = new THnSparseF("hDCAPtAll","hDCAPtAll",3, binsDCAxyDCAzPt, minDCAxyDCAzPt, maxDCAxyDCAzPt);
377 hDCAPtAccepted = new THnSparseF("hDCAPtAccepted","hDCAPtAccepted",3, binsDCAxyDCAzPt, minDCAxyDCAzPt, maxDCAxyDCAzPt);
378 hMCDCAPtSecondary = new THnSparseF("hMCDCAPtSecondary","hMCDCAPtSecondary",3, binsDCAxyDCAzPt, minDCAxyDCAzPt, maxDCAxyDCAzPt);
379 hMCDCAPtPrimary = new THnSparseF("hMCDCAPtPrimary","hMCDCAPtPrimary",3, binsDCAxyDCAzPt, minDCAxyDCAzPt, maxDCAxyDCAzPt);
d0483ba3 380
381 hDCAPtAll->SetBinEdges(2, fBinsPt);
382 hDCAPtAccepted->SetBinEdges(2, fBinsPt);
383 hMCDCAPtSecondary->SetBinEdges(2, fBinsPt);
384 hMCDCAPtPrimary->SetBinEdges(2, fBinsPt);
d25bcbe6 385
fad9b70b 386 hDCAPtAll->Sumw2();
387 hDCAPtAccepted->Sumw2();
388 hMCDCAPtSecondary->Sumw2();
389 hMCDCAPtPrimary->Sumw2();
390
391 hDCAPtAll->GetAxis(0)->SetTitle("DCA_{xy} (cm)");
392 hDCAPtAll->GetAxis(1)->SetTitle("DCA_{z} (cm)");
393 hDCAPtAll->GetAxis(2)->SetTitle("p_{T} (GeV/c)");
394
395 hDCAPtAccepted->GetAxis(0)->SetTitle("DCA_{xy} (cm)");
396 hDCAPtAccepted->GetAxis(1)->SetTitle("DCA_{z} (cm)");
397 hDCAPtAccepted->GetAxis(2)->SetTitle("p_{T} (GeV/c)");
398
399 hMCDCAPtSecondary->GetAxis(0)->SetTitle("DCA_{xy} (cm)");
400 hMCDCAPtSecondary->GetAxis(1)->SetTitle("DCA_{z} (cm)");
401 hMCDCAPtSecondary->GetAxis(2)->SetTitle("p_{T} (GeV/c)");
402
403 hMCDCAPtPrimary->GetAxis(0)->SetTitle("DCA_{xy} (cm)");
404 hMCDCAPtPrimary->GetAxis(1)->SetTitle("DCA_{z} (cm)");
405 hMCDCAPtPrimary->GetAxis(2)->SetTitle("p_{T} (GeV/c)");
406
d25bcbe6 407 // Add Histos, Profiles etc to List
408 fOutputList->Add(hnZvPtEtaCent);
409 fOutputList->Add(hPt);
410 fOutputList->Add(hnMCRecPrimZvPtEtaCent);
411 fOutputList->Add(hnMCGenZvPtEtaCent);
5747f2c6 412 fOutputList->Add(hnMCRecSecZvPtEtaCent);
d25bcbe6 413 fOutputList->Add(hMCPt);
414 fOutputList->Add(hEventStatistics);
415 fOutputList->Add(hEventStatisticsCentrality);
416 fOutputList->Add(hAllEventStatisticsCentrality);
fad9b70b 417 fOutputList->Add(hEventStatisticsCentralityTrigger);
d25bcbe6 418 fOutputList->Add(hnZvMultCent);
419 fOutputList->Add(hTriggerStatistics);
420 fOutputList->Add(hMCTrackPdgCode);
421 fOutputList->Add(hMCTrackStatusCode);
422 fOutputList->Add(hCharge);
423 fOutputList->Add(hMCCharge);
424 fOutputList->Add(hMCPdgPt);
425 fOutputList->Add(hMCHijingPrim);
426 fOutputList->Add(hAccNclsTPC);
95f26ffa 427 fOutputList->Add(hAllCrossedRowsTPC);
428 fOutputList->Add(hFilterCrossedRowsTPC);
d25bcbe6 429 fOutputList->Add(hAccCrossedRowsTPC);
5747f2c6 430 fOutputList->Add(hDCAPtAll);
431 fOutputList->Add(hDCAPtAccepted);
432 fOutputList->Add(hMCDCAPtSecondary);
433 fOutputList->Add(hMCDCAPtPrimary);
fad9b70b 434
d25bcbe6 435
436 PostData(1, fOutputList);
437}
438
439void AlidNdPtAnalysisPbPbAOD::UserExec(Option_t *option)
440{
f856053f 441
d25bcbe6 442 // Main Loop
443 // called for each event
444 hEventStatistics->Fill("all events",1);
445
446 // set ZERO pointers:
447 AliInputEventHandler *inputHandler = NULL;
448 AliAODTrack *track = NULL;
449 AliAODMCParticle *mcPart = NULL;
450 AliAODMCHeader *mcHdr = NULL;
451 AliGenHijingEventHeader *genHijingHeader = NULL;
452 AliGenPythiaEventHeader *genPythiaHeader = NULL;
453
454 Bool_t bIsEventSelectedMB = kFALSE;
455 Bool_t bIsEventSelectedSemi = kFALSE;
456 Bool_t bIsEventSelectedCentral = kFALSE;
457 Bool_t bIsEventSelected = kFALSE;
458 Bool_t bIsPrimary = kFALSE;
459 Bool_t bIsHijingParticle = kFALSE;
460 Bool_t bIsPythiaParticle = kFALSE;
461 Bool_t bEventHasATrack = kFALSE;
462 Bool_t bEventHasATrackInRange = kFALSE;
463 Int_t nTriggerFired = 0;
464
465
466 Double_t dMCTrackZvPtEtaCent[4] = {0};
467 Double_t dTrackZvPtEtaCent[4] = {0};
468
469 Double_t dMCEventZv = -100;
470 Double_t dEventZv = -100;
471 Int_t iAcceptedMultiplicity = 0;
472
473 bIsMonteCarlo = kFALSE;
474
475 AliAODEvent *eventAOD = 0x0;
476 eventAOD = dynamic_cast<AliAODEvent*>( InputEvent() );
477 if (!eventAOD) {
478 AliWarning("ERROR: eventAOD not available \n");
479 return;
480 }
481
482 // check, which trigger has been fired
483 inputHandler = (AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
484 bIsEventSelectedMB = ( inputHandler->IsEventSelected() & AliVEvent::kMB);
485 bIsEventSelectedSemi = ( inputHandler->IsEventSelected() & AliVEvent::kSemiCentral);
486 bIsEventSelectedCentral = ( inputHandler->IsEventSelected() & AliVEvent::kCentral);
487
488 if(bIsEventSelectedMB || bIsEventSelectedSemi || bIsEventSelectedCentral) hTriggerStatistics->Fill("all triggered events",1);
489 if(bIsEventSelectedMB) { hTriggerStatistics->Fill("MB trigger",1); nTriggerFired++; }
490 if(bIsEventSelectedSemi) { hTriggerStatistics->Fill("SemiCentral trigger",1); nTriggerFired++; }
491 if(bIsEventSelectedCentral) { hTriggerStatistics->Fill("Central trigger",1); nTriggerFired++; }
492 if(nTriggerFired == 0) { hTriggerStatistics->Fill("No trigger",1); }
493
494 bIsEventSelected = ( inputHandler->IsEventSelected() & GetCollisionCandidates() );
495
496 // only take tracks of events, which are triggered
497 if(nTriggerFired == 0) { return; }
498
499 // if( !bIsEventSelected || nTriggerFired>1 ) return;
500
501 // hEventStatistics->Fill("events with only coll. cand.", 1);
502
503
504
505 // check if there is a stack, if yes, then do MC loop
506 TList *list = eventAOD->GetList();
507 TClonesArray *stack = 0x0;
508 stack = (TClonesArray*)list->FindObject(AliAODMCParticle::StdBranchName());
509
510 if( stack )
511 {
512 bIsMonteCarlo = kTRUE;
513
514 mcHdr = (AliAODMCHeader*)list->FindObject(AliAODMCHeader::StdBranchName());
515
516 genHijingHeader = GetHijingEventHeader(mcHdr);
517 // genPythiaHeader = GetPythiaEventHeader(mcHdr);
518
519 if(!genHijingHeader) { return; }
520
521 // if(!genPythiaHeader) { return; }
522
523 dMCEventZv = mcHdr->GetVtxZ();
524 dMCTrackZvPtEtaCent[0] = dMCEventZv;
525 hEventStatistics->Fill("MC all events",1);
526 }
527
528 AliCentrality* aCentrality = eventAOD->GetCentrality();
529 Double_t dCentrality = aCentrality->GetCentralityPercentile("V0M");
530
531 if( dCentrality < 0 ) return;
532 hEventStatistics->Fill("after centrality selection",1);
533
534
535
536 // start with MC truth analysis
537 if(bIsMonteCarlo)
538 {
539
540 if( dMCEventZv > dCutMaxZVertex ) { return; }
541
542 dMCTrackZvPtEtaCent[0] = dMCEventZv;
543
544 hEventStatistics->Fill("MC afterZv cut",1);
545
546 for(Int_t iMCtrack = 0; iMCtrack < stack->GetEntriesFast(); iMCtrack++)
547 {
548 mcPart =(AliAODMCParticle*)stack->At(iMCtrack);
549
550 // check for charge
551 if( !(IsMCTrackAccepted(mcPart)) ) continue;
552
553 if(!IsHijingParticle(mcPart, genHijingHeader)) { continue; }
554
555 if(mcPart->IsPhysicalPrimary() )
556 {
557 hMCHijingPrim->Fill("IsPhysicalPrimary",1);
558 }
559 else
560 {
561 hMCHijingPrim->Fill("NOT a primary",1);
562 continue;
563 }
564
d0483ba3 565
d25bcbe6 566 //
567 // ======================== fill histograms ========================
568 dMCTrackZvPtEtaCent[1] = mcPart->Pt();
569 dMCTrackZvPtEtaCent[2] = mcPart->Eta();
570 dMCTrackZvPtEtaCent[3] = dCentrality;
571
572 bEventHasATrack = kTRUE;
573
574 hnMCGenZvPtEtaCent->Fill(dMCTrackZvPtEtaCent);
575
576 if( (dMCTrackZvPtEtaCent[1] > dCutPtMin) &&
577 (dMCTrackZvPtEtaCent[1] < dCutPtMax) &&
578 (dMCTrackZvPtEtaCent[2] > dCutEtaMin) &&
579 (dMCTrackZvPtEtaCent[2] < dCutEtaMax) )
580 {
581 hMCPt->Fill(mcPart->Pt());
582 hMCCharge->Fill(mcPart->Charge()/3.);
583 bEventHasATrackInRange = kTRUE;
584 }
585
586 }
587 } // isMonteCarlo
588 if(bEventHasATrack) { hEventStatistics->Fill("MC events with tracks",1); }
589 if(bEventHasATrackInRange) { hEventStatistics->Fill("MC events with tracks in range",1); }
590 bEventHasATrack = kFALSE;
591 bEventHasATrackInRange = kFALSE;
592
593
594
595 // Loop over recontructed tracks
596
597 dEventZv = eventAOD->GetPrimaryVertex()->GetZ();
598 if( TMath::Abs(dEventZv) > dCutMaxZVertex ) return;
599
600 hAllEventStatisticsCentrality->Fill(dCentrality/*, nTriggerFired*/);
601
602 hEventStatistics->Fill("after Zv cut",1);
603
604 dTrackZvPtEtaCent[0] = dEventZv;
605
606 for(Int_t itrack = 0; itrack < eventAOD->GetNumberOfTracks(); itrack++)
607 {
608 track = eventAOD->GetTrack(itrack);
609 if(!track) continue;
610
611 mcPart = NULL;
612 dMCTrackZvPtEtaCent[1] = 0;
613 dMCTrackZvPtEtaCent[2] = 0;
614 dMCTrackZvPtEtaCent[3] = 0;
615
616 bIsPrimary = kFALSE;
617
d0483ba3 618 Double_t dDCAxyDCAzPt[3] = { GetDCAxy(track, eventAOD), GetDCAz(track, eventAOD), track->Pt() };
619
620 hDCAPtAll->Fill(dDCAxyDCAzPt);
621
d25bcbe6 622 if( !(IsTrackAccepted(track)) ) continue;
623
5747f2c6 624 dTrackZvPtEtaCent[1] = track->Pt();
625 dTrackZvPtEtaCent[2] = track->Eta();
626 dTrackZvPtEtaCent[3] = dCentrality;
627
d25bcbe6 628 if( bIsMonteCarlo )
629 {
630 mcPart = (AliAODMCParticle*)stack->At(TMath::Abs(track->GetLabel()));
631 if( !mcPart ) { continue; }
632
633 // check for charge
634 if( !(IsMCTrackAccepted(mcPart)) ) { continue; }
635
636 bIsHijingParticle = IsHijingParticle(mcPart, genHijingHeader);
637 // bIsPythiaParticle = IsPythiaParticle(mcPart, genPythiaHeader);
638
639 // if(!bIsHijingParticle) continue; // only take real tracks, not injected ones
640
641 bIsPrimary = mcPart->IsPhysicalPrimary();
642
643 dMCTrackZvPtEtaCent[1] = mcPart->Pt();
644 dMCTrackZvPtEtaCent[2] = mcPart->Eta();
645 dMCTrackZvPtEtaCent[3] = dCentrality;
646
647 if(bIsPrimary && bIsHijingParticle)
648 {
d0483ba3 649 hnMCRecPrimZvPtEtaCent->Fill(dMCTrackZvPtEtaCent);
650 hMCDCAPtPrimary->Fill(dDCAxyDCAzPt);
d25bcbe6 651 }
652
653 if(!bIsPrimary /*&& !bIsHijingParticle*/)
654 {
655 Int_t indexMoth = mcPart->GetMother();
656 if(indexMoth >= 0)
657 {
658 AliAODMCParticle* moth = (AliAODMCParticle*)stack->At(indexMoth);
659 Bool_t bMotherIsHijingParticle = IsHijingParticle(moth, genHijingHeader);
660
661 if(bMotherIsHijingParticle) // only store secondaries, which come from a not embedded signal!
662 {
663 hMCTrackStatusCode->Fill(Form("%d",mcPart->GetStatus()), 1);
664 if(TMath::Abs(mcPart->Eta()) < 0.8) { hMCPdgPt->Fill(mcPart->Pt(), Form("%s",GetParticleName(mcPart->GetPdgCode())), 1); }
665
d0483ba3 666 hnMCRecSecZvPtEtaCent->Fill(dMCTrackZvPtEtaCent);
667 hMCDCAPtSecondary->Fill(dDCAxyDCAzPt);
d25bcbe6 668 hMCTrackPdgCode->Fill(Form("%s_H%i_H%i",GetParticleName(moth->GetPdgCode()),bMotherIsHijingParticle, bIsHijingParticle), 1);
669 // delete moth;
670 }
671 }
672 }
673 } // end isMonteCarlo
674
675 // ======================== fill histograms ========================
676
d0483ba3 677 // if(bIsMonteCarlo && !bIsHijingParticle)
678 // {
679 // continue; //only store reco tracks, which do not come from embedded signal
680 // }
681
682
683
684 bEventHasATrack = kTRUE;
685
686 hnZvPtEtaCent->Fill(dTrackZvPtEtaCent);
687 hDCAPtAccepted->Fill(dDCAxyDCAzPt);
688
689 if( (dTrackZvPtEtaCent[1] > dCutPtMin) &&
690 (dTrackZvPtEtaCent[1] < dCutPtMax) &&
691 (dTrackZvPtEtaCent[2] > dCutEtaMin) &&
692 (dTrackZvPtEtaCent[2] < dCutEtaMax) )
693 {
694 iAcceptedMultiplicity++;
695 bEventHasATrackInRange = kTRUE;
696 hPt->Fill(track->Pt());
697 hCharge->Fill(track->Charge());
698 }
d25bcbe6 699 } // end track loop
700
701 if(bEventHasATrack) { hEventStatistics->Fill("events with tracks",1); bEventHasATrack = kFALSE; }
702
703 if(bEventHasATrackInRange)
704 {
705 hEventStatistics->Fill("events with tracks in range",1);
706 hEventStatisticsCentrality->Fill(dCentrality);
707 bEventHasATrackInRange = kFALSE;
fad9b70b 708
709 if(bIsEventSelectedMB) hEventStatisticsCentralityTrigger->Fill(dCentrality, 0);
710 if(bIsEventSelectedSemi) hEventStatisticsCentralityTrigger->Fill(dCentrality, 1);
711 if(bIsEventSelectedCentral) hEventStatisticsCentralityTrigger->Fill(dCentrality, 2);
d25bcbe6 712 }
713
714 Double_t dEventZvMultCent[3] = {dEventZv, iAcceptedMultiplicity, dCentrality};
715 hnZvMultCent->Fill(dEventZvMultCent);
716
717
718
719 PostData(1, fOutputList);
720
721}
722
723Bool_t AlidNdPtAnalysisPbPbAOD::IsTrackAccepted(AliAODTrack *tr)
724{
725 if(!tr) return kFALSE;
726
727 if(tr->Charge()==0) { return kFALSE; }
728
d25bcbe6 729 Double_t dNClustersTPC = tr->GetTPCNcls();
730 Double_t dCrossedRowsTPC = tr->GetTPCClusterInfo(2,1);
731
95f26ffa 732 hAllCrossedRowsTPC->Fill(dCrossedRowsTPC);
733
734 if(!(tr->TestFilterBit(AliAODTrack::kTrkGlobal)) ) { return kFALSE; }
735
736 hFilterCrossedRowsTPC->Fill(dCrossedRowsTPC);
737
f856053f 738 if(dCrossedRowsTPC < GetCutMinNCrossedRowsTPC()) { return kFALSE; }
d0483ba3 739
d25bcbe6 740 hAccNclsTPC->Fill(dNClustersTPC);
741 hAccCrossedRowsTPC->Fill(dCrossedRowsTPC);
742 // Double_t dFindableClustersTPC = tr->GetTPCNclsF();
743 // Double_t dChi2PerClusterTPC = (dNClustersTPC>0)?tr->Chi2perNDF()*(dNClustersTPC-5)/dNClustersTPC:-1.; // see AliDielectronVarManager.h
744 //
745 // Bool_t bIsFromKink = kFALSE;
746 // if(tr->GetProdVertex()->GetType() == AliAODVertex::kKink) bIsFromKink = kTRUE;
747 //
748 // // from AliAnalysisTaskPIDqa.cxx
749 // ULong_t uStatus = tr->GetStatus();
750 // Bool_t bHasRefitTPC = kFALSE;
751 // Bool_t bHasRefitITS = kFALSE;
752 //
753 // if ((uStatus & AliVTrack::kTPCrefit) == AliVTrack::kTPCrefit) bHasRefitTPC = kTRUE;
754 // if ((uStatus & AliVTrack::kITSrefit) == AliVTrack::kITSrefit) bHasRefitITS = kTRUE;
755 //
756 // // from AliDielectronVarManager.h
757 // Bool_t bHasHitInSPD = kFALSE;
758 // for (Int_t iC=0; iC<2; iC++)
759 // {
760 // if (((tr->GetITSClusterMap()) & (1<<(iC))) > 0) { bHasHitInSPD = kTRUE; }
761 // }
762 //
763 // Double_t dNClustersITS = tr->GetITSNcls();
764
765 // cuts to be done:
766 // TPC
767 // esdTrackCuts->SetMinNCrossedRowsTPC(70);
768 // esdTrackCuts->SetMinRatioCrossedRowsOverFindableClustersTPC(0.8);
769 //
770 // esdTrackCuts->SetMaxChi2PerClusterTPC(4);
771 // esdTrackCuts->SetAcceptKinkDaughters(kFALSE);
772 // esdTrackCuts->SetRequireTPCRefit(kTRUE);
773 // ITS
774 // esdTrackCuts->SetRequireITSRefit(kTRUE);
775 // esdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD, AliESDtrackCuts::kAny);
776 //
777 // esdTrackCuts->SetMaxDCAToVertexXYPtDep("0.0105+0.0350/pt^1.1");
778 // esdTrackCuts->SetMaxChi2TPCConstrainedGlobal(36);
779 //
780 // esdTrackCuts->SetMaxDCAToVertexZ(2);
781 // esdTrackCuts->SetDCAToVertex2D(kFALSE);
782 // esdTrackCuts->SetRequireSigmaToVertex(kFALSE);
783 //
784 // esdTrackCuts->SetMaxChi2PerClusterITS(36);
785
786
787 return kTRUE;
788}
789
d0483ba3 790Double_t AlidNdPtAnalysisPbPbAOD::GetDCAz(AliAODTrack *track, AliAODEvent *event)
791{
792 return GetDCA(track, event, kTRUE);
793}
794
795Double_t AlidNdPtAnalysisPbPbAOD::GetDCAxy(AliAODTrack *track, AliAODEvent *event)
796{
797 return GetDCA(track, event, kFALSE);
798}
799
800Double_t AlidNdPtAnalysisPbPbAOD::GetDCA(AliAODTrack *tr, AliAODEvent *evt, Bool_t bDCAz)
801{
802 if(!tr) return -999.;
803
804 if(tr->TestBit(AliAODTrack::kIsDCA))
805 {
806 if(bDCAz) return tr->ZAtDCA();
807 else return sqrt(tr->XAtDCA()*tr->XAtDCA() + tr->YAtDCA()*tr->YAtDCA());
808 }
809
810 Bool_t ok=kFALSE;
811 Double_t d0z0[2];
812 if(evt)
813 {
814 Double_t covd0z0[3];
815
816 AliExternalTrackParam etp;
817 etp.CopyFromVTrack(tr);
818
819 Float_t xstart = etp.GetX();
820 if(xstart>3.)
821 {
822 d0z0[0]=-999.;
823 d0z0[1]=-999.;
824 //printf("This method can be used only for propagation inside the beam pipe \n");
825 if(bDCAz) return d0z0[1];
826 else return d0z0[0];
827 }
fad9b70b 828
d0483ba3 829 AliAODVertex *vtx =(AliAODVertex*)(evt->GetPrimaryVertex());
830 Double_t fBzkG = evt->GetMagneticField(); // z componenent of field in kG
831 ok = etp.PropagateToDCA(vtx,fBzkG,kVeryBig,d0z0,covd0z0);
832 }
fad9b70b 833
d0483ba3 834 if(!ok){
835 d0z0[0]=-999.;
836 d0z0[1]=-999.;
837 if(bDCAz) return d0z0[1];
838 else return d0z0[0];
839 }
840 if(bDCAz) return d0z0[1];
841 else return d0z0[0];
842}
843
d25bcbe6 844Bool_t AlidNdPtAnalysisPbPbAOD::IsMCTrackAccepted(AliAODMCParticle *part)
845{
846 if(!part) return kFALSE;
847
848 Double_t charge = part->Charge()/3.;
849 if (TMath::Abs(charge) < 0.001) return kFALSE;
850
851 return kTRUE;
852}
853
854const char * AlidNdPtAnalysisPbPbAOD::GetParticleName(Int_t pdg)
855{
856 TParticlePDG * p1 = TDatabasePDG::Instance()->GetParticle(pdg);
857 if(p1) return p1->GetName();
858 return Form("%d", pdg);
859}
860
861AliGenHijingEventHeader* AlidNdPtAnalysisPbPbAOD::GetHijingEventHeader(AliAODMCHeader *header)
862{
863 //
864 // inspired by PWGJE/AliPWG4HighPtSpectra.cxx
865 //
866
867 if(!header) return 0x0;
868 AliGenHijingEventHeader* hijingGenHeader = NULL;
869
870 TList* headerList = header->GetCocktailHeaders();
871
872 for(Int_t i = 0; i < headerList->GetEntries(); i++)
873 {
874 hijingGenHeader = dynamic_cast<AliGenHijingEventHeader*>(headerList->At(i));
875 if(hijingGenHeader) break;
876 }
877
878 if(!hijingGenHeader) return 0x0;
879
880 return hijingGenHeader;
881}
882
883AliGenPythiaEventHeader* AlidNdPtAnalysisPbPbAOD::GetPythiaEventHeader(AliAODMCHeader *header)
884{
885 //
886 // inspired by PWGJE/AliPWG4HighPtSpectra.cxx
887 //
888
889 if(!header) return 0x0;
890 AliGenPythiaEventHeader* PythiaGenHeader = NULL;
891
892 TList* headerList = header->GetCocktailHeaders();
893
894 for(Int_t i = 0; i < headerList->GetEntries(); i++)
895 {
896 PythiaGenHeader = dynamic_cast<AliGenPythiaEventHeader*>(headerList->At(i));
897 if(PythiaGenHeader) break;
898 }
899
900 if(!PythiaGenHeader) return 0x0;
901
902 return PythiaGenHeader;
903}
904
905//________________________________________________________________________
906Bool_t AlidNdPtAnalysisPbPbAOD::IsHijingParticle(const AliAODMCParticle *part, AliGenHijingEventHeader* hijingGenHeader){
907
908 // Check whether a particle is from Hijing or some injected
909
910 if(part->Label() > (hijingGenHeader->NProduced()-1)) return kFALSE;
911 return kTRUE;
912}
913
914//________________________________________________________________________
915Bool_t AlidNdPtAnalysisPbPbAOD::IsPythiaParticle(const AliAODMCParticle *part, AliGenPythiaEventHeader* pythiaGenHeader){
916
917 // Check whether a particle is from Pythia or some injected
918
919 if(part->Label() > (pythiaGenHeader->NProduced()-1)) return kFALSE;
920 return kTRUE;
921}
922
5747f2c6 923// Int_t AlidNdPtAnalysisPbPbAOD::IsMCSecondary(AliAODMCParticle *part, TClonesArray *arrayMC)
924// {
d0483ba3 925 // //
926 // // adapted from AliAnalysisTaskSpectraAOD.cxx
927 // //
928 // // returns
929 // // -1: no particle
930 // // 0: is primary
931 // // 1: is secondary from weak
932 // // 2: is secondary from material
933 //
934 // // usage for studies, currrently not implemented
935 //
936 // if(!part) return -1;
937 //
938 // if( part->IsPhysicalPrimary() ) return 0;
939 //
940 // Bool_t isSecondaryMaterial = kFALSE;
941 // Bool_t isSecondaryWeak = kFALSE;
942 // Int_t mfl = -999;
943 // Int_t codemoth = -999;
944 // Int_t indexMoth = part->GetMother(); // FIXME ignore fakes? TO BE CHECKED, on ESD is GetFirstMother()
945 // if(indexMoth >= 0)
946 // {
947 // AliAODMCParticle* moth = (AliAODMCParticle*) arrayMC->At(indexMoth);
948 // codemoth = TMath::Abs(moth->GetPdgCode());
949 // mfl = Int_t (codemoth/ TMath::Power(10, Int_t(TMath::Log10(codemoth))));
950 // }
951 // // add if(partMC->GetStatus() & kPDecay)? FIXME
952 // if(mfl==3) isSecondaryWeak = kTRUE;
953 // else isSecondaryMaterial = kTRUE;
954 //
955 // if(isSecondaryWeak) return 1;
956 // if(isSecondaryMaterial) return 2;
957 //
958 // // if( isSecondaryMaterial || isSecondaryWeak ) return kTRUE;
959 //
960 // // return kFALSE; this line will not be reached, as either isSecondaryMaterial or isSecondaryWeak is true!
961 // // removed due to coverity
962 // }
963
964
965
966 void AlidNdPtAnalysisPbPbAOD::Terminate(Option_t *)
967 {
968
969 }
970
971 Double_t* AlidNdPtAnalysisPbPbAOD::GetArrayClone(Int_t n, Double_t* source)
972 {
973 if (!source || n==0) return 0;
974 Double_t* dest = new Double_t[n];
975 for (Int_t i=0; i<n ; i++) { dest[i] = source[i]; }
976 return dest;
977 }
978