bug fix: deletion of pointers, added: store cut settings
[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
ea3cfeda 19// last modified: 17.10.2013
d25bcbe6 20//------------------------------------------------------------------------------
21
22
23#include "AlidNdPtAnalysisPbPbAOD.h"
24
8a4ab847 25#include "AliAnalysisTaskSE.h"
d25bcbe6 26
27using namespace std;
28
29ClassImp(AlidNdPtAnalysisPbPbAOD)
30
d25bcbe6 31AlidNdPtAnalysisPbPbAOD::AlidNdPtAnalysisPbPbAOD(const char *name) : AliAnalysisTaskSE(name),
32fOutputList(0),
33// Histograms
ea3cfeda
PL
34fPt(0),
35fMCPt(0),
36fZvPtEtaCent(0),
bc2a9da9 37fPhiPtEtaCent(0),
a0036e80 38fPtResptCent(0),
ea3cfeda
PL
39fMCRecPrimZvPtEtaCent(0),
40fMCGenZvPtEtaCent(0),
41fMCRecSecZvPtEtaCent(0),
bc2a9da9
PL
42fMCRecPrimPhiPtEtaCent(0),
43fMCGenPhiPtEtaCent(0),
44fMCRecSecPhiPtEtaCent(0),
ea3cfeda
PL
45fEventStatistics(0),
46fEventStatisticsCentrality(0),
47fMCEventStatisticsCentrality(0),
48fAllEventStatisticsCentrality(0),
49fEventStatisticsCentralityTrigger(0),
50fZvMultCent(0),
51fTriggerStatistics(0),
52fMCTrackPdgCode(0),
53fMCTrackStatusCode(0),
54fCharge(0),
55fMCCharge(0),
56fMCPdgPt(0),
57fMCHijingPrim(0),
58fDCAPtAll(0),
59fDCAPtAccepted(0),
60fMCDCAPtSecondary(0),
61fMCDCAPtPrimary(0),
62fCutPercClusters(0),
63fCutPercCrossed(0),
64fCrossCheckRowsLength(0),
65fCrossCheckClusterLength(0),
66fCrossCheckRowsLengthAcc(0),
67fCrossCheckClusterLengthAcc(0),
8a4ab847 68fCutSettings(0),
d25bcbe6 69//global
ea3cfeda 70fIsMonteCarlo(0),
d25bcbe6 71// event cut variables
ea3cfeda 72fCutMaxZVertex(10.),
d25bcbe6 73// track kinematic cut variables
ea3cfeda
PL
74fCutPtMin(0.15),
75fCutPtMax(200.),
76fCutEtaMin(-0.8),
77fCutEtaMax(0.8),
d25bcbe6 78// track quality cut variables
a0036e80 79fFilterBit(AliAODTrack::kTrkGlobal),
ea3cfeda
PL
80fUseRelativeCuts(kFALSE),
81fCutRequireTPCRefit(kTRUE),
a0036e80 82fCutRequireITSRefit(kTRUE),
ea3cfeda
PL
83fCutMinNumberOfClusters(60),
84fCutPercMinNumberOfClusters(0.2),
85fCutMinNumberOfCrossedRows(120.),
86fCutPercMinNumberOfCrossedRows(0.2),
87fCutMinRatioCrossedRowsOverFindableClustersTPC(0.8),
88fCutMaxChi2PerClusterTPC(4.),
89fCutMaxFractionSharedTPCClusters(0.4),
90fCutMaxDCAToVertexZ(3.0),
91fCutMaxDCAToVertexXY(3.0),
ea3cfeda
PL
92fCutMaxChi2PerClusterITS(36.),
93fCutDCAToVertex2D(kFALSE),
94fCutRequireSigmaToVertex(kFALSE),
95fCutMaxDCAToVertexXYPtDepPar0(0.0182),
96fCutMaxDCAToVertexXYPtDepPar1(0.0350),
97fCutMaxDCAToVertexXYPtDepPar2(1.01),
98fCutAcceptKinkDaughters(kFALSE),
99fCutMaxChi2TPCConstrainedGlobal(36.),
100fCutLengthInTPCPtDependent(kFALSE),
101fPrefactorLengthInTPCPtDependent(1),
d25bcbe6 102// binning for THnSparse
103fMultNbins(0),
104fPtNbins(0),
105fPtCorrNbins(0),
72bb4ceb 106fPtCheckNbins(0),
d25bcbe6 107fEtaNbins(0),
72bb4ceb 108fEtaCheckNbins(0),
d25bcbe6 109fZvNbins(0),
110fCentralityNbins(0),
72bb4ceb 111fPhiNbins(0),
d25bcbe6 112fBinsMult(0),
113fBinsPt(0),
114fBinsPtCorr(0),
72bb4ceb 115fBinsPtCheck(0),
d25bcbe6 116fBinsEta(0),
72bb4ceb 117fBinsEtaCheck(0),
d25bcbe6 118fBinsZv(0),
72bb4ceb 119fBinsCentrality(0),
120fBinsPhi(0)
d25bcbe6 121{
72bb4ceb 122
123 for(Int_t i = 0; i < cqMax; i++)
124 {
ea3cfeda
PL
125 fCrossCheckAll[i] = 0;
126 fCrossCheckAcc[i] = 0;
72bb4ceb 127 }
128
d25bcbe6 129 fMultNbins = 0;
130 fPtNbins = 0;
131 fPtCorrNbins = 0;
72bb4ceb 132 fPtCheckNbins = 0;
d25bcbe6 133 fEtaNbins = 0;
72bb4ceb 134 fEtaCheckNbins = 0;
d25bcbe6 135 fZvNbins = 0;
136 fCentralityNbins = 0;
137 fBinsMult = 0;
138 fBinsPt = 0;
139 fBinsPtCorr = 0;
72bb4ceb 140 fBinsPtCheck = 0;
d25bcbe6 141 fBinsEta = 0;
72bb4ceb 142 fBinsEtaCheck = 0;
d25bcbe6 143 fBinsZv = 0;
144 fBinsCentrality = 0;
72bb4ceb 145 fBinsPhi = 0;
d25bcbe6 146
147 DefineOutput(1, TList::Class());
148}
149
150// destructor
151AlidNdPtAnalysisPbPbAOD::~AlidNdPtAnalysisPbPbAOD()
8a4ab847 152{
ea3cfeda 153 if(fZvPtEtaCent) delete fZvPtEtaCent; fZvPtEtaCent = 0;
8a4ab847 154 if(fPhiPtEtaCent) delete fPhiPtEtaCent; fPhiPtEtaCent = 0;
155 if(fPtResptCent) delete fPtResptCent; fPtResptCent = 0;
ea3cfeda 156 if(fPt) delete fPt; fPt = 0;
8a4ab847 157 if(fMCPt) delete fMCPt; fMCPt = 0;
158
ea3cfeda
PL
159 if(fMCRecPrimZvPtEtaCent) delete fMCRecPrimZvPtEtaCent; fMCRecPrimZvPtEtaCent = 0;
160 if(fMCGenZvPtEtaCent) delete fMCGenZvPtEtaCent; fMCGenZvPtEtaCent = 0;
161 if(fMCRecSecZvPtEtaCent) delete fMCRecSecZvPtEtaCent; fMCRecSecZvPtEtaCent = 0;
8a4ab847 162
163 if(fMCRecPrimPhiPtEtaCent) delete fMCRecPrimPhiPtEtaCent; fMCRecPrimPhiPtEtaCent = 0;
164 if(fMCGenPhiPtEtaCent) delete fMCGenPhiPtEtaCent; fMCGenPhiPtEtaCent = 0;
165 if(fMCRecSecPhiPtEtaCent) delete fMCRecSecPhiPtEtaCent; fMCRecSecPhiPtEtaCent = 0;
166
ea3cfeda
PL
167 if(fEventStatistics) delete fEventStatistics; fEventStatistics = 0;
168 if(fEventStatisticsCentrality) delete fEventStatisticsCentrality; fEventStatisticsCentrality = 0;
169 if(fMCEventStatisticsCentrality) delete fMCEventStatisticsCentrality; fMCEventStatisticsCentrality = 0;
170 if(fAllEventStatisticsCentrality) delete fAllEventStatisticsCentrality; fAllEventStatisticsCentrality = 0;
171 if(fEventStatisticsCentralityTrigger) delete fEventStatisticsCentralityTrigger; fEventStatisticsCentralityTrigger = 0;
172 if(fZvMultCent) delete fZvMultCent; fZvMultCent = 0;
173 if(fTriggerStatistics) delete fTriggerStatistics; fTriggerStatistics = 0;
174 if(fMCTrackPdgCode) delete fMCTrackPdgCode; fMCTrackPdgCode = 0;
175 if(fMCTrackStatusCode) delete fMCTrackStatusCode; fMCTrackStatusCode = 0;
176 if(fCharge) delete fCharge; fCharge = 0;
177 if(fMCCharge) delete fMCCharge; fMCCharge = 0;
178 if(fMCPdgPt) delete fMCPdgPt; fMCPdgPt = 0;
179 if(fMCHijingPrim) delete fMCHijingPrim; fMCHijingPrim = 0;
180 if(fDCAPtAll) delete fDCAPtAll; fDCAPtAll = 0;
181 if(fDCAPtAccepted) delete fDCAPtAccepted; fDCAPtAccepted = 0;
182 if(fMCDCAPtSecondary) delete fMCDCAPtSecondary; fMCDCAPtSecondary = 0;
183 if(fMCDCAPtPrimary) delete fMCDCAPtPrimary; fMCDCAPtPrimary = 0;
72bb4ceb 184
185 for(Int_t i = 0; i < cqMax; i++)
186 {
ea3cfeda
PL
187 if(fCrossCheckAll[i]) delete fCrossCheckAll[i]; fCrossCheckAll[i] = 0;
188 if(fCrossCheckAcc[i]) delete fCrossCheckAcc[i]; fCrossCheckAcc[i] = 0;
72bb4ceb 189 }
190
8a4ab847 191 if(fCutPercClusters) delete fCutPercClusters; fCutPercClusters = 0;
192 if(fCutPercCrossed) delete fCutPercCrossed; fCutPercCrossed = 0;
193 if(fCrossCheckRowsLength) delete fCrossCheckRowsLength; fCrossCheckRowsLength = 0;
194 if(fCrossCheckClusterLength) delete fCrossCheckClusterLength; fCrossCheckClusterLength = 0;
195 if(fCrossCheckRowsLengthAcc) delete fCrossCheckRowsLengthAcc; fCrossCheckRowsLengthAcc = 0;
196 if(fCrossCheckClusterLengthAcc) delete fCrossCheckClusterLengthAcc; fCrossCheckClusterLengthAcc = 0;
197 if(fCutSettings) delete fCutSettings; fCutSettings = 0;
198
199 if(fOutputList)
200 {
201 fOutputList->Clear();
202 delete fOutputList;
203 }
204 fOutputList = 0;
d25bcbe6 205}
206
207void AlidNdPtAnalysisPbPbAOD::UserCreateOutputObjects()
208{
209 // create all output histograms here
210 OpenFile(1, "RECREATE");
211
212 fOutputList = new TList();
213 fOutputList->SetOwner();
214
215 //define default binning
216 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 };
217 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};
72bb4ceb 218 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};
d25bcbe6 219 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};
220 Double_t binsZvDefault[13] = {-30.,-25.,-20.,-15.,-10.,-5.,0.,5.,10.,15.,20.,25.,30.};
221 Double_t binsCentralityDefault[12] = {0., 5., 10., 20., 30., 40., 50., 60., 70., 80., 90., 100.};
222
72bb4ceb 223 Double_t binsPhiDefault[37] = { 0, 0.174533, 0.349066, 0.523599, 0.698132, 0.872665, 1.0472, 1.22173, 1.39626, 1.5708, 1.74533, 1.91986, 2.0944, 2.26893, 2.44346, 2.61799, 2.79253, 2.96706, 3.14159, 3.31613, 3.49066, 3.66519, 3.83972, 4.01426, 4.18879, 4.36332, 4.53786, 4.71239, 4.88692, 5.06145, 5.23599, 5.41052, 5.58505, 5.75959, 5.93412, 6.10865, 2.*TMath::Pi()};
224
225 Double_t binsPtCheckDefault[20] = {0.,0.15,0.5,1.0,2.0,3.0,4.0, 5.0, 10.0, 13.0, 15.0, 20.0, 25.0, 30.0, 40.0, 50.0, 70.0, 100.0, 150.0, 200.0};
226 Double_t binsEtaCheckDefault[7] = {-1.0,-0.8,-0.4,0.,0.4,0.8,1.0};
227
d25bcbe6 228 // if no binning is set, use the default
229 if (!fBinsMult) { SetBinsMult(48,binsMultDefault); }
230 if (!fBinsPt) { SetBinsPt(82,binsPtDefault); }
231 if (!fBinsPtCorr) { SetBinsPtCorr(37,binsPtCorrDefault); }
72bb4ceb 232 if (!fBinsPtCheck) { SetBinsPtCheck(20,binsPtCheckDefault); }
d25bcbe6 233 if (!fBinsEta) { SetBinsEta(31,binsEtaDefault); }
72bb4ceb 234 if (!fBinsEtaCheck) { SetBinsEtaCheck(7,binsEtaCheckDefault); }
d25bcbe6 235 if (!fBinsZv) { SetBinsZv(13,binsZvDefault); }
236 if (!fBinsCentrality) { SetBinsCentrality(12,binsCentralityDefault); }
72bb4ceb 237 if (!fBinsPhi) { SetBinsPhi(37,binsPhiDefault); }
d25bcbe6 238
239 Int_t binsZvPtEtaCent[4]={fZvNbins-1,fPtNbins-1,fEtaNbins-1,fCentralityNbins-1};
bc2a9da9 240 Int_t binsPhiPtEtaCent[4]={fPhiNbins-1,fPtNbins-1,fEtaNbins-1,fCentralityNbins-1};
d25bcbe6 241 Int_t binsZvMultCent[3]={fZvNbins-1,fMultNbins-1,fCentralityNbins-1};
242
a0036e80 243 Int_t binsOneOverPtPtResCent[3]={400,300,11};
244 Double_t minbinsOneOverPtPtResCent[3]={0,0,0};
245 Double_t maxbinsOneOverPtPtResCent[3]={1,0.015,100};
246
d25bcbe6 247 // define Histograms
ea3cfeda
PL
248 fZvPtEtaCent = new THnSparseF("fZvPtEtaCent","Zv:Pt:Eta:Centrality",4,binsZvPtEtaCent);
249 fZvPtEtaCent->SetBinEdges(0,fBinsZv);
250 fZvPtEtaCent->SetBinEdges(1,fBinsPt);
251 fZvPtEtaCent->SetBinEdges(2,fBinsEta);
252 fZvPtEtaCent->SetBinEdges(3,fBinsCentrality);
253 fZvPtEtaCent->GetAxis(0)->SetTitle("Zv (cm)");
254 fZvPtEtaCent->GetAxis(1)->SetTitle("Pt (GeV/c)");
255 fZvPtEtaCent->GetAxis(2)->SetTitle("Eta");
256 fZvPtEtaCent->GetAxis(3)->SetTitle("Centrality");
257 fZvPtEtaCent->Sumw2();
258
bc2a9da9
PL
259 fPhiPtEtaCent = new THnSparseF("fPhiPtEtaCent","Phi:Pt:Eta:Centrality",4,binsPhiPtEtaCent);
260 fPhiPtEtaCent->SetBinEdges(0,fBinsPhi);
261 fPhiPtEtaCent->SetBinEdges(1,fBinsPt);
262 fPhiPtEtaCent->SetBinEdges(2,fBinsEta);
263 fPhiPtEtaCent->SetBinEdges(3,fBinsCentrality);
8a4ab847 264 fPhiPtEtaCent->GetAxis(0)->SetTitle("Phi");
bc2a9da9
PL
265 fPhiPtEtaCent->GetAxis(1)->SetTitle("Pt (GeV/c)");
266 fPhiPtEtaCent->GetAxis(2)->SetTitle("Eta");
267 fPhiPtEtaCent->GetAxis(3)->SetTitle("Centrality");
268 fPhiPtEtaCent->Sumw2();
269
a0036e80 270 fPtResptCent = new THnSparseF("fPtResptCent","OneOverPt:PtRes:Centrality",3,binsOneOverPtPtResCent, minbinsOneOverPtPtResCent, maxbinsOneOverPtPtResCent);
271 fPtResptCent->SetBinEdges(2, fBinsCentrality);
272 fPtResptCent->GetAxis(0)->SetTitle("1/pT (GeV/c)^{-1}");
273 fPtResptCent->GetAxis(1)->SetTitle("#sigma(1/pT)");
274 fPtResptCent->GetAxis(2)->SetTitle("centrality");
275 fPtResptCent->Sumw2();
276
bc2a9da9 277
ea3cfeda
PL
278 fMCRecPrimZvPtEtaCent = new THnSparseF("fMCRecPrimZvPtEtaCent","mcZv:mcPt:mcEta:Centrality",4,binsZvPtEtaCent);
279 fMCRecPrimZvPtEtaCent->SetBinEdges(0,fBinsZv);
280 fMCRecPrimZvPtEtaCent->SetBinEdges(1,fBinsPt);
281 fMCRecPrimZvPtEtaCent->SetBinEdges(2,fBinsEta);
282 fMCRecPrimZvPtEtaCent->SetBinEdges(3,fBinsCentrality);
283 fMCRecPrimZvPtEtaCent->GetAxis(0)->SetTitle("MC Zv (cm)");
284 fMCRecPrimZvPtEtaCent->GetAxis(1)->SetTitle("MC Pt (GeV/c)");
285 fMCRecPrimZvPtEtaCent->GetAxis(2)->SetTitle("MC Eta");
286 fMCRecPrimZvPtEtaCent->GetAxis(3)->SetTitle("Centrality");
287 fMCRecPrimZvPtEtaCent->Sumw2();
288
289 fMCGenZvPtEtaCent = new THnSparseF("fMCGenZvPtEtaCent","mcZv:mcPt:mcEta:Centrality",4,binsZvPtEtaCent);
290 fMCGenZvPtEtaCent->SetBinEdges(0,fBinsZv);
291 fMCGenZvPtEtaCent->SetBinEdges(1,fBinsPt);
292 fMCGenZvPtEtaCent->SetBinEdges(2,fBinsEta);
293 fMCGenZvPtEtaCent->SetBinEdges(3,fBinsCentrality);
294 fMCGenZvPtEtaCent->GetAxis(0)->SetTitle("MC Zv (cm)");
295 fMCGenZvPtEtaCent->GetAxis(1)->SetTitle("MC Pt (GeV/c)");
296 fMCGenZvPtEtaCent->GetAxis(2)->SetTitle("MC Eta");
297 fMCGenZvPtEtaCent->GetAxis(3)->SetTitle("Centrality");
298 fMCGenZvPtEtaCent->Sumw2();
299
300 fMCRecSecZvPtEtaCent = new THnSparseF("fMCRecSecZvPtEtaCent","mcZv:mcPt:mcEta:Centrality",4,binsZvPtEtaCent);
301 fMCRecSecZvPtEtaCent->SetBinEdges(0,fBinsZv);
302 fMCRecSecZvPtEtaCent->SetBinEdges(1,fBinsPt);
303 fMCRecSecZvPtEtaCent->SetBinEdges(2,fBinsEta);
304 fMCRecSecZvPtEtaCent->SetBinEdges(3,fBinsCentrality);
305 fMCRecSecZvPtEtaCent->GetAxis(0)->SetTitle("MC Sec Zv (cm)");
306 fMCRecSecZvPtEtaCent->GetAxis(1)->SetTitle("MC Sec Pt (GeV/c)");
307 fMCRecSecZvPtEtaCent->GetAxis(2)->SetTitle("MC Sec Eta");
308 fMCRecSecZvPtEtaCent->GetAxis(3)->SetTitle("Centrality");
309 fMCRecSecZvPtEtaCent->Sumw2();
310
bc2a9da9
PL
311 fMCRecPrimPhiPtEtaCent = new THnSparseF("fMCRecPrimPhiPtEtaCent","mcPhi:mcPt:mcEta:Centrality",4,binsPhiPtEtaCent);
312 fMCRecPrimPhiPtEtaCent->SetBinEdges(0,fBinsPhi);
313 fMCRecPrimPhiPtEtaCent->SetBinEdges(1,fBinsPt);
314 fMCRecPrimPhiPtEtaCent->SetBinEdges(2,fBinsEta);
315 fMCRecPrimPhiPtEtaCent->SetBinEdges(3,fBinsCentrality);
8a4ab847 316 fMCRecPrimPhiPtEtaCent->GetAxis(0)->SetTitle("MC Phi");
bc2a9da9
PL
317 fMCRecPrimPhiPtEtaCent->GetAxis(1)->SetTitle("MC Pt (GeV/c)");
318 fMCRecPrimPhiPtEtaCent->GetAxis(2)->SetTitle("MC Eta");
319 fMCRecPrimPhiPtEtaCent->GetAxis(3)->SetTitle("Centrality");
320 fMCRecPrimPhiPtEtaCent->Sumw2();
321
322 fMCGenPhiPtEtaCent = new THnSparseF("fMCGenPhiPtEtaCent","mcPhi:mcPt:mcEta:Centrality",4,binsPhiPtEtaCent);
323 fMCGenPhiPtEtaCent->SetBinEdges(0,fBinsPhi);
324 fMCGenPhiPtEtaCent->SetBinEdges(1,fBinsPt);
325 fMCGenPhiPtEtaCent->SetBinEdges(2,fBinsEta);
326 fMCGenPhiPtEtaCent->SetBinEdges(3,fBinsCentrality);
8a4ab847 327 fMCGenPhiPtEtaCent->GetAxis(0)->SetTitle("MC Phi");
bc2a9da9
PL
328 fMCGenPhiPtEtaCent->GetAxis(1)->SetTitle("MC Pt (GeV/c)");
329 fMCGenPhiPtEtaCent->GetAxis(2)->SetTitle("MC Eta");
330 fMCGenPhiPtEtaCent->GetAxis(3)->SetTitle("Centrality");
331 fMCGenPhiPtEtaCent->Sumw2();
332
333 fMCRecSecPhiPtEtaCent = new THnSparseF("fMCRecSecPhiPtEtaCent","mcPhi:mcPt:mcEta:Centrality",4,binsPhiPtEtaCent);
334 fMCRecSecPhiPtEtaCent->SetBinEdges(0,fBinsPhi);
335 fMCRecSecPhiPtEtaCent->SetBinEdges(1,fBinsPt);
336 fMCRecSecPhiPtEtaCent->SetBinEdges(2,fBinsEta);
337 fMCRecSecPhiPtEtaCent->SetBinEdges(3,fBinsCentrality);
8a4ab847 338 fMCRecSecPhiPtEtaCent->GetAxis(0)->SetTitle("MC Sec Phi");
bc2a9da9
PL
339 fMCRecSecPhiPtEtaCent->GetAxis(1)->SetTitle("MC Sec Pt (GeV/c)");
340 fMCRecSecPhiPtEtaCent->GetAxis(2)->SetTitle("MC Sec Eta");
341 fMCRecSecPhiPtEtaCent->GetAxis(3)->SetTitle("Centrality");
342 fMCRecSecPhiPtEtaCent->Sumw2();
343
ea3cfeda
PL
344 fPt = new TH1F("fPt","fPt",2000,0,200);
345 fPt->GetXaxis()->SetTitle("p_{T} (GeV/c)");
346 fPt->GetYaxis()->SetTitle("dN/dp_{T}");
347 fPt->Sumw2();
348
349 fMCPt = new TH1F("fMCPt","fMCPt",2000,0,200);
350 fMCPt->GetXaxis()->SetTitle("MC p_{T} (GeV/c)");
351 fMCPt->GetYaxis()->SetTitle("dN/dp_{T}");
352 fMCPt->Sumw2();
353
354 fEventStatistics = new TH1F("fEventStatistics","fEventStatistics",10,0,10);
355 fEventStatistics->GetYaxis()->SetTitle("number of events");
356 fEventStatistics->SetBit(TH1::kCanRebin);
357
358 fEventStatisticsCentrality = new TH1F("fEventStatisticsCentrality","fEventStatisticsCentrality",fCentralityNbins-1, fBinsCentrality);
359 fEventStatisticsCentrality->GetYaxis()->SetTitle("number of events");
360
361 fMCEventStatisticsCentrality = new TH1F("fMCEventStatisticsCentrality","fMCEventStatisticsCentrality",fCentralityNbins-1, fBinsCentrality);
362 fMCEventStatisticsCentrality->GetYaxis()->SetTitle("number of MC events");
363
364 fAllEventStatisticsCentrality = new TH1F("fAllEventStatisticsCentrality","fAllEventStatisticsCentrality",fCentralityNbins-1, fBinsCentrality);
365 fAllEventStatisticsCentrality->GetYaxis()->SetTitle("number of events");
366
367 fEventStatisticsCentralityTrigger = new TH2F("fEventStatisticsCentralityTrigger","fEventStatisticsCentralityTrigger;centrality;trigger",100,0,100,3,0,3);
368 fEventStatisticsCentralityTrigger->Sumw2();
369
370 fZvMultCent = new THnSparseF("fZvMultCent","Zv:mult:Centrality",3,binsZvMultCent);
371 fZvMultCent->SetBinEdges(0,fBinsZv);
372 fZvMultCent->SetBinEdges(1,fBinsMult);
373 fZvMultCent->SetBinEdges(2,fBinsCentrality);
374 fZvMultCent->GetAxis(0)->SetTitle("Zv (cm)");
375 fZvMultCent->GetAxis(1)->SetTitle("N_{acc}");
376 fZvMultCent->GetAxis(2)->SetTitle("Centrality");
377 fZvMultCent->Sumw2();
378
379 fTriggerStatistics = new TH1F("fTriggerStatistics","fTriggerStatistics",10,0,10);
380 fTriggerStatistics->GetYaxis()->SetTitle("number of events");
381
382 fMCTrackPdgCode = new TH1F("fMCTrackPdgCode","fMCTrackPdgCode",100,0,10);
383 fMCTrackPdgCode->GetYaxis()->SetTitle("number of tracks");
384 fMCTrackPdgCode->SetBit(TH1::kCanRebin);
385
386 fMCTrackStatusCode = new TH1F("fMCTrackStatusCode","fMCTrackStatusCode",100,0,10);
387 fMCTrackStatusCode->GetYaxis()->SetTitle("number of tracks");
388 fMCTrackStatusCode->SetBit(TH1::kCanRebin);
389
390 fCharge = new TH1F("fCharge","fCharge",30, -5, 5);
391 fCharge->GetXaxis()->SetTitle("Charge");
392 fCharge->GetYaxis()->SetTitle("number of tracks");
393
394 fMCCharge = new TH1F("fMCCharge","fMCCharge",30, -5, 5);
395 fMCCharge->GetXaxis()->SetTitle("MC Charge");
396 fMCCharge->GetYaxis()->SetTitle("number of tracks");
397
398 fMCPdgPt = new TH2F("fMCPdgPt","fMCPdgPt",fPtNbins-1, fBinsPt, 100,0,100);
399 fMCPdgPt->GetYaxis()->SetTitle("particle");
400 fMCPdgPt->GetXaxis()->SetTitle("Pt (GeV/c)");
401
402 fMCHijingPrim = new TH1F("fMCHijingPrim","fMCHijingPrim",2,0,2);
403 fMCHijingPrim->GetYaxis()->SetTitle("number of particles");
d25bcbe6 404
95f26ffa 405
95f26ffa 406
8a4ab847 407 Int_t binsDCAxyDCAzPtEtaPhi[6] = { 10,10, fPtCheckNbins-1, fEtaCheckNbins-1, 18, fCentralityNbins-1};
9db7eb94 408 Double_t minDCAxyDCAzPtEtaPhi[6] = { -5, -5, 0, -1.5, 0., 0, };
409 Double_t maxDCAxyDCAzPtEtaPhi[6] = { 5., 5., 100, 1.5, 2.*TMath::Pi(), 100};
d25bcbe6 410
ea3cfeda
PL
411 fDCAPtAll = new THnSparseF("fDCAPtAll","fDCAPtAll",6, binsDCAxyDCAzPtEtaPhi, minDCAxyDCAzPtEtaPhi, maxDCAxyDCAzPtEtaPhi);
412 fDCAPtAccepted = new THnSparseF("fDCAPtAccepted","fDCAPtAccepted",6, binsDCAxyDCAzPtEtaPhi, minDCAxyDCAzPtEtaPhi, maxDCAxyDCAzPtEtaPhi);
413 fMCDCAPtSecondary = new THnSparseF("fMCDCAPtSecondary","fMCDCAPtSecondary",6, binsDCAxyDCAzPtEtaPhi, minDCAxyDCAzPtEtaPhi, maxDCAxyDCAzPtEtaPhi);
414 fMCDCAPtPrimary = new THnSparseF("fMCDCAPtPrimary","fMCDCAPtPrimary",6, binsDCAxyDCAzPtEtaPhi, minDCAxyDCAzPtEtaPhi, maxDCAxyDCAzPtEtaPhi);
415
416 fDCAPtAll->SetBinEdges(2, fBinsPtCheck);
417 fDCAPtAccepted->SetBinEdges(2, fBinsPtCheck);
418 fMCDCAPtSecondary->SetBinEdges(2, fBinsPtCheck);
419 fMCDCAPtPrimary->SetBinEdges(2, fBinsPtCheck);
420
421 fDCAPtAll->SetBinEdges(3, fBinsEtaCheck);
422 fDCAPtAccepted->SetBinEdges(3, fBinsEtaCheck);
423 fMCDCAPtSecondary->SetBinEdges(3, fBinsEtaCheck);
424 fMCDCAPtPrimary->SetBinEdges(3, fBinsEtaCheck);
425
426 fDCAPtAll->SetBinEdges(5, fBinsCentrality);
427 fDCAPtAccepted->SetBinEdges(5, fBinsCentrality);
428 fMCDCAPtSecondary->SetBinEdges(5, fBinsCentrality);
429 fMCDCAPtPrimary->SetBinEdges(5, fBinsCentrality);
430
431 fDCAPtAll->Sumw2();
432 fDCAPtAccepted->Sumw2();
433 fMCDCAPtSecondary->Sumw2();
434 fMCDCAPtPrimary->Sumw2();
435
436 fDCAPtAll->GetAxis(0)->SetTitle("DCA_{xy} (cm)");
437 fDCAPtAll->GetAxis(1)->SetTitle("DCA_{z} (cm)");
438 fDCAPtAll->GetAxis(2)->SetTitle("p_{T} (GeV/c)");
439 fDCAPtAll->GetAxis(3)->SetTitle("#eta");
440 fDCAPtAll->GetAxis(4)->SetTitle("#phi");
441 fDCAPtAll->GetAxis(5)->SetTitle("Centrality");
442
443 fDCAPtAccepted->GetAxis(0)->SetTitle("DCA_{xy} (cm)");
444 fDCAPtAccepted->GetAxis(1)->SetTitle("DCA_{z} (cm)");
445 fDCAPtAccepted->GetAxis(2)->SetTitle("p_{T} (GeV/c)");
446 fDCAPtAccepted->GetAxis(3)->SetTitle("#eta");
447 fDCAPtAccepted->GetAxis(4)->SetTitle("#phi");
448 fDCAPtAccepted->GetAxis(5)->SetTitle("Centrality");
449
450 fMCDCAPtSecondary->GetAxis(0)->SetTitle("DCA_{xy} (cm)");
451 fMCDCAPtSecondary->GetAxis(1)->SetTitle("DCA_{z} (cm)");
452 fMCDCAPtSecondary->GetAxis(2)->SetTitle("p_{T} (GeV/c)");
453 fMCDCAPtSecondary->GetAxis(3)->SetTitle("#eta");
454 fMCDCAPtSecondary->GetAxis(4)->SetTitle("#phi");
455 fMCDCAPtSecondary->GetAxis(5)->SetTitle("Centrality");
456
457 fMCDCAPtPrimary->GetAxis(0)->SetTitle("DCA_{xy} (cm)");
458 fMCDCAPtPrimary->GetAxis(1)->SetTitle("DCA_{z} (cm)");
459 fMCDCAPtPrimary->GetAxis(2)->SetTitle("p_{T} (GeV/c)");
460 fMCDCAPtPrimary->GetAxis(3)->SetTitle("#eta");
461 fMCDCAPtPrimary->GetAxis(4)->SetTitle("#phi");
462 fMCDCAPtPrimary->GetAxis(5)->SetTitle("Centrality");
9db7eb94 463
72bb4ceb 464
465 char cFullTempTitle[255];
466 char cTempTitleAxis0All[255];
467 char cTempTitleAxis0Acc[255];
468 // char cTempTitleAxis1[255];
469 char cFullTempName[255];
470 char cTempNameAxis0[255];
471 // char cTempNameAxis1[255];
472 const Int_t iNbinRowsClusters = 21;
473 // Double_t dBinsRowsClusters[iNbinRowsClusters] = {0, 7.95, 15.9, 23.85, 31.8, 39.75, 47.7, 55.65, 63.6, 71.55, 79.5, 87.45, 95.4, 103.35, 111.3, 119.25, 127.2, 135.15, 143.1, 151.05, 159.};
474
475 const Int_t iNbinChi = 51;
ea3cfeda 476 const Int_t iNbinLength = 165;
72bb4ceb 477 // Double_t dBinsChi[iNbinChi] = {0, 0.2, 0.4, 0.6, 0.8, 1, 1.2, 1.4, 1.6, 1.8, 2, 2.2, 2.4, 2.6, 2.8, 3, 3.2, 3.4, 3.6, 3.8, 4, 4.2, 4.4, 4.6, 4.8, 5, 5.2, 5.4, 5.6, 5.8, 6, 6.2, 6.4, 6.6, 6.8, 7, 7.2, 7.4, 7.6, 7.8, 8, 8.2, 8.4, 8.6, 8.8, 9, 9.2, 9.4, 9.6, 9.8,10.};
478
479 Int_t iNbin = 0;
480 // Double_t *dBins = 0x0;
481 Double_t dBinMin = 0;
482 Double_t dBinMax = 0;
483
484 for(Int_t iCheckQuant = 0; iCheckQuant < cqMax; iCheckQuant++)
485 {
486 // iCheckQuant: 0 = CrossedRows, 1 = Nclusters, 2 = Chi^2/clusterTPC
487 if(iCheckQuant == cqCrossedRows)
488 {
bc2a9da9
PL
489 snprintf(cTempTitleAxis0All,255, "NcrossedRows before Cut");
490 snprintf(cTempTitleAxis0Acc,255, "NcrossedRows after Cut");
491 snprintf(cTempNameAxis0,255, "CrossedRows");
72bb4ceb 492 iNbin = iNbinRowsClusters;
493 dBinMin = 0;
494 dBinMax = 159.;
495 }
496 else if(iCheckQuant == cqNcluster)
497 {
bc2a9da9
PL
498 snprintf(cTempTitleAxis0All,255, "Nclusters before Cut");
499 snprintf(cTempTitleAxis0Acc,255, "Nclusters after Cut");
500 snprintf(cTempNameAxis0,255, "Clusters");
72bb4ceb 501 iNbin = iNbinRowsClusters;
502 dBinMin = 0;
503 dBinMax = 159.;
504 }
505 else if(iCheckQuant == cqChi)
506 {
bc2a9da9
PL
507 snprintf(cTempTitleAxis0All,255, "#Chi^{2}/cluster before Cut");
508 snprintf(cTempTitleAxis0Acc,255, "#Chi^{2}/cluster after Cut");
509 snprintf(cTempNameAxis0,255, "Chi");
72bb4ceb 510 iNbin = iNbinChi;
511 dBinMin = 0;
512 dBinMax = 10.;
513 }
ea3cfeda
PL
514 else if(iCheckQuant == cqLength)
515 {
bc2a9da9
PL
516 snprintf(cTempTitleAxis0All,255, "Length in TPC before Cut (cm)");
517 snprintf(cTempTitleAxis0Acc,255, "Length in TPC after Cut (cm)");
518 snprintf(cTempNameAxis0,255, "Length");
ea3cfeda
PL
519 iNbin = iNbinLength;
520 dBinMin = 0;
521 dBinMax = 165.;
522 }
72bb4ceb 523
524 Int_t binsCheckPtEtaPhi[5] = { iNbin, fPtCheckNbins-1, fEtaCheckNbins-1, 18, fCentralityNbins-1};
a0036e80 525 // Int_t binsCheckPtEtaPhi[5] = { iNbin, fPtNbins-1, fEtaCheckNbins-1, 18, fCentralityNbins-1};
72bb4ceb 526 Double_t minCheckPtEtaPhi[5] = { dBinMin, 0, -1.5, 0., 0, };
527 Double_t maxCheckPtEtaPhi[5] = { dBinMax, 100, 1.5, 2.*TMath::Pi(), 100};
528
bc2a9da9
PL
529 snprintf(cFullTempName, 255, "f%sPtEtaPhiAll",cTempNameAxis0);
530 snprintf(cFullTempTitle, 255,"%s;%s;p_{T} (GeV/c);#eta;#phi;Centrality", cFullTempName, cTempTitleAxis0All);
ea3cfeda
PL
531 fCrossCheckAll[iCheckQuant] = new THnF(cFullTempName, cFullTempTitle, 5, binsCheckPtEtaPhi, minCheckPtEtaPhi, maxCheckPtEtaPhi);
532 fCrossCheckAll[iCheckQuant]->SetBinEdges(1, fBinsPtCheck);
533 fCrossCheckAll[iCheckQuant]->SetBinEdges(2, fBinsEtaCheck);
534 fCrossCheckAll[iCheckQuant]->Sumw2();
72bb4ceb 535
bc2a9da9
PL
536 snprintf(cFullTempName, 255, "f%sPtEtaPhiAcc",cTempNameAxis0);
537 snprintf(cFullTempTitle, 255,"%s;%s;p_{T} (GeV/c);#eta;#phi;Centrality", cFullTempName, cTempTitleAxis0Acc);
ea3cfeda
PL
538 fCrossCheckAcc[iCheckQuant] = new THnF(cFullTempName, cFullTempTitle, 5, binsCheckPtEtaPhi, minCheckPtEtaPhi, maxCheckPtEtaPhi);
539 fCrossCheckAcc[iCheckQuant]->SetBinEdges(1, fBinsPtCheck);
540 fCrossCheckAcc[iCheckQuant]->SetBinEdges(2, fBinsEtaCheck);
541 fCrossCheckAcc[iCheckQuant]->Sumw2();
72bb4ceb 542 } // end iCheckQuant
fad9b70b 543
ea3cfeda
PL
544 fCutPercClusters = new TH1F("fCutPercClusters","fCutPercClusters;NclustersTPC;counts",160,0,160);
545 fCutPercClusters->Sumw2();
546 fCutPercCrossed = new TH1F("fCutPercCrossed","fCutPercCrossed;NcrossedRowsTPC;counts",160,0,160);
547 fCutPercCrossed->Sumw2();
548
549 fCrossCheckRowsLength = new TH2F("fCrossCheckRowsLength","fCrossCheckRowsLength;Length in TPC;NcrossedRows",170,0,170,170,0,170);
550 fCrossCheckRowsLength->Sumw2();
551
552 fCrossCheckClusterLength = new TH2F("fCrossCheckClusterLength","fCrossCheckClusterLength;Length in TPC;NClusters",170,0,170,170,0,170);
553 fCrossCheckClusterLength->Sumw2();
554
555 fCrossCheckRowsLengthAcc = new TH2F("fCrossCheckRowsLengthAcc","fCrossCheckRowsLengthAcc;Length in TPC;NcrossedRows",170,0,170,170,0,170);
556 fCrossCheckRowsLengthAcc->Sumw2();
557
558 fCrossCheckClusterLengthAcc = new TH2F("fCrossCheckClusterLengthAcc","fCrossCheckClusterLengthAcc;Length in TPC;NClusters",170,0,170,170,0,170);
559 fCrossCheckClusterLengthAcc->Sumw2();
560
8a4ab847 561 fCutSettings = new TH1F("fCutSettings","fCutSettings",100,0,10);
562 fCutSettings->GetYaxis()->SetTitle("cut value");
563 fCutSettings->SetBit(TH1::kCanRebin);
ea3cfeda 564
d25bcbe6 565 // Add Histos, Profiles etc to List
ea3cfeda 566 fOutputList->Add(fZvPtEtaCent);
bc2a9da9 567 fOutputList->Add(fPhiPtEtaCent);
a0036e80 568 fOutputList->Add(fPtResptCent);
ea3cfeda
PL
569 fOutputList->Add(fPt);
570 fOutputList->Add(fMCRecPrimZvPtEtaCent);
571 fOutputList->Add(fMCGenZvPtEtaCent);
572 fOutputList->Add(fMCRecSecZvPtEtaCent);
bc2a9da9
PL
573 fOutputList->Add(fMCRecPrimPhiPtEtaCent);
574 fOutputList->Add(fMCGenPhiPtEtaCent);
575 fOutputList->Add(fMCRecSecPhiPtEtaCent);
ea3cfeda
PL
576 fOutputList->Add(fMCPt);
577 fOutputList->Add(fEventStatistics);
578 fOutputList->Add(fEventStatisticsCentrality);
579 fOutputList->Add(fMCEventStatisticsCentrality);
580 fOutputList->Add(fAllEventStatisticsCentrality);
581 fOutputList->Add(fEventStatisticsCentralityTrigger);
582 fOutputList->Add(fZvMultCent);
583 fOutputList->Add(fTriggerStatistics);
584 fOutputList->Add(fMCTrackPdgCode);
585 fOutputList->Add(fMCTrackStatusCode);
586 fOutputList->Add(fCharge);
587 fOutputList->Add(fMCCharge);
588 fOutputList->Add(fMCPdgPt);
589 fOutputList->Add(fMCHijingPrim);
590 fOutputList->Add(fDCAPtAll);
591 fOutputList->Add(fDCAPtAccepted);
592 fOutputList->Add(fMCDCAPtSecondary);
593 fOutputList->Add(fMCDCAPtPrimary);
72bb4ceb 594 for(Int_t i = 0; i < cqMax; i++)
595 {
ea3cfeda
PL
596 fOutputList->Add(fCrossCheckAll[i]);
597 fOutputList->Add(fCrossCheckAcc[i]);
72bb4ceb 598 }
ea3cfeda
PL
599 fOutputList->Add(fCutPercClusters);
600 fOutputList->Add(fCutPercCrossed);
601 fOutputList->Add(fCrossCheckRowsLength);
602 fOutputList->Add(fCrossCheckClusterLength);
603 fOutputList->Add(fCrossCheckRowsLengthAcc);
604 fOutputList->Add(fCrossCheckClusterLengthAcc);
8a4ab847 605 fOutputList->Add(fCutSettings);
606
607 StoreCutSettingsToHistogram();
d25bcbe6 608
609 PostData(1, fOutputList);
610}
611
612void AlidNdPtAnalysisPbPbAOD::UserExec(Option_t *option)
613{
f856053f 614
d25bcbe6 615 // Main Loop
616 // called for each event
ea3cfeda 617 fEventStatistics->Fill("all events",1);
d25bcbe6 618
619 // set ZERO pointers:
620 AliInputEventHandler *inputHandler = NULL;
621 AliAODTrack *track = NULL;
622 AliAODMCParticle *mcPart = NULL;
623 AliAODMCHeader *mcHdr = NULL;
624 AliGenHijingEventHeader *genHijingHeader = NULL;
b3341d37 625 //AliGenPythiaEventHeader *genPythiaHeader = NULL;
d25bcbe6 626
627 Bool_t bIsEventSelectedMB = kFALSE;
628 Bool_t bIsEventSelectedSemi = kFALSE;
629 Bool_t bIsEventSelectedCentral = kFALSE;
630 Bool_t bIsEventSelected = kFALSE;
631 Bool_t bIsPrimary = kFALSE;
632 Bool_t bIsHijingParticle = kFALSE;
9db7eb94 633 Bool_t bMotherIsHijingParticle = kFALSE;
b3341d37 634 //Bool_t bIsPythiaParticle = kFALSE;
d25bcbe6 635 Bool_t bEventHasATrack = kFALSE;
636 Bool_t bEventHasATrackInRange = kFALSE;
637 Int_t nTriggerFired = 0;
638
639
640 Double_t dMCTrackZvPtEtaCent[4] = {0};
641 Double_t dTrackZvPtEtaCent[4] = {0};
642
bc2a9da9
PL
643 Double_t dMCTrackPhiPtEtaCent[4] = {0};
644 Double_t dTrackPhiPtEtaCent[4] = {0};
645
b3341d37 646 Double_t dDCA[2] = {0};
647
d25bcbe6 648 Double_t dMCEventZv = -100;
649 Double_t dEventZv = -100;
650 Int_t iAcceptedMultiplicity = 0;
651
ea3cfeda 652 fIsMonteCarlo = kFALSE;
d25bcbe6 653
654 AliAODEvent *eventAOD = 0x0;
655 eventAOD = dynamic_cast<AliAODEvent*>( InputEvent() );
656 if (!eventAOD) {
657 AliWarning("ERROR: eventAOD not available \n");
658 return;
659 }
660
661 // check, which trigger has been fired
662 inputHandler = (AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
663 bIsEventSelectedMB = ( inputHandler->IsEventSelected() & AliVEvent::kMB);
664 bIsEventSelectedSemi = ( inputHandler->IsEventSelected() & AliVEvent::kSemiCentral);
665 bIsEventSelectedCentral = ( inputHandler->IsEventSelected() & AliVEvent::kCentral);
666
ea3cfeda
PL
667 if(bIsEventSelectedMB || bIsEventSelectedSemi || bIsEventSelectedCentral) fTriggerStatistics->Fill("all triggered events",1);
668 if(bIsEventSelectedMB) { fTriggerStatistics->Fill("MB trigger",1); nTriggerFired++; }
669 if(bIsEventSelectedSemi) { fTriggerStatistics->Fill("SemiCentral trigger",1); nTriggerFired++; }
670 if(bIsEventSelectedCentral) { fTriggerStatistics->Fill("Central trigger",1); nTriggerFired++; }
671 if(nTriggerFired == 0) { fTriggerStatistics->Fill("No trigger",1); }
d25bcbe6 672
673 bIsEventSelected = ( inputHandler->IsEventSelected() & GetCollisionCandidates() );
674
675 // only take tracks of events, which are triggered
676 if(nTriggerFired == 0) { return; }
677
678 // if( !bIsEventSelected || nTriggerFired>1 ) return;
679
ea3cfeda 680 // fEventStatistics->Fill("events with only coll. cand.", 1);
d25bcbe6 681
682
683
684 // check if there is a stack, if yes, then do MC loop
685 TList *list = eventAOD->GetList();
686 TClonesArray *stack = 0x0;
687 stack = (TClonesArray*)list->FindObject(AliAODMCParticle::StdBranchName());
688
689 if( stack )
690 {
ea3cfeda 691 fIsMonteCarlo = kTRUE;
d25bcbe6 692
693 mcHdr = (AliAODMCHeader*)list->FindObject(AliAODMCHeader::StdBranchName());
694
695 genHijingHeader = GetHijingEventHeader(mcHdr);
696 // genPythiaHeader = GetPythiaEventHeader(mcHdr);
697
698 if(!genHijingHeader) { return; }
699
700 // if(!genPythiaHeader) { return; }
701
702 dMCEventZv = mcHdr->GetVtxZ();
703 dMCTrackZvPtEtaCent[0] = dMCEventZv;
ea3cfeda 704 fEventStatistics->Fill("MC all events",1);
d25bcbe6 705 }
706
707 AliCentrality* aCentrality = eventAOD->GetCentrality();
708 Double_t dCentrality = aCentrality->GetCentralityPercentile("V0M");
709
710 if( dCentrality < 0 ) return;
ea3cfeda 711 fEventStatistics->Fill("after centrality selection",1);
d25bcbe6 712
713
714
715 // start with MC truth analysis
ea3cfeda 716 if(fIsMonteCarlo)
d25bcbe6 717 {
718
ea3cfeda 719 if( dMCEventZv > GetCutMaxZVertex() ) { return; }
d25bcbe6 720
721 dMCTrackZvPtEtaCent[0] = dMCEventZv;
722
ea3cfeda 723 fEventStatistics->Fill("MC afterZv cut",1);
d25bcbe6 724
725 for(Int_t iMCtrack = 0; iMCtrack < stack->GetEntriesFast(); iMCtrack++)
726 {
727 mcPart =(AliAODMCParticle*)stack->At(iMCtrack);
728
729 // check for charge
730 if( !(IsMCTrackAccepted(mcPart)) ) continue;
731
732 if(!IsHijingParticle(mcPart, genHijingHeader)) { continue; }
733
734 if(mcPart->IsPhysicalPrimary() )
735 {
ea3cfeda 736 fMCHijingPrim->Fill("IsPhysicalPrimary",1);
d25bcbe6 737 }
738 else
739 {
ea3cfeda 740 fMCHijingPrim->Fill("NOT a primary",1);
d25bcbe6 741 continue;
742 }
743
d0483ba3 744
d25bcbe6 745 //
746 // ======================== fill histograms ========================
747 dMCTrackZvPtEtaCent[1] = mcPart->Pt();
748 dMCTrackZvPtEtaCent[2] = mcPart->Eta();
749 dMCTrackZvPtEtaCent[3] = dCentrality;
bc2a9da9
PL
750 fMCGenZvPtEtaCent->Fill(dMCTrackZvPtEtaCent);
751
752 dMCTrackPhiPtEtaCent[0] = mcPart->Phi();
753 dMCTrackPhiPtEtaCent[1] = mcPart->Pt();
754 dMCTrackPhiPtEtaCent[2] = mcPart->Eta();
755 dMCTrackPhiPtEtaCent[3] = dCentrality;
756 fMCGenPhiPtEtaCent->Fill(dMCTrackPhiPtEtaCent);
d25bcbe6 757
758 bEventHasATrack = kTRUE;
759
d25bcbe6 760
ea3cfeda
PL
761 if( (dMCTrackZvPtEtaCent[1] > GetCutPtMin() ) &&
762 (dMCTrackZvPtEtaCent[1] < GetCutPtMax() ) &&
763 (dMCTrackZvPtEtaCent[2] > GetCutEtaMin() ) &&
764 (dMCTrackZvPtEtaCent[2] < GetCutEtaMax() ) )
d25bcbe6 765 {
ea3cfeda
PL
766 fMCPt->Fill(mcPart->Pt());
767 fMCCharge->Fill(mcPart->Charge()/3.);
d25bcbe6 768 bEventHasATrackInRange = kTRUE;
769 }
770
771 }
772 } // isMonteCarlo
ea3cfeda
PL
773
774 if(bEventHasATrack) { fEventStatistics->Fill("MC events with tracks",1); }
9db7eb94 775 if(bEventHasATrackInRange)
776 {
ea3cfeda
PL
777 fEventStatistics->Fill("MC events with tracks in range",1);
778 fMCEventStatisticsCentrality->Fill(dCentrality);
9db7eb94 779 }
d25bcbe6 780 bEventHasATrack = kFALSE;
781 bEventHasATrackInRange = kFALSE;
782
783
784
785 // Loop over recontructed tracks
786
787 dEventZv = eventAOD->GetPrimaryVertex()->GetZ();
ea3cfeda 788 if( TMath::Abs(dEventZv) > GetCutMaxZVertex() ) return;
d25bcbe6 789
ea3cfeda 790 fAllEventStatisticsCentrality->Fill(dCentrality/*, nTriggerFired*/);
d25bcbe6 791
ea3cfeda 792 fEventStatistics->Fill("after Zv cut",1);
d25bcbe6 793
794 dTrackZvPtEtaCent[0] = dEventZv;
795
ea3cfeda
PL
796 if(AreRelativeCutsEnabled())
797 {
798 if(!SetRelativeCuts(eventAOD)) return;
799 }
800
d25bcbe6 801 for(Int_t itrack = 0; itrack < eventAOD->GetNumberOfTracks(); itrack++)
802 {
803 track = eventAOD->GetTrack(itrack);
804 if(!track) continue;
805
806 mcPart = NULL;
807 dMCTrackZvPtEtaCent[1] = 0;
808 dMCTrackZvPtEtaCent[2] = 0;
809 dMCTrackZvPtEtaCent[3] = 0;
810
bc2a9da9
PL
811 dMCTrackPhiPtEtaCent[0] = 0;
812 dMCTrackPhiPtEtaCent[1] = 0;
813 dMCTrackPhiPtEtaCent[2] = 0;
814 dMCTrackPhiPtEtaCent[3] = 0;
815
d25bcbe6 816 bIsPrimary = kFALSE;
817
b3341d37 818 GetDCA(track, eventAOD, dDCA);
819
9db7eb94 820 Double_t dDCAxyDCAzPt[5] = { dDCA[0], dDCA[1], track->Pt(), track->Eta(), track->Phi() };
d0483ba3 821
ea3cfeda 822 fDCAPtAll->Fill(dDCAxyDCAzPt);
d0483ba3 823
ea3cfeda 824 if( !(IsTrackAccepted(track, dCentrality, eventAOD->GetMagneticField())) ) continue;
d25bcbe6 825
5747f2c6 826 dTrackZvPtEtaCent[1] = track->Pt();
827 dTrackZvPtEtaCent[2] = track->Eta();
828 dTrackZvPtEtaCent[3] = dCentrality;
829
bc2a9da9
PL
830 dTrackPhiPtEtaCent[0] = track->Phi();
831 dTrackPhiPtEtaCent[1] = track->Pt();
832 dTrackPhiPtEtaCent[2] = track->Eta();
833 dTrackPhiPtEtaCent[3] = dCentrality;
834
ea3cfeda 835 if( fIsMonteCarlo )
d25bcbe6 836 {
837 mcPart = (AliAODMCParticle*)stack->At(TMath::Abs(track->GetLabel()));
838 if( !mcPart ) { continue; }
839
840 // check for charge
ea3cfeda 841 // if( !(IsMCTrackAccepted(mcPart)) ) { continue; }
d25bcbe6 842
843 bIsHijingParticle = IsHijingParticle(mcPart, genHijingHeader);
844 // bIsPythiaParticle = IsPythiaParticle(mcPart, genPythiaHeader);
845
846 // if(!bIsHijingParticle) continue; // only take real tracks, not injected ones
847
848 bIsPrimary = mcPart->IsPhysicalPrimary();
849
850 dMCTrackZvPtEtaCent[1] = mcPart->Pt();
851 dMCTrackZvPtEtaCent[2] = mcPart->Eta();
852 dMCTrackZvPtEtaCent[3] = dCentrality;
853
bc2a9da9 854 dMCTrackPhiPtEtaCent[0] = mcPart->Phi();
a0036e80 855 dMCTrackPhiPtEtaCent[1] = mcPart->Pt();
856 dMCTrackPhiPtEtaCent[2] = mcPart->Eta();
857 dMCTrackPhiPtEtaCent[3] = dCentrality;
bc2a9da9 858
d25bcbe6 859 if(bIsPrimary && bIsHijingParticle)
860 {
ea3cfeda 861 fMCRecPrimZvPtEtaCent->Fill(dMCTrackZvPtEtaCent);
bc2a9da9 862 fMCRecPrimPhiPtEtaCent->Fill(dMCTrackPhiPtEtaCent);
ea3cfeda
PL
863 if( (TMath::Abs(mcPart->Eta()) < 0.8) && (dCentrality<5.) ) { fMCPdgPt->Fill(mcPart->Pt(), Form("%s",GetParticleName(mcPart->GetPdgCode())), 1); }
864 fMCDCAPtPrimary->Fill(dDCAxyDCAzPt);
d25bcbe6 865 }
866
867 if(!bIsPrimary /*&& !bIsHijingParticle*/)
868 {
869 Int_t indexMoth = mcPart->GetMother();
870 if(indexMoth >= 0)
871 {
872 AliAODMCParticle* moth = (AliAODMCParticle*)stack->At(indexMoth);
9db7eb94 873 bMotherIsHijingParticle = IsHijingParticle(moth, genHijingHeader);
d25bcbe6 874
875 if(bMotherIsHijingParticle) // only store secondaries, which come from a not embedded signal!
876 {
ea3cfeda
PL
877 fMCTrackStatusCode->Fill(Form("%d",mcPart->GetStatus()), 1);
878
d25bcbe6 879
ea3cfeda 880 fMCRecSecZvPtEtaCent->Fill(dMCTrackZvPtEtaCent);
bc2a9da9 881 fMCRecSecPhiPtEtaCent->Fill(dMCTrackPhiPtEtaCent);
ea3cfeda
PL
882 fMCDCAPtSecondary->Fill(dDCAxyDCAzPt);
883 fMCTrackPdgCode->Fill(Form("%s_H%i_H%i",GetParticleName(moth->GetPdgCode()),bMotherIsHijingParticle, bIsHijingParticle), 1);
d25bcbe6 884 // delete moth;
885 }
886 }
887 }
888 } // end isMonteCarlo
889
890 // ======================== fill histograms ========================
891
9db7eb94 892 // only keep prim and sec from not embedded signal
893 Bool_t bKeepMCTrack = kFALSE;
ea3cfeda 894 if(fIsMonteCarlo)
9db7eb94 895 {
896 if( (bIsHijingParticle && bIsPrimary) ^ (bMotherIsHijingParticle && !bIsPrimary) )
d0483ba3 897 {
9db7eb94 898 bKeepMCTrack = kTRUE;
899 }
900 else
901 {
902 continue;
d0483ba3 903 }
9db7eb94 904 }
905
906 bEventHasATrack = kTRUE;
907
ea3cfeda 908 fZvPtEtaCent->Fill(dTrackZvPtEtaCent);
bc2a9da9
PL
909 fPhiPtEtaCent->Fill(dTrackPhiPtEtaCent);
910
ea3cfeda 911 fDCAPtAccepted->Fill(dDCAxyDCAzPt);
9db7eb94 912
ea3cfeda
PL
913 if( (dTrackZvPtEtaCent[1] > GetCutPtMin()) &&
914 (dTrackZvPtEtaCent[1] < GetCutPtMax()) &&
915 (dTrackZvPtEtaCent[2] > GetCutEtaMin()) &&
916 (dTrackZvPtEtaCent[2] < GetCutEtaMax()) )
9db7eb94 917 {
918 iAcceptedMultiplicity++;
919 bEventHasATrackInRange = kTRUE;
ea3cfeda
PL
920 fPt->Fill(track->Pt());
921 fCharge->Fill(track->Charge());
9db7eb94 922 }
d25bcbe6 923 } // end track loop
924
ea3cfeda 925 if(bEventHasATrack) { fEventStatistics->Fill("events with tracks",1); bEventHasATrack = kFALSE; }
d25bcbe6 926
927 if(bEventHasATrackInRange)
928 {
ea3cfeda
PL
929 fEventStatistics->Fill("events with tracks in range",1);
930 fEventStatisticsCentrality->Fill(dCentrality);
d25bcbe6 931 bEventHasATrackInRange = kFALSE;
fad9b70b 932
ea3cfeda
PL
933 if(bIsEventSelectedMB) fEventStatisticsCentralityTrigger->Fill(dCentrality, 0);
934 if(bIsEventSelectedSemi) fEventStatisticsCentralityTrigger->Fill(dCentrality, 1);
935 if(bIsEventSelectedCentral) fEventStatisticsCentralityTrigger->Fill(dCentrality, 2);
d25bcbe6 936 }
937
938 Double_t dEventZvMultCent[3] = {dEventZv, iAcceptedMultiplicity, dCentrality};
ea3cfeda 939 fZvMultCent->Fill(dEventZvMultCent);
d25bcbe6 940
d25bcbe6 941 PostData(1, fOutputList);
942
8a4ab847 943 // delete pointers:
944
d25bcbe6 945}
946
ea3cfeda
PL
947Bool_t AlidNdPtAnalysisPbPbAOD::SetRelativeCuts(AliAODEvent *event)
948{
949 if(!event) return kFALSE;
950
951 AliAODTrack *tr = 0x0;
952 TH1F *hCluster = new TH1F("hCluster","hCluster",160,0,160);
953 TH1F *hCrossed = new TH1F("hCrossed","hCrossed",160,0,160);
954
955 for(Int_t itrack = 0; itrack < event->GetNumberOfTracks(); itrack++)
956 {
957 tr = event->GetTrack(itrack);
958 if(!tr) continue;
959
960 // do some selection already
961 //if(!(tr->TestFilterBit(AliAODTrack::kTrkGlobal)) ) { continue; }
962
963 Double_t dNClustersTPC = tr->GetTPCNcls();
964 Double_t dCrossedRowsTPC = tr->GetTPCClusterInfo(2,1);
965
966 hCluster->Fill(dNClustersTPC);
967 hCrossed->Fill(dCrossedRowsTPC);
968 }
969
970 // loop trough histogram to check, where percentage is reach
971 Double_t dTotIntCluster = hCluster->Integral();
972 Double_t dTotIntCrossed = hCrossed->Integral();
973 Float_t dIntCluster = 0;
974 Float_t dIntCrossed = 0;
975
976 if(dTotIntCluster)
977 {
978 for(Int_t i = 0; i < hCluster->GetNbinsX(); i++)
979 {
980 if(hCluster->GetBinCenter(i) < 0) continue;
981 dIntCluster += hCluster->GetBinContent(i);
982 if(dIntCluster/dTotIntCluster > (1-GetCutPercMinNClustersTPC()))
983 {
984 SetCutMinNClustersTPC(hCluster->GetBinCenter(i));
985 fCutPercClusters->Fill(hCluster->GetBinCenter(i));
986 break;
987 }
988 }
989 }
990
991 if(dTotIntCrossed)
992 {
993 for(Int_t i = 0; i < hCrossed->GetNbinsX(); i++)
994 {
995 if(hCrossed->GetBinCenter(i) < 0) continue;
996 dIntCrossed += hCrossed->GetBinContent(i);
997 if(dIntCrossed/dTotIntCrossed > (1-GetCutPercMinNCrossedRowsTPC()))
998 {
999 SetCutMinNClustersTPC(hCrossed->GetBinCenter(i));
1000 fCutPercCrossed->Fill(hCrossed->GetBinCenter(i));
1001 break;
1002 }
1003 }
1004 }
1005
1006 delete hCrossed;
1007 delete hCluster;
1008 return kTRUE;
1009
1010}
1011
1012Bool_t AlidNdPtAnalysisPbPbAOD::IsTrackAccepted(AliAODTrack *tr, Double_t dCentrality, Double_t bMagZ)
d25bcbe6 1013{
1014 if(!tr) return kFALSE;
1015
1016 if(tr->Charge()==0) { return kFALSE; }
1017
ea3cfeda
PL
1018 //
1019 // as done in AliAnalysisTaskFragmentationFunction
1020 //
1021
1022 Short_t sign = tr->Charge();
1023 Double_t xyz[50];
1024 Double_t pxpypz[50];
a0036e80 1025 Double_t cv[21];
ea3cfeda 1026
a0036e80 1027 for(Int_t i = 0; i < 21; i++) cv[i] = 0;
ea3cfeda
PL
1028 for(Int_t i = 0; i < 50; i++) xyz[i] = 0;
1029 for(Int_t i = 0; i < 50; i++) pxpypz[i] = 0;
a0036e80 1030
ea3cfeda
PL
1031 tr->GetXYZ(xyz);
1032 tr->GetPxPyPz(pxpypz);
a0036e80 1033 tr->GetCovarianceXYZPxPyPz(cv);
ea3cfeda
PL
1034
1035 // similar error occured as this one:
1036 // See https://savannah.cern.ch/bugs/?102721
1037 // which is one of the two 11h re-filtering follow-ups:
1038 // Andrea Dainese now first does the beam pipe
1039 // check and then copies from the vtrack (was the other
1040 // way around) to avoid the crash in the etp::Set()
a0036e80 1041
1042 // if(xyz[0]*xyz[0]+xyz[1]*xyz[1] > 3.*3.) { return kFALSE; }
ea3cfeda 1043
8a4ab847 1044 AliExternalTrackParam par(xyz, pxpypz, cv, sign);
1045// AliExternalTrackParam *par = new AliExternalTrackParam(xyz, pxpypz, cv, sign); // high mem consumption!!!!
ea3cfeda 1046 AliESDtrack dummy;
a0036e80 1047 // Double_t dLength = dummy.GetLengthInActiveZone(par,3,236, -5 ,0,0);
1048 // Double_t dLengthInTPC = GetLengthInTPC(tr, 1.8, 220, bMagZ);
1049
8a4ab847 1050 Double_t dLengthInTPC = dummy.GetLengthInActiveZone(&par,3,236, bMagZ ,0,0);
d25bcbe6 1051 Double_t dNClustersTPC = tr->GetTPCNcls();
1052 Double_t dCrossedRowsTPC = tr->GetTPCClusterInfo(2,1);
ea3cfeda 1053 Double_t dFindableClustersTPC = tr->GetTPCNclsF();
9db7eb94 1054 Double_t dChi2PerClusterTPC = (dNClustersTPC>0)?tr->Chi2perNDF()*(dNClustersTPC-5)/dNClustersTPC:-1.; // see AliDielectronVarManager.h
a0036e80 1055 Double_t dOneOverPt = tr->OneOverPt();
8a4ab847 1056 Double_t dSigmaOneOverPt = TMath::Sqrt(par.GetSigma1Pt2());
9db7eb94 1057
1058 // hAllCrossedRowsTPC->Fill(dCrossedRowsTPC);
d25bcbe6 1059
ea3cfeda
PL
1060 Double_t dCheck[cqMax] = {dCrossedRowsTPC, dNClustersTPC, dChi2PerClusterTPC, dLengthInTPC};// = new Double_t[cqMax];
1061 Double_t dKine[kqMax] = {tr->Pt(), tr->Eta(), tr->Phi()};// = new Double_t[kqMax];
1062 // dKine[0] = tr->Pt();
1063 // dKine[1] = tr->Eta();
1064 // dKine[2] = tr->Phi();
1065 //
1066 // dCheck[0] = dCrossedRowsTPC;
1067 // dCheck[1] = dNClustersTPC;
1068 // dCheck[2] = dChi2PerClusterTPC;
72bb4ceb 1069
1070
1071 FillDebugHisto(dCheck, dKine, dCentrality, kFALSE);
1072
ea3cfeda
PL
1073 // first cut on length
1074
1075 if( DoCutLengthInTPCPtDependent() && ( dLengthInTPC < GetPrefactorLengthInTPCPtDependent()*(130-5*TMath::Abs(1./tr->Pt())) ) ) { return kFALSE; }
95f26ffa 1076
9db7eb94 1077 // filter bit 5
a0036e80 1078 // if(!(tr->TestFilterBit(AliAODTrack::kTrkGlobal)) ) { return kFALSE; }
1079 if(!(tr->TestFilterBit(GetFilterBit())) ) { return kFALSE; }
95f26ffa 1080
9db7eb94 1081 // filter bit 4
1082 // if(!(tr->TestFilterBit(AliAODTrack::kTrkGlobalNoDCA)) ) { return kFALSE; }
1083
1084 // hFilterCrossedRowsTPC->Fill(dCrossedRowsTPC);
95f26ffa 1085
d0483ba3 1086
a0036e80 1087 if(dFindableClustersTPC == 0) {return kFALSE; }
1088 if(dCrossedRowsTPC < GetCutMinNCrossedRowsTPC()) { return kFALSE; }
1089 if( (dCrossedRowsTPC/dFindableClustersTPC) < GetCutMinRatioCrossedRowsOverFindableClustersTPC() ) { return kFALSE; }
1090 if(dNClustersTPC < GetCutMinNClustersTPC()) { return kFALSE; }
72bb4ceb 1091
a0036e80 1092 if (IsITSRefitRequired() && !(tr->GetStatus() & AliVTrack::kITSrefit)) { return kFALSE; } // no ITS refit
1093
1094 // do a relativ cut in Nclusters, both time at 80% of mean
1095 // if(fIsMonteCarlo)
1096 // {
1097 // if(dNClustersTPC < 88) { return kFALSE; }
ea3cfeda 1098 // }
a0036e80 1099 // else
ea3cfeda 1100 // {
a0036e80 1101 // if(dNClustersTPC < 76) { return kFALSE; }
ea3cfeda 1102 // }
ea3cfeda 1103
a0036e80 1104 // fill histogram for pT resolution correction
1105 Double_t dPtResolutionHisto[3] = { dOneOverPt, dSigmaOneOverPt, dCentrality };
1106 fPtResptCent->Fill(dPtResolutionHisto);
ea3cfeda 1107
a0036e80 1108 // fill debug histogram for all accepted tracks
1109 FillDebugHisto(dCheck, dKine, dCentrality, kTRUE);
8a4ab847 1110
1111 // delete pointers
1112
ea3cfeda 1113 return kTRUE;
d25bcbe6 1114}
1115
72bb4ceb 1116Bool_t AlidNdPtAnalysisPbPbAOD::FillDebugHisto(Double_t *dCrossCheckVar, Double_t *dKineVar, Double_t dCentrality, Bool_t bIsAccepted)
1117{
1118 if(bIsAccepted)
1119 {
1120 for(Int_t iCrossCheck = 0; iCrossCheck < cqMax; iCrossCheck++)
1121 {
1122 Double_t dFillIt[5] = {dCrossCheckVar[iCrossCheck], dKineVar[0], dKineVar[1], dKineVar[2], dCentrality};
ea3cfeda 1123 fCrossCheckAcc[iCrossCheck]->Fill(dFillIt);
72bb4ceb 1124 }
ea3cfeda
PL
1125
1126 fCrossCheckRowsLengthAcc->Fill(dCrossCheckVar[cqLength], dCrossCheckVar[cqCrossedRows]);
1127 fCrossCheckClusterLengthAcc->Fill(dCrossCheckVar[cqLength], dCrossCheckVar[cqNcluster]);
72bb4ceb 1128 }
1129 else
1130 {
1131 for(Int_t iCrossCheck = 0; iCrossCheck < cqMax; iCrossCheck++)
1132 {
1133 Double_t dFillIt[5] = {dCrossCheckVar[iCrossCheck], dKineVar[0], dKineVar[1], dKineVar[2], dCentrality};
ea3cfeda 1134 fCrossCheckAll[iCrossCheck]->Fill(dFillIt);
72bb4ceb 1135 }
ea3cfeda
PL
1136
1137 fCrossCheckRowsLength->Fill(dCrossCheckVar[cqLength], dCrossCheckVar[cqCrossedRows]);
1138 fCrossCheckClusterLength->Fill(dCrossCheckVar[cqLength], dCrossCheckVar[cqNcluster]);
72bb4ceb 1139 }
1140
1141 return kTRUE;
1142
1143}
1144
8a4ab847 1145void AlidNdPtAnalysisPbPbAOD::StoreCutSettingsToHistogram()
1146{
1147 //
1148 // this function stores all cut settings to a histograms
1149 //
1150
1151 fCutSettings->Fill("IsMonteCarlo",fIsMonteCarlo);
1152
1153 fCutSettings->Fill("fCutMaxZVertex", fCutMaxZVertex);
1154
1155 // kinematic cuts
1156 fCutSettings->Fill("fCutPtMin", fCutPtMin);
1157 fCutSettings->Fill("fCutPtMax", fCutPtMax);
1158 fCutSettings->Fill("fCutEtaMin", fCutEtaMin);
1159 fCutSettings->Fill("fCutEtaMax", fCutEtaMax);
1160
1161 // track quality cut variables
1162 fCutSettings->Fill("fFilterBit", fFilterBit);
1163 if(fUseRelativeCuts) fCutSettings->Fill("fUseRelativeCuts", 1);
1164 if(fCutRequireTPCRefit) fCutSettings->Fill("fCutRequireTPCRefit", 1);
1165 if(fCutRequireITSRefit) fCutSettings->Fill("fCutRequireITSRefit", 1);
1166
1167 fCutSettings->Fill("fCutMinNumberOfClusters", fCutMinNumberOfClusters);
1168 fCutSettings->Fill("fCutPercMinNumberOfClusters", fCutPercMinNumberOfClusters);
1169 fCutSettings->Fill("fCutMinNumberOfCrossedRows", fCutMinNumberOfCrossedRows);
1170 fCutSettings->Fill("fCutPercMinNumberOfCrossedRows", fCutPercMinNumberOfCrossedRows);
1171
1172 fCutSettings->Fill("fCutMinRatioCrossedRowsOverFindableClustersTPC", fCutMinRatioCrossedRowsOverFindableClustersTPC);
1173 fCutSettings->Fill("fCutMaxFractionSharedTPCClusters", fCutMaxFractionSharedTPCClusters);
1174 fCutSettings->Fill("fCutMaxDCAToVertexXY", fCutMaxDCAToVertexXY);
1175 fCutSettings->Fill("fCutMaxChi2PerClusterITS", fCutMaxChi2PerClusterITS);
1176
1177 if(fCutDCAToVertex2D) fCutSettings->Fill("fCutDCAToVertex2D", 1);
1178 if(fCutRequireSigmaToVertex) fCutSettings->Fill("fCutRequireSigmaToVertex",1);
1179 fCutSettings->Fill("fCutMaxDCAToVertexXYPtDepPar0", fCutMaxDCAToVertexXYPtDepPar0);
1180 fCutSettings->Fill("fCutMaxDCAToVertexXYPtDepPar1", fCutMaxDCAToVertexXYPtDepPar1);
1181 fCutSettings->Fill("fCutMaxDCAToVertexXYPtDepPar2", fCutMaxDCAToVertexXYPtDepPar2);
1182
1183 if(fCutAcceptKinkDaughters) fCutSettings->Fill("fCutAcceptKinkDaughters", 1);
1184 fCutSettings->Fill("fCutMaxChi2TPCConstrainedGlobal", fCutMaxChi2TPCConstrainedGlobal);
1185 if(fCutLengthInTPCPtDependent) fCutSettings->Fill("fCutLengthInTPCPtDependent", 1);
1186 fCutSettings->Fill("fPrefactorLengthInTPCPtDependent", fPrefactorLengthInTPCPtDependent);
1187}
1188
b3341d37 1189Bool_t AlidNdPtAnalysisPbPbAOD::GetDCA(const AliAODTrack *track, AliAODEvent *evt, Double_t d0z0[2])
d0483ba3 1190{
b3341d37 1191 // function adapted from AliDielectronVarManager.h
d0483ba3 1192
b3341d37 1193 if(track->TestBit(AliAODTrack::kIsDCA)){
1194 d0z0[0]=track->DCA();
1195 d0z0[1]=track->ZAtDCA();
1196 return kTRUE;
d0483ba3 1197 }
1198
1199 Bool_t ok=kFALSE;
b3341d37 1200 if(evt) {
d0483ba3 1201 Double_t covd0z0[3];
b3341d37 1202 //AliAODTrack copy(*track);
1203 AliExternalTrackParam etp; etp.CopyFromVTrack(track);
d0483ba3 1204
1205 Float_t xstart = etp.GetX();
b3341d37 1206 if(xstart>3.) {
d0483ba3 1207 d0z0[0]=-999.;
1208 d0z0[1]=-999.;
1209 //printf("This method can be used only for propagation inside the beam pipe \n");
b3341d37 1210 return kFALSE;
d0483ba3 1211 }
b3341d37 1212
1213
d0483ba3 1214 AliAODVertex *vtx =(AliAODVertex*)(evt->GetPrimaryVertex());
1215 Double_t fBzkG = evt->GetMagneticField(); // z componenent of field in kG
1216 ok = etp.PropagateToDCA(vtx,fBzkG,kVeryBig,d0z0,covd0z0);
b3341d37 1217 //ok = copy.PropagateToDCA(vtx,fBzkG,kVeryBig,d0z0,covd0z0);
d0483ba3 1218 }
1219 if(!ok){
1220 d0z0[0]=-999.;
1221 d0z0[1]=-999.;
d0483ba3 1222 }
b3341d37 1223 return ok;
d0483ba3 1224}
1225
b3341d37 1226
d25bcbe6 1227Bool_t AlidNdPtAnalysisPbPbAOD::IsMCTrackAccepted(AliAODMCParticle *part)
1228{
1229 if(!part) return kFALSE;
1230
1231 Double_t charge = part->Charge()/3.;
1232 if (TMath::Abs(charge) < 0.001) return kFALSE;
1233
1234 return kTRUE;
1235}
1236
1237const char * AlidNdPtAnalysisPbPbAOD::GetParticleName(Int_t pdg)
1238{
1239 TParticlePDG * p1 = TDatabasePDG::Instance()->GetParticle(pdg);
1240 if(p1) return p1->GetName();
1241 return Form("%d", pdg);
1242}
1243
1244AliGenHijingEventHeader* AlidNdPtAnalysisPbPbAOD::GetHijingEventHeader(AliAODMCHeader *header)
1245{
1246 //
1247 // inspired by PWGJE/AliPWG4HighPtSpectra.cxx
1248 //
1249
1250 if(!header) return 0x0;
1251 AliGenHijingEventHeader* hijingGenHeader = NULL;
1252
1253 TList* headerList = header->GetCocktailHeaders();
1254
1255 for(Int_t i = 0; i < headerList->GetEntries(); i++)
1256 {
1257 hijingGenHeader = dynamic_cast<AliGenHijingEventHeader*>(headerList->At(i));
1258 if(hijingGenHeader) break;
1259 }
1260
1261 if(!hijingGenHeader) return 0x0;
1262
1263 return hijingGenHeader;
1264}
1265
1266AliGenPythiaEventHeader* AlidNdPtAnalysisPbPbAOD::GetPythiaEventHeader(AliAODMCHeader *header)
1267{
1268 //
1269 // inspired by PWGJE/AliPWG4HighPtSpectra.cxx
1270 //
1271
1272 if(!header) return 0x0;
1273 AliGenPythiaEventHeader* PythiaGenHeader = NULL;
1274
1275 TList* headerList = header->GetCocktailHeaders();
1276
1277 for(Int_t i = 0; i < headerList->GetEntries(); i++)
1278 {
1279 PythiaGenHeader = dynamic_cast<AliGenPythiaEventHeader*>(headerList->At(i));
1280 if(PythiaGenHeader) break;
1281 }
1282
1283 if(!PythiaGenHeader) return 0x0;
1284
1285 return PythiaGenHeader;
1286}
1287
1288//________________________________________________________________________
1289Bool_t AlidNdPtAnalysisPbPbAOD::IsHijingParticle(const AliAODMCParticle *part, AliGenHijingEventHeader* hijingGenHeader){
1290
1291 // Check whether a particle is from Hijing or some injected
ea3cfeda 1292 // returns kFALSE if particle is injected
d25bcbe6 1293
1294 if(part->Label() > (hijingGenHeader->NProduced()-1)) return kFALSE;
1295 return kTRUE;
1296}
1297
1298//________________________________________________________________________
1299Bool_t AlidNdPtAnalysisPbPbAOD::IsPythiaParticle(const AliAODMCParticle *part, AliGenPythiaEventHeader* pythiaGenHeader){
1300
1301 // Check whether a particle is from Pythia or some injected
1302
1303 if(part->Label() > (pythiaGenHeader->NProduced()-1)) return kFALSE;
1304 return kTRUE;
b3341d37 1305}
1306
1307Double_t* AlidNdPtAnalysisPbPbAOD::GetArrayClone(Int_t n, Double_t* source)
1308{
1309 if (!source || n==0) return 0;
1310 Double_t* dest = new Double_t[n];
1311 for (Int_t i=0; i<n ; i++) { dest[i] = source[i]; }
1312 return dest;
d25bcbe6 1313}
1314
b3341d37 1315void AlidNdPtAnalysisPbPbAOD::Terminate(Option_t *)
1316{
1317
1318}
1319
1320