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