]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGLF/SPECTRA/ChargedHadrons/dNdPt/AlidNdPtAnalysisPbPbAOD.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
added correction in phi, modified AddTask acoordingly
[u/mrichter/AliRoot.git] / PWGLF / SPECTRA / ChargedHadrons / dNdPt / AlidNdPtAnalysisPbPbAOD.cxx
CommitLineData
d25bcbe6 1/**************************************************************************
ab1375ce 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**************************************************************************/
d25bcbe6 15//------------------------------------------------------------------------------
16// AlidNdPtAnalysisPbPbAOD class.
17//
18// Author: P. Luettig, 15.05.2013
ab1375ce 19// last modified: 10.06.2014
d25bcbe6 20//------------------------------------------------------------------------------
3dd0b8f4 21/*
ab1375ce 22* This task analysis measured data in PbPb collisions stored in AODs and extract
23* transverse momentum spectra for unidentified charged hadrons vs. centrality.
24* Based on MC the efficiency and secondary contamination are determined,
25* to correct the measured pT distribution.
26* Histograms for the pT resolution correction are also filled.
27*
28*/
bc80e684 29
d25bcbe6 30
31#include "AlidNdPtAnalysisPbPbAOD.h"
32
8a4ab847 33#include "AliAnalysisTaskSE.h"
d25bcbe6 34
35using namespace std;
36
37ClassImp(AlidNdPtAnalysisPbPbAOD)
38
d25bcbe6 39AlidNdPtAnalysisPbPbAOD::AlidNdPtAnalysisPbPbAOD(const char *name) : AliAnalysisTaskSE(name),
40fOutputList(0),
41// Histograms
ea3cfeda
PL		 42fPt(0),
43fMCPt(0),
44fZvPtEtaCent(0),
4ffd6c6e 45fDeltaphiPtEtaPhiCent(0),
a0036e80 46fPtResptCent(0),
ea3cfeda
PL		 47fMCRecPrimZvPtEtaCent(0),
48fMCGenZvPtEtaCent(0),
49fMCRecSecZvPtEtaCent(0),
4ffd6c6e 50fMCRecPrimDeltaphiPtEtaPhiCent(0),
51fMCGenDeltaphiPtEtaPhiCent(0),
52fMCRecSecDeltaphiPtEtaPhiCent(0),
ea3cfeda
PL		 53fEventStatistics(0),
54fEventStatisticsCentrality(0),
55fMCEventStatisticsCentrality(0),
56fAllEventStatisticsCentrality(0),
57fEventStatisticsCentralityTrigger(0),
58fZvMultCent(0),
59fTriggerStatistics(0),
ea3cfeda
PL		 60fCharge(0),
61fMCCharge(0),
ea3cfeda
PL		 62fDCAPtAll(0),
63fDCAPtAccepted(0),
64fMCDCAPtSecondary(0),
65fMCDCAPtPrimary(0),
66fCutPercClusters(0),
67fCutPercCrossed(0),
68fCrossCheckRowsLength(0),
69fCrossCheckClusterLength(0),
70fCrossCheckRowsLengthAcc(0),
71fCrossCheckClusterLengthAcc(0),
aa7eca65 72fCrossCheckPtresLength(0),
73fCrossCheckPtresRows(0),
8a4ab847 74fCutSettings(0),
bc80e684 75fEventplaneDist(0),
89f04ae7 76fEventplaneRunDist(0),
bc80e684 77fMCEventplaneDist(0),
1444967d 78fCorrelEventplaneMCDATA(0),
96ebdea7 79fCorrelEventplaneDefaultCorrected(0),
80fEventplaneSubtractedPercentage(0),
ab1375ce 81// cross check for event plane resolution
82fEPDistCent(0),
83fPhiCent(0),
84fPcosEPCent(0),
85fPsinEPCent(0),
86fPcosPhiCent(0),
87fPsinPhiCent(0),
82a24e4d 88// cross check for event plane determination
89fDeltaPhiCent(0),
b7741813 90fCrossCheckFilterBitPhiCent(0),
d25bcbe6 91//global
ea3cfeda 92fIsMonteCarlo(0),
0d3e3f7e 93fEPselector("Q"),
d25bcbe6 94// event cut variables
ea3cfeda 95fCutMaxZVertex(10.),
d25bcbe6 96// track kinematic cut variables
ea3cfeda
PL		 97fCutPtMin(0.15),
98fCutPtMax(200.),
99fCutEtaMin(-0.8),
100fCutEtaMax(0.8),
d25bcbe6 101// track quality cut variables
a0036e80 102fFilterBit(AliAODTrack::kTrkGlobal),
ea3cfeda
PL		 103fUseRelativeCuts(kFALSE),
104fCutRequireTPCRefit(kTRUE),
a0036e80 105fCutRequireITSRefit(kTRUE),
ea3cfeda
PL		 106fCutMinNumberOfClusters(60),
107fCutPercMinNumberOfClusters(0.2),
108fCutMinNumberOfCrossedRows(120.),
109fCutPercMinNumberOfCrossedRows(0.2),
110fCutMinRatioCrossedRowsOverFindableClustersTPC(0.8),
111fCutMaxChi2PerClusterTPC(4.),
112fCutMaxFractionSharedTPCClusters(0.4),
113fCutMaxDCAToVertexZ(3.0),
114fCutMaxDCAToVertexXY(3.0),
ea3cfeda
PL		 115fCutMaxChi2PerClusterITS(36.),
116fCutDCAToVertex2D(kFALSE),
117fCutRequireSigmaToVertex(kFALSE),
118fCutMaxDCAToVertexXYPtDepPar0(0.0182),
119fCutMaxDCAToVertexXYPtDepPar1(0.0350),
120fCutMaxDCAToVertexXYPtDepPar2(1.01),
121fCutAcceptKinkDaughters(kFALSE),
122fCutMaxChi2TPCConstrainedGlobal(36.),
123fCutLengthInTPCPtDependent(kFALSE),
124fPrefactorLengthInTPCPtDependent(1),
d25bcbe6 125// binning for THnSparse
126fMultNbins(0),
127fPtNbins(0),
128fPtCorrNbins(0),
72bb4ceb 129fPtCheckNbins(0),
d25bcbe6 130fEtaNbins(0),
72bb4ceb 131fEtaCheckNbins(0),
d25bcbe6 132fZvNbins(0),
133fCentralityNbins(0),
72bb4ceb 134fPhiNbins(0),
4ffd6c6e 135fDeltaphiNbins(0),
89f04ae7 136fRunNumberNbins(0),
d25bcbe6 137fBinsMult(0),
138fBinsPt(0),
139fBinsPtCorr(0),
72bb4ceb 140fBinsPtCheck(0),
d25bcbe6 141fBinsEta(0),
72bb4ceb 142fBinsEtaCheck(0),
d25bcbe6 143fBinsZv(0),
72bb4ceb 144fBinsCentrality(0),
89f04ae7 145fBinsPhi(0),
4ffd6c6e 146fBinsDeltaphi(0),
89f04ae7 147fBinsRunNumber(0)
d25bcbe6 148{
72bb4ceb 149
150 for(Int_t i = 0; i < cqMax; i++)
151 {
3dd0b8f4 152 fCrossCheckAll[i] = 0;
153 fCrossCheckAcc[i] = 0;
72bb4ceb 154 }
155
d25bcbe6 156 fMultNbins = 0;
157 fPtNbins = 0;
158 fPtCorrNbins = 0;
72bb4ceb 159 fPtCheckNbins = 0;
d25bcbe6 160 fEtaNbins = 0;
72bb4ceb 161 fEtaCheckNbins = 0;
d25bcbe6 162 fZvNbins = 0;
163 fCentralityNbins = 0;
89f04ae7 164 fPhiNbins = 0;
4ffd6c6e 165 fDeltaphiNbins = 0;
89f04ae7 166 fRunNumberNbins = 0;
d25bcbe6 167 fBinsMult = 0;
168 fBinsPt = 0;
169 fBinsPtCorr = 0;
72bb4ceb 170 fBinsPtCheck = 0;
d25bcbe6 171 fBinsEta = 0;
72bb4ceb 172 fBinsEtaCheck = 0;
d25bcbe6 173 fBinsZv = 0;
174 fBinsCentrality = 0;
72bb4ceb 175 fBinsPhi = 0;
4ffd6c6e 176 fBinsDeltaphi = 0;
89f04ae7 177 fBinsRunNumber = 0;
d25bcbe6 178
179 DefineOutput(1, TList::Class());
180}
181
182// destructor
183AlidNdPtAnalysisPbPbAOD::~AlidNdPtAnalysisPbPbAOD()
3dd0b8f4 184{
185 //
186 // because task is owner of the output list, all objects are deleted, when list->Clear() is called
187 //
8a4ab847 188 if(fOutputList)
189 {
3dd0b8f4 190 fOutputList->Clear();
191 delete fOutputList;
8a4ab847 192 }
193 fOutputList = 0;
d25bcbe6 194}
195
196void AlidNdPtAnalysisPbPbAOD::UserCreateOutputObjects()
197{
198 // create all output histograms here
199 OpenFile(1, "RECREATE");
200
201 fOutputList = new TList();
202 fOutputList->SetOwner();
203
204 //define default binning
205 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 };
206 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};
4ffd6c6e 207 Double_t binsPtCorrDefault[51] = {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, 200.0};
d25bcbe6 208 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};
3dd0b8f4 209 Double_t binsZvDefault[7] = {-30.,-10.,-5.,0.,5.,10.,30.};
d25bcbe6 210 Double_t binsCentralityDefault[12] = {0., 5., 10., 20., 30., 40., 50., 60., 70., 80., 90., 100.};
211
3dd0b8f4 212 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()};
72bb4ceb 213
4ffd6c6e 214 Double_t binsDeltaphiDefault[9] = { 0, 1./16.*TMath::Pi(), 2./16.*TMath::Pi(), 3./16.*TMath::Pi(), 4./16.*TMath::Pi(), 5./16.*TMath::Pi(), 6./16.*TMath::Pi(), 7./16.*TMath::Pi(), 8./16.*TMath::Pi()};
215
216
72bb4ceb 217 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};
218 Double_t binsEtaCheckDefault[7] = {-1.0,-0.8,-0.4,0.,0.4,0.8,1.0};
219
89f04ae7 220 Double_t binsRunNumbers2011[186] = {
221 167693, 167706, 167711, 167712, 167713, 167806, 167807, 167808, 167813, 167814, 167818, 167841, 167842, 167844, 167846, 167902, 167903, 167909, 167915, 167920, 167921, 167985, 167986, 167987, 167988, 168066, 168068, 168069, 168076, 168103, 168104, 168105, 168107, 168108, 168115, 168171, 168172, 168173, 168175, 168177, 168181, 168203, 168204, 168205, 168206, 168207, 168208, 168212, 168213, 168310, 168311, 168318, 168322, 168325, 168341, 168342, 168356, 168361, 168362, 168458, 168460, 168461, 168464, 168467, 168511, 168512, 168514, 168644, 168777, 168826, 168984, 168988, 168992, 169035, 169040, 169044, 169045, 169091, 169094, 169099, 169138, 169143, 169144, 169145, 169148, 169156, 169160, 169167, 169236, 169238, 169377, 169382, 169411, 169415, 169417, 169418, 169419, 169420, 169475, 169498, 169504, 169506, 169512, 169515, 169550, 169553, 169554, 169555, 169557, 169584, 169586, 169587, 169588, 169590, 169591, 169628, 169683, 169835, 169837, 169838, 169846, 169855, 169858, 169859, 169914, 169918, 169919, 169920, 169922, 169923, 169924, 169926, 169956, 169961, 169965, 169969, 169975, 169981, 170027, 170036, 170038, 170040, 170081, 170083, 170084, 170085, 170088, 170089, 170091, 170152, 170155, 170159, 170162, 170163, 170193, 170195, 170203, 170204, 170205, 170207, 170208, 170228, 170230, 170264, 170267, 170268, 170269, 170270, 170306, 170308, 170309, 170311, 170312, 170313, 170315, 170374, 170387, 170388, 170389, 170390, 170546, 170552, 170556, 170572, 170593, 170593+1
222 };
223
d25bcbe6 224 // if no binning is set, use the default
4ffd6c6e 225 if (!fBinsMult) { SetBinsMult(48,binsMultDefault); }
226 if (!fBinsPt) { SetBinsPt(82,binsPtDefault); }
227 if (!fBinsPtCorr) { SetBinsPtCorr(51,binsPtCorrDefault); }
72bb4ceb 228 if (!fBinsPtCheck) { SetBinsPtCheck(20,binsPtCheckDefault); }
4ffd6c6e 229 if (!fBinsEta) { SetBinsEta(31,binsEtaDefault); }
72bb4ceb 230 if (!fBinsEtaCheck) { SetBinsEtaCheck(7,binsEtaCheckDefault); }
4ffd6c6e 231 if (!fBinsZv) { SetBinsZv(7,binsZvDefault); }
d25bcbe6 232 if (!fBinsCentrality) { SetBinsCentrality(12,binsCentralityDefault); }
4ffd6c6e 233 if (!fBinsPhi) { SetBinsPhi(37,binsPhiDefault); }
234 if (!fBinsDeltaphi) { SetBinsDeltaphi(9,binsDeltaphiDefault); }
235 if (!fBinsRunNumber) { SetBinsRunNumber(186, binsRunNumbers2011); }
d25bcbe6 236
237 Int_t binsZvPtEtaCent[4]={fZvNbins-1,fPtNbins-1,fEtaNbins-1,fCentralityNbins-1};
4ffd6c6e 238 Int_t binsPhiPtEtaCent[5]={fDeltaphiNbins-1,fPtNbins-1,fEtaNbins-1,fPhiNbins-1,fCentralityNbins-1};
239 Int_t binsPhiPtCorrEtaCent[5]={fDeltaphiNbins-1,fPtCorrNbins-1,fEtaNbins-1,fPhiNbins-1,fCentralityNbins-1};
d25bcbe6 240 Int_t binsZvMultCent[3]={fZvNbins-1,fMultNbins-1,fCentralityNbins-1};
241
a0036e80 242 Int_t binsOneOverPtPtResCent[3]={400,300,11};
243 Double_t minbinsOneOverPtPtResCent[3]={0,0,0};
244 Double_t maxbinsOneOverPtPtResCent[3]={1,0.015,100};
245
d25bcbe6 246 // define Histograms
ea3cfeda
PL		 247 fZvPtEtaCent = new THnSparseF("fZvPtEtaCent","Zv:Pt:Eta:Centrality",4,binsZvPtEtaCent);
248 fZvPtEtaCent->SetBinEdges(0,fBinsZv);
249 fZvPtEtaCent->SetBinEdges(1,fBinsPt);
250 fZvPtEtaCent->SetBinEdges(2,fBinsEta);
251 fZvPtEtaCent->SetBinEdges(3,fBinsCentrality);
252 fZvPtEtaCent->GetAxis(0)->SetTitle("Zv (cm)");
253 fZvPtEtaCent->GetAxis(1)->SetTitle("Pt (GeV/c)");
254 fZvPtEtaCent->GetAxis(2)->SetTitle("Eta");
255 fZvPtEtaCent->GetAxis(3)->SetTitle("Centrality");
256 fZvPtEtaCent->Sumw2();
257
4ffd6c6e 258 fDeltaphiPtEtaPhiCent = new THnSparseF("fDeltaphiPtEtaPhiCent","Deltaphi:Pt:Eta:Phi:Centrality",5,binsPhiPtEtaCent);
259 fDeltaphiPtEtaPhiCent->SetBinEdges(0,fBinsDeltaphi);
260 fDeltaphiPtEtaPhiCent->SetBinEdges(1,fBinsPt);
261 fDeltaphiPtEtaPhiCent->SetBinEdges(2,fBinsEta);
262 fDeltaphiPtEtaPhiCent->SetBinEdges(3,fBinsPhi);
263 fDeltaphiPtEtaPhiCent->SetBinEdges(4,fBinsCentrality);
264 fDeltaphiPtEtaPhiCent->GetAxis(0)->SetTitle("#Delta phi to ep");
265 fDeltaphiPtEtaPhiCent->GetAxis(1)->SetTitle("Pt (GeV/c)");
266 fDeltaphiPtEtaPhiCent->GetAxis(2)->SetTitle("Eta");
267 fDeltaphiPtEtaPhiCent->GetAxis(3)->SetTitle("Phi");
268 fDeltaphiPtEtaPhiCent->GetAxis(4)->SetTitle("Centrality");
269 fDeltaphiPtEtaPhiCent->Sumw2();
bc2a9da9 270
a0036e80 271 fPtResptCent = new THnSparseF("fPtResptCent","OneOverPt:PtRes:Centrality",3,binsOneOverPtPtResCent, minbinsOneOverPtPtResCent, maxbinsOneOverPtPtResCent);
272 fPtResptCent->SetBinEdges(2, fBinsCentrality);
273 fPtResptCent->GetAxis(0)->SetTitle("1/pT (GeV/c)^{-1}");
274 fPtResptCent->GetAxis(1)->SetTitle("#sigma(1/pT)");
275 fPtResptCent->GetAxis(2)->SetTitle("centrality");
276 fPtResptCent->Sumw2();
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
4ffd6c6e 311 fMCRecPrimDeltaphiPtEtaPhiCent = new THnSparseF("fMCRecPrimDeltaphiPtEtaPhiCent","mcDeltaphi:mcPt:mcEta:mcPhi:Centrality",5,binsPhiPtCorrEtaCent);
312 fMCRecPrimDeltaphiPtEtaPhiCent->SetBinEdges(0,fBinsDeltaphi);
313 fMCRecPrimDeltaphiPtEtaPhiCent->SetBinEdges(1,fBinsPtCorr);
314 fMCRecPrimDeltaphiPtEtaPhiCent->SetBinEdges(2,fBinsEta);
315 fMCRecPrimDeltaphiPtEtaPhiCent->SetBinEdges(3,fBinsPhi);
316 fMCRecPrimDeltaphiPtEtaPhiCent->SetBinEdges(4,fBinsCentrality);
317 fMCRecPrimDeltaphiPtEtaPhiCent->GetAxis(0)->SetTitle("MC #Delta phi to rp");
318 fMCRecPrimDeltaphiPtEtaPhiCent->GetAxis(1)->SetTitle("MC Pt (GeV/c)");
319 fMCRecPrimDeltaphiPtEtaPhiCent->GetAxis(2)->SetTitle("MC Eta");
320 fMCRecPrimDeltaphiPtEtaPhiCent->GetAxis(3)->SetTitle("MC Phi");
321 fMCRecPrimDeltaphiPtEtaPhiCent->GetAxis(4)->SetTitle("Centrality");
322 fMCRecPrimDeltaphiPtEtaPhiCent->Sumw2();
323
324 fMCGenDeltaphiPtEtaPhiCent = new THnSparseF("fMCGenDeltaphiPtEtaPhiCent","mcDeltaphi:mcPt:mcEta:mcPhi:Centrality",5,binsPhiPtCorrEtaCent);
325 fMCGenDeltaphiPtEtaPhiCent->SetBinEdges(0,fBinsDeltaphi);
326 fMCGenDeltaphiPtEtaPhiCent->SetBinEdges(1,fBinsPtCorr);
327 fMCGenDeltaphiPtEtaPhiCent->SetBinEdges(2,fBinsEta);
328 fMCGenDeltaphiPtEtaPhiCent->SetBinEdges(3,fBinsPhi);
329 fMCGenDeltaphiPtEtaPhiCent->SetBinEdges(4,fBinsCentrality);
330 fMCGenDeltaphiPtEtaPhiCent->GetAxis(0)->SetTitle("MC #Delta phi to rp");
331 fMCGenDeltaphiPtEtaPhiCent->GetAxis(1)->SetTitle("MC Pt (GeV/c)");
332 fMCGenDeltaphiPtEtaPhiCent->GetAxis(2)->SetTitle("MC Eta");
333 fMCGenDeltaphiPtEtaPhiCent->GetAxis(3)->SetTitle("MC Phi");
334 fMCGenDeltaphiPtEtaPhiCent->GetAxis(4)->SetTitle("Centrality");
335 fMCGenDeltaphiPtEtaPhiCent->Sumw2();
336
337 fMCRecSecDeltaphiPtEtaPhiCent = new THnSparseF("fMCRecSecDeltaphiPtEtaPhiCent","mcDeltaphi:mcPt:mcEta:mcPhi:Centrality",5,binsPhiPtCorrEtaCent);
338 fMCRecSecDeltaphiPtEtaPhiCent->SetBinEdges(0,fBinsDeltaphi);
339 fMCRecSecDeltaphiPtEtaPhiCent->SetBinEdges(1,fBinsPtCorr);
340 fMCRecSecDeltaphiPtEtaPhiCent->SetBinEdges(2,fBinsEta);
341 fMCRecSecDeltaphiPtEtaPhiCent->SetBinEdges(3,fBinsPhi);
342 fMCRecSecDeltaphiPtEtaPhiCent->SetBinEdges(4,fBinsCentrality);
343 fMCRecSecDeltaphiPtEtaPhiCent->GetAxis(0)->SetTitle("MC Sec #Delta phi to rp");
344 fMCRecSecDeltaphiPtEtaPhiCent->GetAxis(1)->SetTitle("MC Sec Pt (GeV/c)");
345 fMCRecSecDeltaphiPtEtaPhiCent->GetAxis(2)->SetTitle("MC Sec Eta");
346 fMCRecSecDeltaphiPtEtaPhiCent->GetAxis(3)->SetTitle("MC Phi");
347 fMCRecSecDeltaphiPtEtaPhiCent->GetAxis(4)->SetTitle("Centrality");
348 fMCRecSecDeltaphiPtEtaPhiCent->Sumw2();
bc2a9da9 349
ea3cfeda
PL		 350 fPt = new TH1F("fPt","fPt",2000,0,200);
351 fPt->GetXaxis()->SetTitle("p_{T} (GeV/c)");
352 fPt->GetYaxis()->SetTitle("dN/dp_{T}");
353 fPt->Sumw2();
354
355 fMCPt = new TH1F("fMCPt","fMCPt",2000,0,200);
356 fMCPt->GetXaxis()->SetTitle("MC p_{T} (GeV/c)");
357 fMCPt->GetYaxis()->SetTitle("dN/dp_{T}");
358 fMCPt->Sumw2();
359
360 fEventStatistics = new TH1F("fEventStatistics","fEventStatistics",10,0,10);
361 fEventStatistics->GetYaxis()->SetTitle("number of events");
362 fEventStatistics->SetBit(TH1::kCanRebin);
363
364 fEventStatisticsCentrality = new TH1F("fEventStatisticsCentrality","fEventStatisticsCentrality",fCentralityNbins-1, fBinsCentrality);
365 fEventStatisticsCentrality->GetYaxis()->SetTitle("number of events");
366
367 fMCEventStatisticsCentrality = new TH1F("fMCEventStatisticsCentrality","fMCEventStatisticsCentrality",fCentralityNbins-1, fBinsCentrality);
368 fMCEventStatisticsCentrality->GetYaxis()->SetTitle("number of MC events");
369
370 fAllEventStatisticsCentrality = new TH1F("fAllEventStatisticsCentrality","fAllEventStatisticsCentrality",fCentralityNbins-1, fBinsCentrality);
371 fAllEventStatisticsCentrality->GetYaxis()->SetTitle("number of events");
372
373 fEventStatisticsCentralityTrigger = new TH2F("fEventStatisticsCentralityTrigger","fEventStatisticsCentralityTrigger;centrality;trigger",100,0,100,3,0,3);
374 fEventStatisticsCentralityTrigger->Sumw2();
375
376 fZvMultCent = new THnSparseF("fZvMultCent","Zv:mult:Centrality",3,binsZvMultCent);
377 fZvMultCent->SetBinEdges(0,fBinsZv);
378 fZvMultCent->SetBinEdges(1,fBinsMult);
379 fZvMultCent->SetBinEdges(2,fBinsCentrality);
380 fZvMultCent->GetAxis(0)->SetTitle("Zv (cm)");
381 fZvMultCent->GetAxis(1)->SetTitle("N_{acc}");
382 fZvMultCent->GetAxis(2)->SetTitle("Centrality");
383 fZvMultCent->Sumw2();
384
385 fTriggerStatistics = new TH1F("fTriggerStatistics","fTriggerStatistics",10,0,10);
386 fTriggerStatistics->GetYaxis()->SetTitle("number of events");
387
ea3cfeda
PL		 388 fCharge = new TH1F("fCharge","fCharge",30, -5, 5);
389 fCharge->GetXaxis()->SetTitle("Charge");
390 fCharge->GetYaxis()->SetTitle("number of tracks");
391
392 fMCCharge = new TH1F("fMCCharge","fMCCharge",30, -5, 5);
393 fMCCharge->GetXaxis()->SetTitle("MC Charge");
3dd0b8f4 394 fMCCharge->GetYaxis()->SetTitle("number of tracks");
95f26ffa 395
3dd0b8f4 396 Int_t binsDCAxyDCAzPtEtaPhi[6] = { 10 , 10 , fPtCheckNbins-1, fEtaCheckNbins-1, 18, fCentralityNbins-1 };
397 Double_t minDCAxyDCAzPtEtaPhi[6] = { -5 , -5 , 0, -1.5, 0., 0 };
398 Double_t maxDCAxyDCAzPtEtaPhi[6] = { 5., 5., 100, 1.5, 2.*TMath::Pi(), 100 };
d25bcbe6 399
ea3cfeda
PL		 400 fDCAPtAll = new THnSparseF("fDCAPtAll","fDCAPtAll",6, binsDCAxyDCAzPtEtaPhi, minDCAxyDCAzPtEtaPhi, maxDCAxyDCAzPtEtaPhi);
401 fDCAPtAccepted = new THnSparseF("fDCAPtAccepted","fDCAPtAccepted",6, binsDCAxyDCAzPtEtaPhi, minDCAxyDCAzPtEtaPhi, maxDCAxyDCAzPtEtaPhi);
402 fMCDCAPtSecondary = new THnSparseF("fMCDCAPtSecondary","fMCDCAPtSecondary",6, binsDCAxyDCAzPtEtaPhi, minDCAxyDCAzPtEtaPhi, maxDCAxyDCAzPtEtaPhi);
403 fMCDCAPtPrimary = new THnSparseF("fMCDCAPtPrimary","fMCDCAPtPrimary",6, binsDCAxyDCAzPtEtaPhi, minDCAxyDCAzPtEtaPhi, maxDCAxyDCAzPtEtaPhi);
404
405 fDCAPtAll->SetBinEdges(2, fBinsPtCheck);
406 fDCAPtAccepted->SetBinEdges(2, fBinsPtCheck);
407 fMCDCAPtSecondary->SetBinEdges(2, fBinsPtCheck);
408 fMCDCAPtPrimary->SetBinEdges(2, fBinsPtCheck);
409
410 fDCAPtAll->SetBinEdges(3, fBinsEtaCheck);
411 fDCAPtAccepted->SetBinEdges(3, fBinsEtaCheck);
412 fMCDCAPtSecondary->SetBinEdges(3, fBinsEtaCheck);
413 fMCDCAPtPrimary->SetBinEdges(3, fBinsEtaCheck);
414
415 fDCAPtAll->SetBinEdges(5, fBinsCentrality);
416 fDCAPtAccepted->SetBinEdges(5, fBinsCentrality);
417 fMCDCAPtSecondary->SetBinEdges(5, fBinsCentrality);
418 fMCDCAPtPrimary->SetBinEdges(5, fBinsCentrality);
419
420 fDCAPtAll->Sumw2();
421 fDCAPtAccepted->Sumw2();
422 fMCDCAPtSecondary->Sumw2();
423 fMCDCAPtPrimary->Sumw2();
424
425 fDCAPtAll->GetAxis(0)->SetTitle("DCA_{xy} (cm)");
426 fDCAPtAll->GetAxis(1)->SetTitle("DCA_{z} (cm)");
427 fDCAPtAll->GetAxis(2)->SetTitle("p_{T} (GeV/c)");
428 fDCAPtAll->GetAxis(3)->SetTitle("#eta");
429 fDCAPtAll->GetAxis(4)->SetTitle("#phi");
430 fDCAPtAll->GetAxis(5)->SetTitle("Centrality");
431
432 fDCAPtAccepted->GetAxis(0)->SetTitle("DCA_{xy} (cm)");
433 fDCAPtAccepted->GetAxis(1)->SetTitle("DCA_{z} (cm)");
434 fDCAPtAccepted->GetAxis(2)->SetTitle("p_{T} (GeV/c)");
435 fDCAPtAccepted->GetAxis(3)->SetTitle("#eta");
436 fDCAPtAccepted->GetAxis(4)->SetTitle("#phi");
437 fDCAPtAccepted->GetAxis(5)->SetTitle("Centrality");
438
439 fMCDCAPtSecondary->GetAxis(0)->SetTitle("DCA_{xy} (cm)");
440 fMCDCAPtSecondary->GetAxis(1)->SetTitle("DCA_{z} (cm)");
441 fMCDCAPtSecondary->GetAxis(2)->SetTitle("p_{T} (GeV/c)");
442 fMCDCAPtSecondary->GetAxis(3)->SetTitle("#eta");
443 fMCDCAPtSecondary->GetAxis(4)->SetTitle("#phi");
444 fMCDCAPtSecondary->GetAxis(5)->SetTitle("Centrality");
445
446 fMCDCAPtPrimary->GetAxis(0)->SetTitle("DCA_{xy} (cm)");
447 fMCDCAPtPrimary->GetAxis(1)->SetTitle("DCA_{z} (cm)");
448 fMCDCAPtPrimary->GetAxis(2)->SetTitle("p_{T} (GeV/c)");
449 fMCDCAPtPrimary->GetAxis(3)->SetTitle("#eta");
450 fMCDCAPtPrimary->GetAxis(4)->SetTitle("#phi");
451 fMCDCAPtPrimary->GetAxis(5)->SetTitle("Centrality");
9db7eb94 452
72bb4ceb 453
454 char cFullTempTitle[255];
455 char cTempTitleAxis0All[255];
456 char cTempTitleAxis0Acc[255];
457 // char cTempTitleAxis1[255];
458 char cFullTempName[255];
459 char cTempNameAxis0[255];
460 // char cTempNameAxis1[255];
461 const Int_t iNbinRowsClusters = 21;
462 // 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.};
463
464 const Int_t iNbinChi = 51;
ea3cfeda 465 const Int_t iNbinLength = 165;
bc80e684 466 const Int_t iNbinRowsOverClusters = 60;
72bb4ceb 467 // 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.};
468
469 Int_t iNbin = 0;
470 // Double_t *dBins = 0x0;
471 Double_t dBinMin = 0;
472 Double_t dBinMax = 0;
473
474 for(Int_t iCheckQuant = 0; iCheckQuant < cqMax; iCheckQuant++)
475 {
3dd0b8f4 476 // iCheckQuant: 0 = CrossedRows, 1 = Nclusters, 2 = Chi^2/clusterTPC
477 if(iCheckQuant == cqCrossedRows)
478 {
479 snprintf(cTempTitleAxis0All,255, "NcrossedRows before Cut");
480 snprintf(cTempTitleAxis0Acc,255, "NcrossedRows after Cut");
481 snprintf(cTempNameAxis0,255, "CrossedRows");
482 iNbin = iNbinRowsClusters;
483 dBinMin = 0;
484 dBinMax = 159.;
485 }
486 else if(iCheckQuant == cqNcluster)
487 {
488 snprintf(cTempTitleAxis0All,255, "Nclusters before Cut");
489 snprintf(cTempTitleAxis0Acc,255, "Nclusters after Cut");
490 snprintf(cTempNameAxis0,255, "Clusters");
491 iNbin = iNbinRowsClusters;
492 dBinMin = 0;
493 dBinMax = 159.;
494 }
495 else if(iCheckQuant == cqChi)
496 {
497 snprintf(cTempTitleAxis0All,255, "#Chi^{2}/cluster before Cut");
498 snprintf(cTempTitleAxis0Acc,255, "#Chi^{2}/cluster after Cut");
499 snprintf(cTempNameAxis0,255, "Chi");
500 iNbin = iNbinChi;
501 dBinMin = 0;
502 dBinMax = 10.;
503 }
504 else if(iCheckQuant == cqLength)
505 {
506 snprintf(cTempTitleAxis0All,255, "Length in TPC before Cut (cm)");
507 snprintf(cTempTitleAxis0Acc,255, "Length in TPC after Cut (cm)");
508 snprintf(cTempNameAxis0,255, "Length");
509 iNbin = iNbinLength;
510 dBinMin = 0;
511 dBinMax = 165.;
512 }
bc80e684 513 else if(iCheckQuant == cqRowsOverFindable)
514 {
515 snprintf(cTempTitleAxis0All,255, "Number of Crossed Rows / Number of Findable Clusters before Cut");
516 snprintf(cTempTitleAxis0Acc,255, "Number of Crossed Rows / Number of Findable Clusters before Cut");
517 snprintf(cTempNameAxis0,255, "RowsOverFindable");
518 iNbin = iNbinRowsOverClusters;
519 dBinMin = 0.6;
520 dBinMax = 1.2;
521 }
522
3dd0b8f4 523
524 Int_t binsCheckPtEtaPhi[5] = { iNbin, fPtCheckNbins-1, fEtaCheckNbins-1, 18, fCentralityNbins-1};
525 // Int_t binsCheckPtEtaPhi[5] = { iNbin, fPtNbins-1, fEtaCheckNbins-1, 18, fCentralityNbins-1};
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
529 snprintf(cFullTempName, 255, "f%sPtEtaPhiAll",cTempNameAxis0);
530 snprintf(cFullTempTitle, 255,"%s;%s;p_{T} (GeV/c);#eta;#phi;Centrality", cFullTempName, cTempTitleAxis0All);
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();
535
536 snprintf(cFullTempName, 255, "f%sPtEtaPhiAcc",cTempNameAxis0);
537 snprintf(cFullTempTitle, 255,"%s;%s;p_{T} (GeV/c);#eta;#phi;Centrality", cFullTempName, cTempTitleAxis0Acc);
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
aa7eca65 561 fCrossCheckPtresLength = new TH2F("fCrossCheckPtresLength","fCrossCheckPtresLength;Length in TPC;#sigma(1/pT)*pT",170,0,170,100,0,1);
562 fCrossCheckPtresLength->Sumw2();
563
564 fCrossCheckPtresRows = new TH2F("fCrossCheckPtresRows","fCrossCheckPtresRows;NcrossedRows;#sigma(1/pT)*pT",170,0,170,100,0,1);
565 fCrossCheckPtresRows->Sumw2();
566
8a4ab847 567 fCutSettings = new TH1F("fCutSettings","fCutSettings",100,0,10);
568 fCutSettings->GetYaxis()->SetTitle("cut value");
569 fCutSettings->SetBit(TH1::kCanRebin);
ea3cfeda 570
227a6341 571 fEventplaneDist = new TH1F("fEventplaneDist","fEventplaneDist",200, 0, 2.*TMath::Pi());
bc80e684 572 fEventplaneDist->GetXaxis()->SetTitle("#phi (event plane)");
573 fEventplaneDist->Sumw2();
574
227a6341 575 fEventplaneRunDist = new TH2F("fEventplaneRunDist","fEventplaneRunDist",200, 0, 2.*TMath::Pi(),fRunNumberNbins-1, fBinsRunNumber );
89f04ae7 576 fEventplaneRunDist->GetXaxis()->SetTitle("#phi (event plane)");
577 fEventplaneRunDist->GetYaxis()->SetTitle("runnumber");
578 fEventplaneRunDist->Sumw2();
579
bc80e684 580 fMCEventplaneDist = new TH1F("fMCEventplaneDist","fMCEventplaneDist",20, -1.*TMath::Pi(), TMath::Pi());
581 fMCEventplaneDist->GetXaxis()->SetTitle("#phi (MC event plane)");
582 fMCEventplaneDist->Sumw2();
583
1444967d 584 fCorrelEventplaneMCDATA = new TH2F("fCorrelEventplaneMCDATA","fCorrelEventplaneMCDATA",40, -2.*TMath::Pi(), 2.*TMath::Pi(), 40, -2.*TMath::Pi(), 2.*TMath::Pi());
585 fCorrelEventplaneMCDATA->GetXaxis()->SetTitle("#phi (event plane)");
586 fCorrelEventplaneMCDATA->GetYaxis()->SetTitle("#phi (MC event plane)");
587 fCorrelEventplaneMCDATA->Sumw2();
588
96ebdea7 589 Int_t binsCorrelPhiPhiCent[3] = { 40, 40, 10};
590 Double_t minCorrelPhiPhiCent[3] = { -2.*TMath::Pi(), -2.*TMath::Pi(), 0};
591 Double_t maxCorrelPhiPhiCent[3] = { 2.*TMath::Pi(), 2.*TMath::Pi(), 100};
592
593 fCorrelEventplaneDefaultCorrected = new THnSparseF("fCorrelEventplaneDefaultCorrected","fCorrelEventplaneDefaultCorrected",3,binsCorrelPhiPhiCent, minCorrelPhiPhiCent, maxCorrelPhiPhiCent);
594 fCorrelEventplaneDefaultCorrected->SetBinEdges(2, fBinsCentrality);
595 fCorrelEventplaneDefaultCorrected->GetAxis(0)->SetTitle("#phi (event plane)");
596 fCorrelEventplaneDefaultCorrected->GetAxis(1)->SetTitle("#phi (corrected event plane)");
597 fCorrelEventplaneDefaultCorrected->GetAxis(2)->SetTitle("centrality");
598 fCorrelEventplaneDefaultCorrected->Sumw2();
599
600 fEventplaneSubtractedPercentage = new TH2F("fEventplaneSubtractedPercentage","fEventplaneSubtractedPercentage",100, 0,1, fCentralityNbins-1, fBinsCentrality);
601 fEventplaneSubtractedPercentage->GetXaxis()->SetTitle("percentage of tracks, which have been subtracted during analysis");
602 fEventplaneSubtractedPercentage->GetYaxis()->SetTitle("centrality");
603 fEventplaneSubtractedPercentage->Sumw2();
604
ab1375ce 605 // cross check for event plane resolution
227a6341 606 fEPDistCent = new TH2F("fEPDistCent","fEPDistCent",20, -2.*TMath::Pi(), 2.*TMath::Pi(), fCentralityNbins-1, fBinsCentrality);
ab1375ce 607 fEPDistCent->GetXaxis()->SetTitle("#phi (#Psi_{EP})");
608 fEPDistCent->GetYaxis()->SetTitle("Centrality");
609 fEPDistCent->Sumw2();
610
611 fPhiCent = new TH2F("fPhiCent","fPhiCent",200, -2.*TMath::Pi(), 2.*TMath::Pi(), fCentralityNbins-1, fBinsCentrality);
612 fPhiCent->GetXaxis()->SetTitle("#phi");
613 fPhiCent->GetYaxis()->SetTitle("Centrality");
614 fPhiCent->Sumw2();
615
60134a87 616 fPcosEPCent = new TProfile("fPcosEPCent","fPcosEPCent", 100,0,100);
ab1375ce 617 fPcosEPCent->GetXaxis()->SetTitle("Centrality");
618 fPcosEPCent->GetYaxis()->SetTitle("#LT cos 2 #Psi_{EP} #GT");
619 fPcosEPCent->Sumw2();
620
60134a87 621 fPsinEPCent = new TProfile("fPsinEPCent","fPsinEPCent", 100,0,100);
ab1375ce 622 fPsinEPCent->GetXaxis()->SetTitle("Centrality");
623 fPsinEPCent->GetYaxis()->SetTitle("#LT sin 2 #Psi_{EP} #GT");
624 fPsinEPCent->Sumw2();
625
60134a87 626 fPcosPhiCent = new TProfile("fPcosPhiCent","fPcosPhiCent", 100,0,100);
ab1375ce 627 fPcosPhiCent->GetXaxis()->SetTitle("Centrality");
628 fPcosPhiCent->GetYaxis()->SetTitle("#LT cos 2 #phi #GT");
629 fPcosPhiCent->Sumw2();
630
60134a87 631 fPsinPhiCent = new TProfile("fPsinPhiCent","fPsinPhiCent", 100,0,100);
ab1375ce 632 fPsinPhiCent->GetXaxis()->SetTitle("Centrality");
633 fPsinPhiCent->GetYaxis()->SetTitle("#LT sin 2 #phi #GT");
634 fPsinPhiCent->Sumw2();
635
82a24e4d 636 fDeltaPhiCent = new TH2F("fDeltaPhiCent","fDeltaPhiCent",200, -2.*TMath::Pi(), 2.*TMath::Pi(), fCentralityNbins-1, fBinsCentrality);
637 fDeltaPhiCent->GetXaxis()->SetTitle("#Delta #phi");
638 fDeltaPhiCent->GetYaxis()->SetTitle("Centrality");
639 fDeltaPhiCent->Sumw2();
640
b7741813 641 Int_t binsFilterBitPhiCent[3]={3,200,fCentralityNbins-1};
642 Double_t minbinsFilterBitPhiCent[3]={0,-2.*TMath::Pi(),0};
643 Double_t maxbinsFilterBitPhiCent[3]={3,2.*TMath::Pi(),100};
644
645 fCrossCheckFilterBitPhiCent = new THnSparseF("fCrossCheckFilterBitPhiCent","fCrossCheckFilterBitPhiCent",3, binsFilterBitPhiCent, minbinsFilterBitPhiCent, maxbinsFilterBitPhiCent);
646 fCrossCheckFilterBitPhiCent->SetBinEdges(2,fBinsCentrality);
647 fCrossCheckFilterBitPhiCent->GetAxis(0)->SetTitle("FilterBit");
648 fCrossCheckFilterBitPhiCent->GetAxis(1)->SetTitle("#phi");
649 fCrossCheckFilterBitPhiCent->GetAxis(2)->SetTitle("Centrality");
650 fCrossCheckFilterBitPhiCent->Sumw2();
651
d25bcbe6 652 // Add Histos, Profiles etc to List
ea3cfeda 653 fOutputList->Add(fZvPtEtaCent);
4ffd6c6e 654 fOutputList->Add(fDeltaphiPtEtaPhiCent);
a0036e80 655 fOutputList->Add(fPtResptCent);
ea3cfeda
PL		 656 fOutputList->Add(fPt);
657 fOutputList->Add(fMCRecPrimZvPtEtaCent);
658 fOutputList->Add(fMCGenZvPtEtaCent);
659 fOutputList->Add(fMCRecSecZvPtEtaCent);
4ffd6c6e 660 fOutputList->Add(fMCRecPrimDeltaphiPtEtaPhiCent);
661 fOutputList->Add(fMCGenDeltaphiPtEtaPhiCent);
662 fOutputList->Add(fMCRecSecDeltaphiPtEtaPhiCent);
ea3cfeda
PL		 663 fOutputList->Add(fMCPt);
664 fOutputList->Add(fEventStatistics);
665 fOutputList->Add(fEventStatisticsCentrality);
666 fOutputList->Add(fMCEventStatisticsCentrality);
667 fOutputList->Add(fAllEventStatisticsCentrality);
668 fOutputList->Add(fEventStatisticsCentralityTrigger);
669 fOutputList->Add(fZvMultCent);
670 fOutputList->Add(fTriggerStatistics);
ea3cfeda
PL		 671 fOutputList->Add(fCharge);
672 fOutputList->Add(fMCCharge);
ea3cfeda
PL		 673 fOutputList->Add(fDCAPtAll);
674 fOutputList->Add(fDCAPtAccepted);
675 fOutputList->Add(fMCDCAPtSecondary);
676 fOutputList->Add(fMCDCAPtPrimary);
72bb4ceb 677 for(Int_t i = 0; i < cqMax; i++)
678 {
3dd0b8f4 679 fOutputList->Add(fCrossCheckAll[i]);
680 fOutputList->Add(fCrossCheckAcc[i]);
72bb4ceb 681 }
ea3cfeda
PL		 682 fOutputList->Add(fCutPercClusters);
683 fOutputList->Add(fCutPercCrossed);
684 fOutputList->Add(fCrossCheckRowsLength);
685 fOutputList->Add(fCrossCheckClusterLength);
686 fOutputList->Add(fCrossCheckRowsLengthAcc);
687 fOutputList->Add(fCrossCheckClusterLengthAcc);
aa7eca65 688 fOutputList->Add(fCrossCheckPtresLength);
689 fOutputList->Add(fCrossCheckPtresRows);
8a4ab847 690 fOutputList->Add(fCutSettings);
bc80e684 691 fOutputList->Add(fEventplaneDist);
89f04ae7 692 fOutputList->Add(fEventplaneRunDist);
bc80e684 693 fOutputList->Add(fMCEventplaneDist);
1444967d 694 fOutputList->Add(fCorrelEventplaneMCDATA);
96ebdea7 695 fOutputList->Add(fCorrelEventplaneDefaultCorrected);
696 fOutputList->Add(fEventplaneSubtractedPercentage);
8a4ab847 697
ab1375ce 698 fOutputList->Add(fEPDistCent);
699 fOutputList->Add(fPhiCent);
700 fOutputList->Add(fPcosEPCent);
701 fOutputList->Add(fPsinEPCent);
702 fOutputList->Add(fPcosPhiCent);
703 fOutputList->Add(fPsinPhiCent);
704
82a24e4d 705 fOutputList->Add(fDeltaPhiCent);
b7741813 706
707 fOutputList->Add(fCrossCheckFilterBitPhiCent);
708
8a4ab847 709 StoreCutSettingsToHistogram();
d25bcbe6 710
711 PostData(1, fOutputList);
712}
713
714void AlidNdPtAnalysisPbPbAOD::UserExec(Option_t *option)
715{
3dd0b8f4 716 //
d25bcbe6 717 // Main Loop
718 // called for each event
3dd0b8f4 719 //
4ffd6c6e 720 cout << fBinsPhi[fPhiNbins-1] << endl;
ea3cfeda 721 fEventStatistics->Fill("all events",1);
d25bcbe6 722
723 // set ZERO pointers:
724 AliInputEventHandler *inputHandler = NULL;
725 AliAODTrack *track = NULL;
726 AliAODMCParticle *mcPart = NULL;
727 AliAODMCHeader *mcHdr = NULL;
728 AliGenHijingEventHeader *genHijingHeader = NULL;
b3341d37 729 //AliGenPythiaEventHeader *genPythiaHeader = NULL;
bc80e684 730 AliEventplane *ep = NULL;
d25bcbe6 731
3b07fd51 732 TVector2 *epQvector = NULL;
733
d25bcbe6 734 Bool_t bIsEventSelectedMB = kFALSE;
735 Bool_t bIsEventSelectedSemi = kFALSE;
736 Bool_t bIsEventSelectedCentral = kFALSE;
737 Bool_t bIsEventSelected = kFALSE;
738 Bool_t bIsPrimary = kFALSE;
739 Bool_t bIsHijingParticle = kFALSE;
9db7eb94 740 Bool_t bMotherIsHijingParticle = kFALSE;
b3341d37 741 //Bool_t bIsPythiaParticle = kFALSE;
d25bcbe6 742 Bool_t bEventHasATrack = kFALSE;
743 Bool_t bEventHasATrackInRange = kFALSE;
744 Int_t nTriggerFired = 0;
745
746
747 Double_t dMCTrackZvPtEtaCent[4] = {0};
748 Double_t dTrackZvPtEtaCent[4] = {0};
749
4ffd6c6e 750 Double_t dMCTrackDeltaphiPtEtaPhiCent[4] = {0};
751 Double_t dTrackDeltaphiPtEtaPhiCent[4] = {0};
bc2a9da9 752
b3341d37 753 Double_t dDCA[2] = {0};
754
d25bcbe6 755 Double_t dMCEventZv = -100;
756 Double_t dEventZv = -100;
757 Int_t iAcceptedMultiplicity = 0;
bc80e684 758 Double_t dEventplaneAngle = -10;
3b07fd51 759 Double_t dEventplaneAngleCorrected = -10; // event plane angle, where tracks contributing to this angle have been subtracted
bc80e684 760 Double_t dMCEventplaneAngle = -10;
d25bcbe6 761
ea3cfeda 762 fIsMonteCarlo = kFALSE;
d25bcbe6 763
764 AliAODEvent *eventAOD = 0x0;
765 eventAOD = dynamic_cast<AliAODEvent*>( InputEvent() );
766 if (!eventAOD) {
3dd0b8f4 767 AliWarning("ERROR: eventAOD not available \n");
768 return;
d25bcbe6 769 }
770
771 // check, which trigger has been fired
772 inputHandler = (AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
773 bIsEventSelectedMB = ( inputHandler->IsEventSelected() & AliVEvent::kMB);
774 bIsEventSelectedSemi = ( inputHandler->IsEventSelected() & AliVEvent::kSemiCentral);
775 bIsEventSelectedCentral = ( inputHandler->IsEventSelected() & AliVEvent::kCentral);
776
ea3cfeda
PL		 777 if(bIsEventSelectedMB || bIsEventSelectedSemi || bIsEventSelectedCentral) fTriggerStatistics->Fill("all triggered events",1);
778 if(bIsEventSelectedMB) { fTriggerStatistics->Fill("MB trigger",1); nTriggerFired++; }
779 if(bIsEventSelectedSemi) { fTriggerStatistics->Fill("SemiCentral trigger",1); nTriggerFired++; }
780 if(bIsEventSelectedCentral) { fTriggerStatistics->Fill("Central trigger",1); nTriggerFired++; }
781 if(nTriggerFired == 0) { fTriggerStatistics->Fill("No trigger",1); }
d25bcbe6 782
783 bIsEventSelected = ( inputHandler->IsEventSelected() & GetCollisionCandidates() );
784
785 // only take tracks of events, which are triggered
786 if(nTriggerFired == 0) { return; }
787
d68506b8 788
d25bcbe6 789 // if( !bIsEventSelected || nTriggerFired>1 ) return;
790
ea3cfeda 791 // fEventStatistics->Fill("events with only coll. cand.", 1);
d25bcbe6 792
793
794
795 // check if there is a stack, if yes, then do MC loop
796 TList *list = eventAOD->GetList();
797 TClonesArray *stack = 0x0;
798 stack = (TClonesArray*)list->FindObject(AliAODMCParticle::StdBranchName());
799
800 if( stack )
801 {
3dd0b8f4 802 fIsMonteCarlo = kTRUE;
803
804 mcHdr = (AliAODMCHeader*)list->FindObject(AliAODMCHeader::StdBranchName());
805
806 genHijingHeader = GetHijingEventHeader(mcHdr);
807 // genPythiaHeader = GetPythiaEventHeader(mcHdr);
808
809 if(!genHijingHeader) { return; }
810
811 // if(!genPythiaHeader) { return; }
812
bc80e684 813
3dd0b8f4 814 dMCEventZv = mcHdr->GetVtxZ();
815 dMCTrackZvPtEtaCent[0] = dMCEventZv;
aa7eca65 816 dMCEventplaneAngle = genHijingHeader->ReactionPlaneAngle();//MoveEventplane(genHijingHeader->ReactionPlaneAngle());
3dd0b8f4 817 fEventStatistics->Fill("MC all events",1);
bc80e684 818 fMCEventplaneDist->Fill(dMCEventplaneAngle);
d25bcbe6 819 }
820
821 AliCentrality* aCentrality = eventAOD->GetCentrality();
822 Double_t dCentrality = aCentrality->GetCentralityPercentile("V0M");
823
824 if( dCentrality < 0 ) return;
d68506b8 825
826 // protection for bias on pt spectra if all triggers selected
ab1375ce 827 if( (bIsEventSelectedCentral) && (dCentrality > 10) ) return;
828 if( (bIsEventSelectedSemi) && ((dCentrality < 20) || (dCentrality > 50))) return;
d68506b8 829
ea3cfeda 830 fEventStatistics->Fill("after centrality selection",1);
d25bcbe6 831
1444967d 832 // get event plane Angle from AODHeader, default is Q
833 ep = const_cast<AliAODEvent*>(eventAOD)->GetEventplane();
834 if(ep) {
aa7eca65 835 dEventplaneAngle = ep->GetEventplane(GetEventplaneSelector().Data(),eventAOD);//MoveEventplane(ep->GetEventplane(GetEventplaneSelector().Data(),eventAOD));
3b07fd51 836 if(GetEventplaneSelector().CompareTo("Q") == 0)
837 {
838 epQvector = ep->GetQVector();
aa7eca65 839 if(epQvector) dEventplaneAngle = epQvector->Phi()/2.;//MoveEventplane(epQvector->Phi());
3b07fd51 840 }
841 }
842
843 if( (GetEventplaneSelector().CompareTo("Q") == 0) && !epQvector )
844 {
845 AliWarning("ERROR: epQvector not available \n");
846 return;
1444967d 847 }
d25bcbe6 848
1444967d 849 // cout << dEventplaneAngle << endl;
850 fEventplaneDist->Fill(dEventplaneAngle);
89f04ae7 851 fEventplaneRunDist->Fill(dEventplaneAngle, (Double_t)eventAOD->GetRunNumber());
d25bcbe6 852
ab1375ce 853 // fill crosscheck histos
854 fEPDistCent->Fill(dEventplaneAngle, dCentrality);
855 fPcosEPCent->Fill(dCentrality, TMath::Cos(2.*dEventplaneAngle));
856 fPsinEPCent->Fill(dCentrality, TMath::Sin(2.*dEventplaneAngle));
857
d25bcbe6 858 // start with MC truth analysis
ea3cfeda 859 if(fIsMonteCarlo)
d25bcbe6 860 {
3dd0b8f4 861
862 if( dMCEventZv > GetCutMaxZVertex() ) { return; }
863
864 dMCTrackZvPtEtaCent[0] = dMCEventZv;
865
866 fEventStatistics->Fill("MC afterZv cut",1);
867
868 for(Int_t iMCtrack = 0; iMCtrack < stack->GetEntriesFast(); iMCtrack++)
869 {
870 mcPart =(AliAODMCParticle*)stack->At(iMCtrack);
871
872 // check for charge
873 if( !(IsMCTrackAccepted(mcPart)) ) continue;
874
875 if(!IsHijingParticle(mcPart, genHijingHeader)) { continue; }
876
877 if(mcPart->IsPhysicalPrimary() )
878 {
879 // fMCHijingPrim->Fill("IsPhysicalPrimary",1);
880 }
881 else
882 {
883 // fMCHijingPrim->Fill("NOT a primary",1);
884 continue;
885 }
886
887
888 //
889 // ======================== fill histograms ========================
890 dMCTrackZvPtEtaCent[1] = mcPart->Pt();
891 dMCTrackZvPtEtaCent[2] = mcPart->Eta();
892 dMCTrackZvPtEtaCent[3] = dCentrality;
893 fMCGenZvPtEtaCent->Fill(dMCTrackZvPtEtaCent);
894
4ffd6c6e 895 dMCTrackDeltaphiPtEtaPhiCent[0] = RotatePhi(mcPart->Phi(), dEventplaneAngle, fBinsPhi[fPhiNbins-1]); // use eventplane and not reactionplan, similar to centrality vs impact paramter
896 // if( dMCTrackDeltaphiPtEtaPhiCent[0] < 0) dMCTrackDeltaphiPtEtaPhiCent[0] += 2.*TMath::Pi();
897 // else if( dMCTrackDeltaphiPtEtaPhiCent[0] > 2.*TMath::Pi()) dMCTrackDeltaphiPtEtaPhiCent[0] -= 2.*TMath::Pi();
898 dMCTrackDeltaphiPtEtaPhiCent[1] = mcPart->Pt();
899 dMCTrackDeltaphiPtEtaPhiCent[2] = mcPart->Eta();
900 dMCTrackDeltaphiPtEtaPhiCent[3] = mcPart->Phi();
901 dMCTrackDeltaphiPtEtaPhiCent[4] = dCentrality;
902 fMCGenDeltaphiPtEtaPhiCent->Fill(dMCTrackDeltaphiPtEtaPhiCent);
3dd0b8f4 903
904 bEventHasATrack = kTRUE;
905
906
907 if( (dMCTrackZvPtEtaCent[1] > GetCutPtMin() ) &&
908 (dMCTrackZvPtEtaCent[1] < GetCutPtMax() ) &&
909 (dMCTrackZvPtEtaCent[2] > GetCutEtaMin() ) &&
910 (dMCTrackZvPtEtaCent[2] < GetCutEtaMax() ) )
911 {
912 fMCPt->Fill(mcPart->Pt());
913 fMCCharge->Fill(mcPart->Charge()/3.);
914 bEventHasATrackInRange = kTRUE;
915 }
916
917 }
d25bcbe6 918 } // isMonteCarlo
ea3cfeda
PL		 919
920 if(bEventHasATrack) { fEventStatistics->Fill("MC events with tracks",1); }
9db7eb94 921 if(bEventHasATrackInRange)
922 {
3dd0b8f4 923 fEventStatistics->Fill("MC events with tracks in range",1);
924 fMCEventStatisticsCentrality->Fill(dCentrality);
9db7eb94 925 }
d25bcbe6 926 bEventHasATrack = kFALSE;
927 bEventHasATrackInRange = kFALSE;
928
929
3dd0b8f4 930 //
d25bcbe6 931 // Loop over recontructed tracks
3dd0b8f4 932 //
d25bcbe6 933
934 dEventZv = eventAOD->GetPrimaryVertex()->GetZ();
ea3cfeda 935 if( TMath::Abs(dEventZv) > GetCutMaxZVertex() ) return;
d25bcbe6 936
3dd0b8f4 937 // count all events, which are within zv distribution
ea3cfeda 938 fAllEventStatisticsCentrality->Fill(dCentrality/*, nTriggerFired*/);
d25bcbe6 939
ea3cfeda 940 fEventStatistics->Fill("after Zv cut",1);
d25bcbe6 941
942 dTrackZvPtEtaCent[0] = dEventZv;
943
bc80e684 944
bc80e684 945
ea3cfeda
PL		 946 if(AreRelativeCutsEnabled())
947 {
3dd0b8f4 948 if(!SetRelativeCuts(eventAOD)) return;
ea3cfeda
PL		 949 }
950
96ebdea7 951 Int_t iSubtractedTracks = 0;
952
d25bcbe6 953 for(Int_t itrack = 0; itrack < eventAOD->GetNumberOfTracks(); itrack++)
954 {
f15c1f69 955 track = dynamic_cast<AliAODTrack*>(eventAOD->GetTrack(itrack));
956 if(!track) AliFatal("Not a standard AOD");
b7741813 957
3dd0b8f4 958 if(!track) continue;
959
960 mcPart = NULL;
961 dMCTrackZvPtEtaCent[1] = 0;
962 dMCTrackZvPtEtaCent[2] = 0;
963 dMCTrackZvPtEtaCent[3] = 0;
964
4ffd6c6e 965 dMCTrackDeltaphiPtEtaPhiCent[0] = 0;
966 dMCTrackDeltaphiPtEtaPhiCent[1] = 0;
967 dMCTrackDeltaphiPtEtaPhiCent[2] = 0;
968 dMCTrackDeltaphiPtEtaPhiCent[3] = 0;
3dd0b8f4 969
970 bIsPrimary = kFALSE;
971
972 GetDCA(track, eventAOD, dDCA);
973
974 Double_t dDCAxyDCAzPt[5] = { dDCA[0], dDCA[1], track->Pt(), track->Eta(), track->Phi() };
975
976 fDCAPtAll->Fill(dDCAxyDCAzPt);
977
978 if( !(IsTrackAccepted(track, dCentrality, eventAOD->GetMagneticField())) ) continue;
979
980 dTrackZvPtEtaCent[1] = track->Pt();
981 dTrackZvPtEtaCent[2] = track->Eta();
982 dTrackZvPtEtaCent[3] = dCentrality;
983
3b07fd51 984 if(GetEventplaneSelector().CompareTo("Q") == 0)
985 {
986 // subtract track contribution from eventplane
96ebdea7 987 Double_t dX = -1000;
988 Double_t dY = -1000;
3b07fd51 989
990 dX = epQvector->X();
991 dY = epQvector->Y();
96ebdea7 992 if( (dX>-1000) && (dY>-1000) ) // only subtract, if not default!
993 {
994 dX -= ep->GetQContributionX(track);
995 dY -= ep->GetQContributionY(track);
996 iSubtractedTracks++;
997 }
3b07fd51 998 TVector2 epCorrected(dX, dY);
aa7eca65 999 dEventplaneAngleCorrected = epCorrected.Phi()/2.; // see AlEPSelectionTask.cxx:354
3b07fd51 1000 }
1001 else
1002 {
1003 dEventplaneAngleCorrected = dEventplaneAngle;
1004 }
1005
96ebdea7 1006 Double_t dFillEPCorrectionCheck[] = {dEventplaneAngle, dEventplaneAngleCorrected, dCentrality};
1007 fCorrelEventplaneDefaultCorrected->Fill(dFillEPCorrectionCheck);
1008
1009
4ffd6c6e 1010 dTrackDeltaphiPtEtaPhiCent[0] = RotatePhi(track->Phi(), dEventplaneAngleCorrected, fBinsPhi[fPhiNbins-1]);
1444967d 1011
4ffd6c6e 1012 dTrackDeltaphiPtEtaPhiCent[1] = track->Pt();
1013 dTrackDeltaphiPtEtaPhiCent[2] = track->Eta();
1014 dTrackDeltaphiPtEtaPhiCent[3] = track->Phi();
1015 dTrackDeltaphiPtEtaPhiCent[4] = dCentrality;
3dd0b8f4 1016
82a24e4d 1017
3dd0b8f4 1018 if( fIsMonteCarlo )
d25bcbe6 1019 {
3dd0b8f4 1020 mcPart = (AliAODMCParticle*)stack->At(TMath::Abs(track->GetLabel()));
1021 if( !mcPart ) { continue; }
1022
1023 // check for charge
1024 // if( !(IsMCTrackAccepted(mcPart)) ) { continue; }
1025
1026 bIsHijingParticle = IsHijingParticle(mcPart, genHijingHeader);
1027 // bIsPythiaParticle = IsPythiaParticle(mcPart, genPythiaHeader);
1028
1029 bIsPrimary = mcPart->IsPhysicalPrimary();
1030
1031 dMCTrackZvPtEtaCent[1] = mcPart->Pt();
1032 dMCTrackZvPtEtaCent[2] = mcPart->Eta();
1033 dMCTrackZvPtEtaCent[3] = dCentrality;
d25bcbe6 1034
4ffd6c6e 1035 dMCTrackDeltaphiPtEtaPhiCent[0] = RotatePhi(mcPart->Phi(), dEventplaneAngle, fBinsPhi[fPhiNbins-1]); // use eventplane and not reactionplan, similar to centrality vs impact paramter
1444967d 1036
4ffd6c6e 1037 dMCTrackDeltaphiPtEtaPhiCent[1] = mcPart->Pt();
1038 dMCTrackDeltaphiPtEtaPhiCent[2] = mcPart->Eta();
1039 dMCTrackDeltaphiPtEtaPhiCent[3] = mcPart->Phi();
1040 dMCTrackDeltaphiPtEtaPhiCent[4] = dCentrality;
3dd0b8f4 1041
1042 if(bIsPrimary && bIsHijingParticle)
1043 {
1044 fMCRecPrimZvPtEtaCent->Fill(dMCTrackZvPtEtaCent);
4ffd6c6e 1045 fMCRecPrimDeltaphiPtEtaPhiCent->Fill(dMCTrackDeltaphiPtEtaPhiCent);
3dd0b8f4 1046 fMCDCAPtPrimary->Fill(dDCAxyDCAzPt);
1047 }
1048
1049 if(!bIsPrimary /*&& !bIsHijingParticle*/)
1050 {
1051 Int_t indexMoth = mcPart->GetMother();
1052 if(indexMoth >= 0)
1053 {
1054 AliAODMCParticle* moth = (AliAODMCParticle*)stack->At(indexMoth);
1055 bMotherIsHijingParticle = IsHijingParticle(moth, genHijingHeader);
1056
1057 if(bMotherIsHijingParticle) // only store secondaries, which come from a not embedded signal!
1058 {
1059 fMCRecSecZvPtEtaCent->Fill(dMCTrackZvPtEtaCent);
4ffd6c6e 1060 fMCRecSecDeltaphiPtEtaPhiCent->Fill(dMCTrackDeltaphiPtEtaPhiCent);
3dd0b8f4 1061 fMCDCAPtSecondary->Fill(dDCAxyDCAzPt);
1062 // delete moth;
1063 }
1064 }
1065 }
1066 } // end isMonteCarlo
1067
1068 // ======================== fill histograms ========================
1069
1070 // only keep prim and sec from not embedded signal
1071 Bool_t bKeepMCTrack = kFALSE;
1072 if(fIsMonteCarlo)
1073 {
1074 if( (bIsHijingParticle && bIsPrimary) ^ (bMotherIsHijingParticle && !bIsPrimary) )
1075 {
1076 bKeepMCTrack = kTRUE;
1077 }
1078 else
d25bcbe6 1079 {
3dd0b8f4 1080 continue;
d25bcbe6 1081 }
3dd0b8f4 1082 }
1083
1084 bEventHasATrack = kTRUE;
1085
1086 fZvPtEtaCent->Fill(dTrackZvPtEtaCent);
4ffd6c6e 1087 fDeltaphiPtEtaPhiCent->Fill(dTrackDeltaphiPtEtaPhiCent);
3dd0b8f4 1088
1089 fDCAPtAccepted->Fill(dDCAxyDCAzPt);
1090
1091 if( (dTrackZvPtEtaCent[1] > GetCutPtMin()) &&
1092 (dTrackZvPtEtaCent[1] < GetCutPtMax()) &&
1093 (dTrackZvPtEtaCent[2] > GetCutEtaMin()) &&
1094 (dTrackZvPtEtaCent[2] < GetCutEtaMax()) )
1095 {
1096 iAcceptedMultiplicity++;
1097 bEventHasATrackInRange = kTRUE;
1098 fPt->Fill(track->Pt());
1099 fCharge->Fill(track->Charge());
ab1375ce 1100
1101 fPhiCent->Fill(track->Phi(), dCentrality);
1102 fPcosPhiCent->Fill(dCentrality, TMath::Cos(2.*track->Phi()));
1103 fPsinPhiCent->Fill(dCentrality, TMath::Sin(2.*track->Phi()));
82a24e4d 1104
1105 Double_t deltaphi = track->Phi() - dEventplaneAngleCorrected;
b7741813 1106 // if(deltaphi > TMath::Pi()) deltaphi -= 2.*TMath::Pi();
82a24e4d 1107
1108 fDeltaPhiCent->Fill(deltaphi, dCentrality);
3dd0b8f4 1109 }
d25bcbe6 1110 } // end track loop
1111
96ebdea7 1112 Int_t iContributorsQVector = ep->GetQContributionXArray()->GetSize();
1113 if(iContributorsQVector) fEventplaneSubtractedPercentage->Fill((Double_t)iSubtractedTracks/(Double_t)iContributorsQVector, dCentrality);
1114
ea3cfeda 1115 if(bEventHasATrack) { fEventStatistics->Fill("events with tracks",1); bEventHasATrack = kFALSE; }
d25bcbe6 1116
1117 if(bEventHasATrackInRange)
1118 {
3dd0b8f4 1119 fEventStatistics->Fill("events with tracks in range",1);
1120 fEventStatisticsCentrality->Fill(dCentrality);
1121
1122 bEventHasATrackInRange = kFALSE;
d25bcbe6 1123 }
1124
3dd0b8f4 1125 if(bIsEventSelectedMB) fEventStatisticsCentralityTrigger->Fill(dCentrality, 0);
1126 if(bIsEventSelectedSemi) fEventStatisticsCentralityTrigger->Fill(dCentrality, 1);
1127 if(bIsEventSelectedCentral) fEventStatisticsCentralityTrigger->Fill(dCentrality, 2);
1128
2942f542 1129 Double_t dEventZvMultCent[3] = {dEventZv, static_cast<Double_t>(iAcceptedMultiplicity), dCentrality};
ea3cfeda 1130 fZvMultCent->Fill(dEventZvMultCent);
d25bcbe6 1131
1444967d 1132 // store correlation between data and MC eventplane
1133 if(fIsMonteCarlo) fCorrelEventplaneMCDATA->Fill(dEventplaneAngle, dMCEventplaneAngle);
1134
d25bcbe6 1135 PostData(1, fOutputList);
1136
8a4ab847 1137 // delete pointers:
b7741813 1138 // delete [] iIndexAcceptedTracks;
d25bcbe6 1139}
1140
1444967d 1141
4ffd6c6e 1142Double_t AlidNdPtAnalysisPbPbAOD::RotatePhi(Double_t phiTrack, Double_t phiEP, Double_t dMaxDeltaPhi)
1444967d 1143{
1144 Double_t dPhi = 0;
4b6bf723 1145 dPhi = TMath::Abs(phiTrack - phiEP);
1444967d 1146
4ffd6c6e 1147// if( dPhi <= TMath::Pi() )
1148// {
1149// return dPhi;
1150// }
901067a9 1151 if( dPhi > TMath::Pi() )
1444967d 1152 {
aa7eca65 1153 dPhi = 2.*TMath::Pi() - dPhi;
4ffd6c6e 1154// return dPhi;
1155 }
1156
1157 if( dPhi > dMaxDeltaPhi)
1158 {
1159 dPhi = 2.*dMaxDeltaPhi - dPhi;
1444967d 1160 }
1161
ab1375ce 1162 // Printf("[E] dphi = %.4f , phiTrack = %.4f, phiEP = %.4f", dPhi, phiTrack, phiEP);
1444967d 1163
4ffd6c6e 1164// return -9999.;
1165
1166 return dPhi;
1444967d 1167}
1168
ea3cfeda
PL		 1169Bool_t AlidNdPtAnalysisPbPbAOD::SetRelativeCuts(AliAODEvent *event)
1170{
3dd0b8f4 1171 //
1172 // this function determines the absolute cut event-by-event based on the
1173 // the percentage given from outside
1174 // - cut set on Nclusters and NcrossedRows
1175 //
1176
ea3cfeda
PL		 1177 if(!event) return kFALSE;
1178
1179 AliAODTrack *tr = 0x0;
1180 TH1F *hCluster = new TH1F("hCluster","hCluster",160,0,160);
1181 TH1F *hCrossed = new TH1F("hCrossed","hCrossed",160,0,160);
1182
1183 for(Int_t itrack = 0; itrack < event->GetNumberOfTracks(); itrack++)
1184 {
f15c1f69 1185 tr = dynamic_cast<AliAODTrack*>(event->GetTrack(itrack));
1186 if(!tr) AliFatal("Not a standard AOD");
3dd0b8f4 1187 if(!tr) continue;
1188
1189 // do some selection already
1190 //if(!(tr->TestFilterBit(AliAODTrack::kTrkGlobal)) ) { continue; }
1191
1192 Double_t dNClustersTPC = tr->GetTPCNcls();
1193 Double_t dCrossedRowsTPC = tr->GetTPCClusterInfo(2,1);
1194
1195 hCluster->Fill(dNClustersTPC);
1196 hCrossed->Fill(dCrossedRowsTPC);
ea3cfeda
PL		 1197 }
1198
1199 // loop trough histogram to check, where percentage is reach
1200 Double_t dTotIntCluster = hCluster->Integral();
1201 Double_t dTotIntCrossed = hCrossed->Integral();
1202 Float_t dIntCluster = 0;
1203 Float_t dIntCrossed = 0;
1204
1205 if(dTotIntCluster)
1206 {
3dd0b8f4 1207 for(Int_t i = 0; i < hCluster->GetNbinsX(); i++)
1208 {
1209 if(hCluster->GetBinCenter(i) < 0) continue;
1210 dIntCluster += hCluster->GetBinContent(i);
1211 if(dIntCluster/dTotIntCluster > (1-GetCutPercMinNClustersTPC()))
1212 {
1213 SetCutMinNClustersTPC(hCluster->GetBinCenter(i));
1214 fCutPercClusters->Fill(hCluster->GetBinCenter(i));
1215 break;
1216 }
1217 }
ea3cfeda
PL		 1218 }
1219
1220 if(dTotIntCrossed)
1221 {
3dd0b8f4 1222 for(Int_t i = 0; i < hCrossed->GetNbinsX(); i++)
1223 {
1224 if(hCrossed->GetBinCenter(i) < 0) continue;
1225 dIntCrossed += hCrossed->GetBinContent(i);
1226 if(dIntCrossed/dTotIntCrossed > (1-GetCutPercMinNCrossedRowsTPC()))
1227 {
1228 SetCutMinNClustersTPC(hCrossed->GetBinCenter(i));
1229 fCutPercCrossed->Fill(hCrossed->GetBinCenter(i));
1230 break;
1231 }
1232 }
ea3cfeda
PL		 1233 }
1234
1235 delete hCrossed;
1236 delete hCluster;
1237 return kTRUE;
1238
1239}
1240
1241Bool_t AlidNdPtAnalysisPbPbAOD::IsTrackAccepted(AliAODTrack *tr, Double_t dCentrality, Double_t bMagZ)
d25bcbe6 1242{
3dd0b8f4 1243 //
1244 // this function checks the track parameters for quality
1245 // returns kTRUE if track is accepted
1246 //
1247 // - debug histograms (cuts vs pt,eta,phi) are filled in this function
1248 // - histogram for pt resolution correction are filled here as well
1249 //
1250
d25bcbe6 1251 if(!tr) return kFALSE;
b7741813 1252
d25bcbe6 1253 if(tr->Charge()==0) { return kFALSE; }
1254
ea3cfeda
PL		 1255 //
1256 // as done in AliAnalysisTaskFragmentationFunction
1257 //
1258
1259 Short_t sign = tr->Charge();
1260 Double_t xyz[50];
1261 Double_t pxpypz[50];
a0036e80 1262 Double_t cv[21];
ea3cfeda 1263
a0036e80 1264 for(Int_t i = 0; i < 21; i++) cv[i] = 0;
ea3cfeda
PL		 1265 for(Int_t i = 0; i < 50; i++) xyz[i] = 0;
1266 for(Int_t i = 0; i < 50; i++) pxpypz[i] = 0;
a0036e80 1267
ea3cfeda
PL		 1268 tr->GetXYZ(xyz);
1269 tr->GetPxPyPz(pxpypz);
a0036e80 1270 tr->GetCovarianceXYZPxPyPz(cv);
ea3cfeda
PL		 1271
1272 // similar error occured as this one:
1273 // See https://savannah.cern.ch/bugs/?102721
1274 // which is one of the two 11h re-filtering follow-ups:
1275 // Andrea Dainese now first does the beam pipe
1276 // check and then copies from the vtrack (was the other
1277 // way around) to avoid the crash in the etp::Set()
a0036e80 1278
1279 // if(xyz[0]*xyz[0]+xyz[1]*xyz[1] > 3.*3.) { return kFALSE; }
ea3cfeda 1280
8a4ab847 1281 AliExternalTrackParam par(xyz, pxpypz, cv, sign);
3dd0b8f4 1282 // AliExternalTrackParam *par = new AliExternalTrackParam(xyz, pxpypz, cv, sign); // high mem consumption!!!!
1283 static AliESDtrack dummy;
a0036e80 1284 // Double_t dLength = dummy.GetLengthInActiveZone(par,3,236, -5 ,0,0);
1285 // Double_t dLengthInTPC = GetLengthInTPC(tr, 1.8, 220, bMagZ);
1286
bc80e684 1287 Double_t dLengthInTPC = 0;
1288 if ( DoCutLengthInTPCPtDependent() ) { dLengthInTPC = dummy.GetLengthInActiveZone(&par,3,236, bMagZ ,0,0); }
b7741813 1289
d25bcbe6 1290 Double_t dNClustersTPC = tr->GetTPCNcls();
ba9a71a2 1291 Double_t dCrossedRowsTPC = tr->GetTPCNCrossedRows();//GetTPCClusterInfo(2,1);
ea3cfeda 1292 Double_t dFindableClustersTPC = tr->GetTPCNclsF();
9db7eb94 1293 Double_t dChi2PerClusterTPC = (dNClustersTPC>0)?tr->Chi2perNDF()*(dNClustersTPC-5)/dNClustersTPC:-1.; // see AliDielectronVarManager.h
a0036e80 1294 Double_t dOneOverPt = tr->OneOverPt();
8a4ab847 1295 Double_t dSigmaOneOverPt = TMath::Sqrt(par.GetSigma1Pt2());
9db7eb94 1296
1297 // hAllCrossedRowsTPC->Fill(dCrossedRowsTPC);
d25bcbe6 1298
bc80e684 1299 Double_t dCrossedRowsTPCOverFindableClustersTPC = 0;
1300 if(dFindableClustersTPC) dCrossedRowsTPCOverFindableClustersTPC = dCrossedRowsTPC/dFindableClustersTPC;
1301 Double_t dCheck[cqMax] = {dCrossedRowsTPC, dNClustersTPC, dChi2PerClusterTPC, dLengthInTPC, dCrossedRowsTPCOverFindableClustersTPC};// = new Double_t[cqMax];
ea3cfeda 1302 Double_t dKine[kqMax] = {tr->Pt(), tr->Eta(), tr->Phi()};// = new Double_t[kqMax];
bc80e684 1303
ea3cfeda
PL		 1304 // dKine[0] = tr->Pt();
1305 // dKine[1] = tr->Eta();
1306 // dKine[2] = tr->Phi();
1307 //
1308 // dCheck[0] = dCrossedRowsTPC;
1309 // dCheck[1] = dNClustersTPC;
1310 // dCheck[2] = dChi2PerClusterTPC;
72bb4ceb 1311
1312
1313 FillDebugHisto(dCheck, dKine, dCentrality, kFALSE);
1314
aa7eca65 1315 fCrossCheckPtresLength->Fill(dLengthInTPC, dSigmaOneOverPt*tr->Pt());
1316 fCrossCheckPtresRows->Fill(dCrossedRowsTPC, dSigmaOneOverPt*tr->Pt());
1317
1318
ea3cfeda
PL		 1319 // first cut on length
1320
1321 if( DoCutLengthInTPCPtDependent() && ( dLengthInTPC < GetPrefactorLengthInTPCPtDependent()*(130-5*TMath::Abs(1./tr->Pt())) ) ) { return kFALSE; }
95f26ffa 1322
9db7eb94 1323 // filter bit 5
a0036e80 1324 // if(!(tr->TestFilterBit(AliAODTrack::kTrkGlobal)) ) { return kFALSE; }
1325 if(!(tr->TestFilterBit(GetFilterBit())) ) { return kFALSE; }
95f26ffa 1326
9db7eb94 1327 // filter bit 4
1328 // if(!(tr->TestFilterBit(AliAODTrack::kTrkGlobalNoDCA)) ) { return kFALSE; }
1329
1330 // hFilterCrossedRowsTPC->Fill(dCrossedRowsTPC);
95f26ffa 1331
d0483ba3 1332
a0036e80 1333 if(dFindableClustersTPC == 0) {return kFALSE; }
1334 if(dCrossedRowsTPC < GetCutMinNCrossedRowsTPC()) { return kFALSE; }
bc80e684 1335 if( (dCrossedRowsTPCOverFindableClustersTPC) < GetCutMinRatioCrossedRowsOverFindableClustersTPC() ) { return kFALSE; }
a0036e80 1336 if(dNClustersTPC < GetCutMinNClustersTPC()) { return kFALSE; }
72bb4ceb 1337
a0036e80 1338 if (IsITSRefitRequired() && !(tr->GetStatus() & AliVTrack::kITSrefit)) { return kFALSE; } // no ITS refit
3dd0b8f4 1339
ab1375ce 1340 // do a relativ cut in Nclusters, both time at 80% of mean
1341 // if(fIsMonteCarlo)
1342 // {
1343 // if(dNClustersTPC < 88) { return kFALSE; }
1344 // }
1345 // else
1346 // {
1347 // if(dNClustersTPC < 76) { return kFALSE; }
1348 // }
1349
1350 // fill histogram for pT resolution correction
1351 Double_t dPtResolutionHisto[3] = { dOneOverPt, dSigmaOneOverPt, dCentrality };
1352 fPtResptCent->Fill(dPtResolutionHisto);
1353
1354 // fill debug histogram for all accepted tracks
1355 FillDebugHisto(dCheck, dKine, dCentrality, kTRUE);
1356
b7741813 1357 Double_t dFilterBitPhiCent[3] = {-10, -10, -10};
1358 if(tr->TestFilterBit(AliAODTrack::kTrkGlobal)) dFilterBitPhiCent[0] = 0;
1359 else if(tr->TestFilterBit(AliAODTrack::kTrkGlobalSDD)) dFilterBitPhiCent[0] = 1;
1360
1361 dFilterBitPhiCent[1] = tr->Phi();
1362 dFilterBitPhiCent[2] = dCentrality;
1363 fCrossCheckFilterBitPhiCent->Fill(dFilterBitPhiCent);
1364
ab1375ce 1365 // delete pointers
1366
1367 return kTRUE;
60134a87 1368}
1369
1370Bool_t AlidNdPtAnalysisPbPbAOD::FillDebugHisto(Double_t *dCrossCheckVar, Double_t *dKineVar, Double_t dCentrality, Bool_t bIsAccepted)
1371{
1372 if(bIsAccepted)
1373 {
1374 for(Int_t iCrossCheck = 0; iCrossCheck < cqMax; iCrossCheck++)
1375 {
1376 Double_t dFillIt[5] = {dCrossCheckVar[iCrossCheck], dKineVar[0], dKineVar[1], dKineVar[2], dCentrality};
1377 fCrossCheckAcc[iCrossCheck]->Fill(dFillIt);
1378 }
1379
1380 fCrossCheckRowsLengthAcc->Fill(dCrossCheckVar[cqLength], dCrossCheckVar[cqCrossedRows]);
1381 fCrossCheckClusterLengthAcc->Fill(dCrossCheckVar[cqLength], dCrossCheckVar[cqNcluster]);
1382 }
1383 else
1384 {
1385 for(Int_t iCrossCheck = 0; iCrossCheck < cqMax; iCrossCheck++)
1386 {
1387 Double_t dFillIt[5] = {dCrossCheckVar[iCrossCheck], dKineVar[0], dKineVar[1], dKineVar[2], dCentrality};
1388 fCrossCheckAll[iCrossCheck]->Fill(dFillIt);
1389 }
1390
1391 fCrossCheckRowsLength->Fill(dCrossCheckVar[cqLength], dCrossCheckVar[cqCrossedRows]);
1392 fCrossCheckClusterLength->Fill(dCrossCheckVar[cqLength], dCrossCheckVar[cqNcluster]);
1393 }
1394
1395 return kTRUE;
1396
1397}
1398
1399void AlidNdPtAnalysisPbPbAOD::StoreCutSettingsToHistogram()
1400{
1401 //
1402 // this function stores all cut settings to a histograms
1403 //
1404
1405 fCutSettings->Fill("IsMonteCarlo",fIsMonteCarlo);
1406
1407 fCutSettings->Fill("fCutMaxZVertex", fCutMaxZVertex);
1408
1409 // kinematic cuts
1410 fCutSettings->Fill("fCutPtMin", fCutPtMin);
1411 fCutSettings->Fill("fCutPtMax", fCutPtMax);
1412 fCutSettings->Fill("fCutEtaMin", fCutEtaMin);
1413 fCutSettings->Fill("fCutEtaMax", fCutEtaMax);
1414
1415 // track quality cut variables
1416 fCutSettings->Fill("fFilterBit", fFilterBit);
1417 if(fUseRelativeCuts) fCutSettings->Fill("fUseRelativeCuts", 1);
1418 if(fCutRequireTPCRefit) fCutSettings->Fill("fCutRequireTPCRefit", 1);
1419 if(fCutRequireITSRefit) fCutSettings->Fill("fCutRequireITSRefit", 1);
1420
1421 fCutSettings->Fill("fCutMinNumberOfClusters", fCutMinNumberOfClusters);
1422 fCutSettings->Fill("fCutPercMinNumberOfClusters", fCutPercMinNumberOfClusters);
1423 fCutSettings->Fill("fCutMinNumberOfCrossedRows", fCutMinNumberOfCrossedRows);
1424 fCutSettings->Fill("fCutPercMinNumberOfCrossedRows", fCutPercMinNumberOfCrossedRows);
1425
1426 fCutSettings->Fill("fCutMinRatioCrossedRowsOverFindableClustersTPC", fCutMinRatioCrossedRowsOverFindableClustersTPC);
1427 fCutSettings->Fill("fCutMaxFractionSharedTPCClusters", fCutMaxFractionSharedTPCClusters);
1428 fCutSettings->Fill("fCutMaxDCAToVertexXY", fCutMaxDCAToVertexXY);
1429 fCutSettings->Fill("fCutMaxChi2PerClusterITS", fCutMaxChi2PerClusterITS);
1430
1431 if(fCutDCAToVertex2D) fCutSettings->Fill("fCutDCAToVertex2D", 1);
1432 if(fCutRequireSigmaToVertex) fCutSettings->Fill("fCutRequireSigmaToVertex",1);
1433 fCutSettings->Fill("fCutMaxDCAToVertexXYPtDepPar0", fCutMaxDCAToVertexXYPtDepPar0);
1434 fCutSettings->Fill("fCutMaxDCAToVertexXYPtDepPar1", fCutMaxDCAToVertexXYPtDepPar1);
1435 fCutSettings->Fill("fCutMaxDCAToVertexXYPtDepPar2", fCutMaxDCAToVertexXYPtDepPar2);
1436
1437 if(fCutAcceptKinkDaughters) fCutSettings->Fill("fCutAcceptKinkDaughters", 1);
1438 fCutSettings->Fill("fCutMaxChi2TPCConstrainedGlobal", fCutMaxChi2TPCConstrainedGlobal);
1439 if(fCutLengthInTPCPtDependent) fCutSettings->Fill("fCutLengthInTPCPtDependent", 1);
1440 fCutSettings->Fill("fPrefactorLengthInTPCPtDependent", fPrefactorLengthInTPCPtDependent);
1441 fCutSettings->Fill(Form("EP selector %s", fEPselector.Data()), 1);
1442}
1443
1444Bool_t AlidNdPtAnalysisPbPbAOD::GetDCA(const AliAODTrack *track, AliAODEvent *evt, Double_t d0z0[2])
1445{
1446 // function adapted from AliDielectronVarManager.h
1447
1448 if(track->TestBit(AliAODTrack::kIsDCA)){
1449 d0z0[0]=track->DCA();
1450 d0z0[1]=track->ZAtDCA();
1451 return kTRUE;
1452 }
1453
1454 Bool_t ok=kFALSE;
1455 if(evt) {
1456 Double_t covd0z0[3];
1457 //AliAODTrack copy(*track);
1458 AliExternalTrackParam etp; etp.CopyFromVTrack(track);
1459
1460 Float_t xstart = etp.GetX();
1461 if(xstart>3.) {
1462 d0z0[0]=-999.;
1463 d0z0[1]=-999.;
1464 //printf("This method can be used only for propagation inside the beam pipe \n");
1465 return kFALSE;
1466 }
1467
1468
1469 AliAODVertex *vtx =(AliAODVertex*)(evt->GetPrimaryVertex());
1470 Double_t fBzkG = evt->GetMagneticField(); // z componenent of field in kG
1471 ok = etp.PropagateToDCA(vtx,fBzkG,kVeryBig,d0z0,covd0z0);
1472 //ok = copy.PropagateToDCA(vtx,fBzkG,kVeryBig,d0z0,covd0z0);
1473 }
1474 if(!ok){
1475 d0z0[0]=-999.;
1476 d0z0[1]=-999.;
1477 }
1478 return ok;
1479}
1480
1481
1482Bool_t AlidNdPtAnalysisPbPbAOD::IsMCTrackAccepted(AliAODMCParticle *part)
1483{
1484 if(!part) return kFALSE;
1485
1486 Double_t charge = part->Charge()/3.;
1487 if (TMath::Abs(charge) < 0.001) return kFALSE;
1488
1489 return kTRUE;
1490}
1491
1492const char * AlidNdPtAnalysisPbPbAOD::GetParticleName(Int_t pdg)
1493{
1494 TParticlePDG * p1 = TDatabasePDG::Instance()->GetParticle(pdg);
1495 if(p1) return p1->GetName();
1496 return Form("%d", pdg);
1497}
1498
1499AliGenHijingEventHeader* AlidNdPtAnalysisPbPbAOD::GetHijingEventHeader(AliAODMCHeader *header)
1500{
1501 //
1502 // inspired by PWGJE/AliPWG4HighPtSpectra.cxx
1503 //
1504
1505 if(!header) return 0x0;
1506 AliGenHijingEventHeader* hijingGenHeader = NULL;
1507
1508 TList* headerList = header->GetCocktailHeaders();
1509
1510 for(Int_t i = 0; i < headerList->GetEntries(); i++)
1511 {
1512 hijingGenHeader = dynamic_cast<AliGenHijingEventHeader*>(headerList->At(i));
1513 if(hijingGenHeader) break;
1514 }
1515
1516 if(!hijingGenHeader) return 0x0;
1517
1518 return hijingGenHeader;
1519}
1520
1521AliGenPythiaEventHeader* AlidNdPtAnalysisPbPbAOD::GetPythiaEventHeader(AliAODMCHeader *header)
1522{
1523 //
1524 // inspired by PWGJE/AliPWG4HighPtSpectra.cxx
1525 //
1526
1527 if(!header) return 0x0;
1528 AliGenPythiaEventHeader* PythiaGenHeader = NULL;
1529
1530 TList* headerList = header->GetCocktailHeaders();
1531
1532 for(Int_t i = 0; i < headerList->GetEntries(); i++)
1533 {
1534 PythiaGenHeader = dynamic_cast<AliGenPythiaEventHeader*>(headerList->At(i));
1535 if(PythiaGenHeader) break;
1536 }
1537
1538 if(!PythiaGenHeader) return 0x0;
1539
1540 return PythiaGenHeader;
1541}
1542
1543//________________________________________________________________________
1544Bool_t AlidNdPtAnalysisPbPbAOD::IsHijingParticle(const AliAODMCParticle *part, AliGenHijingEventHeader* hijingGenHeader){
1545
1546 // Check whether a particle is from Hijing or some injected
1547 // returns kFALSE if particle is injected
1548
1549 if(part->Label() > (hijingGenHeader->NProduced()-1)) return kFALSE;
1550 return kTRUE;
1551}
1552
1553//________________________________________________________________________
1554Bool_t AlidNdPtAnalysisPbPbAOD::IsPythiaParticle(const AliAODMCParticle *part, AliGenPythiaEventHeader* pythiaGenHeader){
1555
1556 // Check whether a particle is from Pythia or some injected
1557
1558 if(part->Label() > (pythiaGenHeader->NProduced()-1)) return kFALSE;
1559 return kTRUE;
1560}
1561
1562Double_t* AlidNdPtAnalysisPbPbAOD::GetArrayClone(Int_t n, Double_t* source)
1563{
1564 if (!source || n==0) return 0;
1565 Double_t* dest = new Double_t[n];
1566 for (Int_t i=0; i<n ; i++) { dest[i] = source[i]; }
1567 return dest;
1568}
1569
1570void AlidNdPtAnalysisPbPbAOD::Terminate(Option_t *)
1571{
1572
1573}
1574
1575