]>
Commit | Line | Data |
---|---|---|
70da6c5a | 1 | /************************************************************************* |
2 | * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. * | |
3 | * * | |
4 | * Author: The ALICE Off-line Project. * | |
5 | * Contributors are mentioned in the code where appropriate. * | |
6 | * * | |
7 | * Permission to use, copy, modify and distribute this software and its * | |
8 | * documentation strictly for non-commercial purposes is hereby granted * | |
9 | * without fee, provided that the above copyright notice appears in all * | |
10 | * copies and that both the copyright notice and this permission notice * | |
11 | * appear in the supporting documentation. The authors make no claims * | |
12 | * about the suitability of this software for any purpose. It is * | |
13 | * provided "as is" without express or implied warranty. * | |
14 | **************************************************************************/ | |
15 | // | |
16 | // Class for impact parameter (DCA) of charged particles | |
17 | // Study resolution and pull: prepare for beauty study | |
18 | // | |
19 | // Authors: | |
20 | // Hongyan Yang <hongyan@physi.uni-heidelberg.de> | |
21 | // Carlo Bombonati <carlo.bombonati@cern.ch> | |
22 | // | |
23 | ||
24 | #include "TMath.h" | |
25 | #include "TH1F.h" | |
26 | #include "TList.h" | |
27 | #include <TParticle.h> | |
28 | #include "AliMCParticle.h" | |
29 | #include "AliESDtrack.h" | |
30 | #include "AliESDEvent.h" | |
31 | #include "AliMCEvent.h" | |
faee3b18 | 32 | #include "AliMCVertex.h" |
33 | ||
34 | #include "AliKFParticle.h" | |
35 | #include "AliKFVertex.h" | |
36 | ||
37 | #include "AliESDVertex.h" | |
38 | ||
39 | #include "AliPID.h" | |
70da6c5a | 40 | |
41 | #include "AliHFEdca.h" | |
42 | ||
faee3b18 | 43 | |
70da6c5a | 44 | ClassImp(AliHFEdca) |
45 | ||
46 | //________________________________________________________________________ | |
47 | const Char_t* AliHFEdca::fgkParticles[12] = { | |
faee3b18 | 48 | // particles name |
49 | "electron", "muonMinus","pionMinus", "kaonMinus", "protonMinus", | |
50 | "positron", "muonPlus", "pionPlus", "kaonPlus", "protonPlus", | |
51 | "allNegative", "allPositive" | |
70da6c5a | 52 | }; |
faee3b18 | 53 | |
54 | const Int_t AliHFEdca::fgkPdgParticle[10] = { | |
55 | // 11, 13, -211, -233, -2122, | |
56 | // -11, -13, 211, 233, 2122}; | |
57 | kPDGelectron, kPDGmuon, -kPDGpion, -kPDGkaon, -kPDGproton, | |
58 | -kPDGelectron, -kPDGmuon, kPDGpion, kPDGkaon, kPDGproton}; | |
59 | ||
70da6c5a | 60 | //________________________________________________________________________ |
61 | const Int_t AliHFEdca::fgkColorPart[12] = { | |
faee3b18 | 62 | // colors assigned to particles |
63 | kRed, kBlue, kGreen+2, kYellow+2, kMagenta, | |
64 | kRed+2, kBlue+2, kGreen+4, kYellow+4, kMagenta+2, | |
65 | kBlack, kGray+1 | |
70da6c5a | 66 | }; |
67 | ||
70da6c5a | 68 | //________________________________________________________________________ |
faee3b18 | 69 | const Float_t AliHFEdca::fgkPtIntv[51] = { |
70 | // define pT bins | |
71 | 0.10, 0.15, 0.20, 0.25, 0.30, 0.35, 0.40, 0.45, 0.50, 0.55, | |
72 | 0.60, 0.70, 0.80, 0.90, 1.00, 1.10, 1.20, 1.30, 1.40, 1.50, | |
73 | 1.70, 1.90, 2.10, 2.30, 2.50, 2.70, 2.90, 3.10, 3.30, 3.50, | |
74 | 3.80, 4.10, 4.40, 4.70, 5.00, 5.30, 5.60, 5.90, 6.20, 6.50, | |
75 | 7.00, 7.50, 8.00, 9.00, 10.0, 11.0, 12.0, 13.0, 15.0, 18.0, | |
76 | 20.00}; | |
70da6c5a | 77 | |
78 | //________________________________________________________________________ | |
79 | const Char_t *AliHFEdca::fgkDcaVar[2] = { | |
faee3b18 | 80 | "DcaXY", "DcaZ"}; |
70da6c5a | 81 | |
82 | //________________________________________________________________________ | |
83 | const Char_t *AliHFEdca::fgkDcaVarTitle[2] ={ | |
faee3b18 | 84 | ";dca_{xy} [#mum];counts", ";dca_{z} [#mum];counts"}; |
85 | ||
86 | //________________________________________________________________________ | |
87 | const Char_t *AliHFEdca::fgkVertexVar[3] = { | |
88 | "VertexX", "VertexY", "VertexZ"}; | |
89 | ||
90 | //________________________________________________________________________ | |
91 | const Char_t *AliHFEdca::fgkVertexVarTitle[3] ={ | |
92 | ";vertex_{x} [#mum];counts", ";vertex_{y} [#mum];counts", ";vertex_{z} [#mum];counts"}; | |
93 | ||
94 | //________________________________________________________________________ | |
95 | const Char_t *AliHFEdca::fgkResDcaVar[2] = { | |
96 | "deltaDcaXY", "deltaDcaZ"}; | |
70da6c5a | 97 | |
faee3b18 | 98 | //________________________________________________________________________ |
99 | const Char_t *AliHFEdca::fgkResDcaVarTitle[2] ={ | |
100 | ";residual #Delta(d_{xy}) [#mum];counts", ";residual #Delta(d_{z}) [#mum];counts"}; | |
70da6c5a | 101 | |
102 | //________________________________________________________________________ | |
103 | const Char_t *AliHFEdca::fgkPullDcaVar[2] = { | |
faee3b18 | 104 | "pullDcaXY", "pullDcaZ" |
70da6c5a | 105 | }; |
106 | ||
107 | //________________________________________________________________________ | |
108 | const Char_t *AliHFEdca::fgkPullDcaVarTitle[2] = { | |
faee3b18 | 109 | ";residual dca_{xy}/(error dca_{xy});counts", |
110 | ";residual dca_{z}/(error dca_{z});counts" | |
70da6c5a | 111 | }; |
112 | ||
113 | //________________________________________________________________________ | |
faee3b18 | 114 | const Char_t *AliHFEdca::fgkPullDataDcaVarTitle[2] = { |
115 | ";dca_{xy}^{data}/error dca_{xy};counts", | |
116 | ";dca_{z}^{data}/error dca_{z};counts" | |
117 | }; | |
70da6c5a | 118 | |
faee3b18 | 119 | //________________________________________________________________________ |
120 | AliHFEdca::AliHFEdca(): | |
121 | TObject(), | |
122 | fStat(NULL) | |
70da6c5a | 123 | { |
faee3b18 | 124 | // default constructor |
125 | ||
69ac0e6f | 126 | for(Int_t j=0; j<kNParticles; j++){ |
127 | fHistMcPid[j] = new TH1F(); | |
128 | fHistEsdPid[j] = new TH1F(); | |
129 | fHistDataEsdPid[j] = new TH1F(); | |
130 | } | |
131 | ||
132 | for(Int_t i=0; i<3; i++){ | |
133 | fHistMCvertex[i] = new TH1F(); | |
134 | fHistESDvertex[i] = new TH1F(); | |
135 | fHistDatavertex[i] = new TH1F(); | |
136 | } | |
137 | ||
138 | for(Int_t iEle=0; iEle<2; iEle++){ | |
139 | fHistDataHfePid[iEle] = new TH1F(); | |
140 | } | |
141 | ||
70da6c5a | 142 | } |
143 | ||
144 | //________________________________________________________________________ | |
faee3b18 | 145 | AliHFEdca::AliHFEdca(const AliHFEdca &ref): |
146 | TObject(ref), | |
147 | fStat(ref.fStat) | |
70da6c5a | 148 | { |
faee3b18 | 149 | // copy constructor |
70da6c5a | 150 | |
69ac0e6f | 151 | for(Int_t j=0; j<kNParticles; j++){ |
152 | fHistMcPid[j] = ref.fHistMcPid[j]; | |
153 | fHistEsdPid[j] = ref.fHistEsdPid[j]; | |
154 | fHistDataEsdPid[j] = ref.fHistDataEsdPid[j]; | |
155 | } | |
156 | ||
157 | for(Int_t i=0; i<3; i++){ | |
158 | fHistMCvertex[i] = ref.fHistMCvertex[i]; | |
159 | fHistESDvertex[i] = ref.fHistESDvertex[i]; | |
160 | fHistDatavertex[i] = ref.fHistDatavertex[i]; | |
161 | } | |
162 | ||
163 | for(Int_t iEle=0; iEle<2; iEle++){ | |
164 | fHistDataHfePid[iEle] = ref.fHistDataHfePid[iEle]; | |
165 | } | |
166 | ||
70da6c5a | 167 | } |
168 | //_______________________________________________________________________________________________ | |
faee3b18 | 169 | AliHFEdca&AliHFEdca::operator=(const AliHFEdca &ref) |
70da6c5a | 170 | { |
faee3b18 | 171 | // |
172 | // Assignment operator | |
173 | // | |
174 | ||
faee3b18 | 175 | if(this == &ref) return *this; |
176 | TObject::operator=(ref); | |
177 | return *this; | |
69ac0e6f | 178 | |
70da6c5a | 179 | } |
180 | ||
181 | //________________________________________________________________________ | |
182 | AliHFEdca::~AliHFEdca() | |
183 | { | |
faee3b18 | 184 | // default destructor |
185 | ||
186 | for(Int_t j=0; j<kNParticles; j++){ | |
187 | for(Int_t i=0; i<kNPtBins; i++){ | |
188 | if(fHistDcaXYRes[j][i]) delete fHistDcaXYRes[j][i]; | |
189 | if(fHistDcaZRes[j][i]) delete fHistDcaZRes[j][i]; | |
190 | ||
191 | if(fHistDcaXYPull[j][i]) delete fHistDcaXYPull[j][i]; | |
192 | if(fHistDcaZPull[j][i]) delete fHistDcaZPull[j][i]; | |
193 | ||
194 | if(fHistDcaXY[j][i]) delete fHistDcaXY[j][i]; | |
195 | if(fHistDcaZ[j][i]) delete fHistDcaZ[j][i]; | |
70da6c5a | 196 | |
faee3b18 | 197 | if(j<(kNParticles-2)){ |
198 | if(fHistEPDcaXYRes[j][i]) delete fHistEPDcaXYRes[j][i]; | |
199 | if(fHistEPDcaZRes[j][i]) delete fHistEPDcaZRes[j][i]; | |
200 | ||
201 | if(fHistEPDcaXYPull[j][i]) delete fHistEPDcaXYPull[j][i]; | |
202 | if(fHistEPDcaZPull[j][i]) delete fHistEPDcaZPull[j][i]; | |
203 | ||
204 | if(fHistEPDcaXY[j][i]) delete fHistEPDcaXY[j][i]; | |
205 | if(fHistEPDcaZ[j][i]) delete fHistEPDcaZ[j][i]; | |
206 | } | |
207 | ||
208 | if(fHistKFDcaXY[j][i]) delete fHistKFDcaXY[j][i]; | |
209 | if(fHistKFDcaZ[j][i]) delete fHistKFDcaZ[j][i]; | |
210 | ||
211 | if(fHistDataDcaXY[j][i]) delete fHistDataDcaXY[j][i]; | |
212 | if(fHistDataDcaZ[j][i]) delete fHistDataDcaZ[j][i]; | |
213 | if(fHistDataWoDcaXY[j][i]) delete fHistDataWoDcaXY[j][i]; | |
214 | if(fHistDataWoDcaZ[j][i]) delete fHistDataWoDcaZ[j][i]; | |
215 | ||
216 | if(fHistDataDcaXYPull[j][i]) delete fHistDataDcaXYPull[j][i]; | |
217 | if(fHistDataDcaZPull[j][i]) delete fHistDataDcaZPull[j][i]; | |
218 | if(fHistDataWoDcaXYPull[j][i]) delete fHistDataWoDcaXYPull[j][i]; | |
219 | if(fHistDataWoDcaZPull[j][i]) delete fHistDataWoDcaZPull[j][i]; | |
220 | } | |
221 | ||
222 | if(fHistMcPid[j]) delete fHistMcPid[j]; | |
223 | if(fHistEsdPid[j]) delete fHistEsdPid[j]; | |
224 | if(fHistDataEsdPid[j]) delete fHistDataEsdPid[j]; | |
225 | } | |
226 | ||
227 | for(Int_t i=0; i<3; i++){ | |
228 | if(fHistMCvertex[i]) delete fHistMCvertex[i]; | |
229 | if(fHistESDvertex[i]) delete fHistESDvertex[i]; | |
230 | if(fHistDatavertex[i]) delete fHistDatavertex[i]; | |
231 | } | |
232 | ||
233 | // for the HFEpid | |
234 | for(Int_t iEle=0; iEle<2; iEle++){ | |
235 | for(Int_t iPt=0; iPt<kNPtBins; iPt++){ | |
236 | if(fHistHPDcaXYRes[iEle][iPt]) delete fHistHPDcaXYRes[iEle][iPt]; | |
237 | if(fHistHPDcaZRes[iEle][iPt]) delete fHistHPDcaZRes[iEle][iPt]; | |
238 | if(fHistHPDcaXYPull[iEle][iPt]) delete fHistHPDcaXYPull[iEle][iPt]; | |
239 | if(fHistHPDcaZPull[iEle][iPt]) delete fHistHPDcaZPull[iEle][iPt]; | |
240 | if(fHistHPDcaXY[iEle][iPt]) delete fHistHPDcaXY[iEle][iPt]; | |
241 | if(fHistHPDcaZ[iEle][iPt]) delete fHistHPDcaZ[iEle][iPt]; | |
69ac0e6f | 242 | |
243 | ||
244 | // Data | |
faee3b18 | 245 | if(fHistHPDataDcaXY[iEle][iPt]) delete fHistHPDataDcaXY[iEle][iPt]; |
246 | if(fHistHPDataDcaZ[iEle][iPt]) delete fHistHPDataDcaZ[iEle][iPt]; | |
247 | if(fHistHPDataDcaXYPull[iEle][iPt]) delete fHistHPDataDcaXYPull[iEle][iPt]; | |
248 | if(fHistHPDataDcaZPull[iEle][iPt]) delete fHistHPDataDcaZPull[iEle][iPt]; | |
249 | ||
250 | } | |
251 | for(Int_t i=0; i<2; i++) | |
252 | if(fHistHfePid[iEle][i]) delete fHistHfePid[iEle][i]; | |
253 | ||
254 | if(fHistDataHfePid[iEle]) delete fHistDataHfePid[iEle]; | |
255 | ||
256 | } | |
70da6c5a | 257 | |
faee3b18 | 258 | if(fStat) delete fStat; |
70da6c5a | 259 | |
6555e2ad | 260 | //Printf("analysis done\n"); |
69ac0e6f | 261 | |
70da6c5a | 262 | } |
263 | ||
264 | //________________________________________________________________________ | |
c2690925 | 265 | void AliHFEdca::InitAnalysis()const{ |
faee3b18 | 266 | |
6555e2ad | 267 | //Printf("initialize analysis\n"); |
70da6c5a | 268 | |
269 | } | |
270 | ||
271 | ||
272 | //________________________________________________________________________ | |
273 | void AliHFEdca::PostAnalysis() const | |
274 | { | |
faee3b18 | 275 | // do fit |
276 | // moved to dcaPostAnalysis.C | |
70da6c5a | 277 | |
278 | } | |
faee3b18 | 279 | |
280 | ||
70da6c5a | 281 | //________________________________________________________________________ |
282 | void AliHFEdca::CreateHistogramsResidual(TList *residualList){ | |
faee3b18 | 283 | // define histogram |
284 | // 1. residual | |
285 | ||
286 | // for residuals | |
69ac0e6f | 287 | // fHistDcaXYRes[kNParticles][kNPtBins]=0x0; |
288 | // fHistDcaZRes[kNParticles][kNPtBins]=0x0; | |
289 | // fHistEPDcaXYRes[kNParticles-2][kNPtBins]=0x0; | |
290 | // fHistEPDcaZRes[kNParticles-2][kNPtBins]=0x0; | |
faee3b18 | 291 | |
292 | const Int_t nBins = 1000; | |
293 | const Float_t maxXYBin = 1000.; | |
294 | const Float_t maxZBin = 1000.; | |
295 | ||
296 | ||
297 | for(Int_t k=0; k<kNDcaVar; k++){ | |
298 | TString histTitle((const char*)fgkResDcaVarTitle[k]); | |
299 | ||
300 | for(Int_t j=0; j<kNParticles; j++){ | |
301 | for(Int_t i=0; i<kNPtBins; i++){ | |
70da6c5a | 302 | |
303 | TString histName((const char*)fgkParticles[j]); | |
faee3b18 | 304 | histName += Form("_MCpid_%s_pT-%.2f-%.2f", (const char*)fgkResDcaVar[k], fgkPtIntv[i], fgkPtIntv[i+1]); |
305 | TString histEPName((const char*)fgkParticles[j]); | |
306 | histEPName += Form("_ESDpid_%s_pT-%.2f-%.2f", (const char*)fgkResDcaVar[k], fgkPtIntv[i], fgkPtIntv[i+1]); | |
70da6c5a | 307 | |
308 | if(k==0){ | |
309 | fHistDcaXYRes[j][i] = new TH1F((const char*)histName, (const char*)histTitle, nBins, -maxXYBin, maxXYBin); | |
6555e2ad | 310 | fHistDcaXYRes[j][i]->SetLineColor((int)fgkColorPart[j]); |
faee3b18 | 311 | if(j<(kNParticles-2)){ |
312 | fHistEPDcaXYRes[j][i] = new TH1F((const char*)histEPName, (const char*)histTitle, nBins, -maxXYBin, maxXYBin); | |
69ac0e6f | 313 | fHistEPDcaXYRes[j][i]->SetLineColor((int)fgkColorPart[j]);} |
70da6c5a | 314 | } |
315 | if(k==1){ | |
316 | fHistDcaZRes[j][i] = new TH1F((const char*)histName, (const char*)histTitle, nBins, -maxZBin, maxZBin); | |
69ac0e6f | 317 | fHistDcaZRes[j][i]->SetLineColor((int)fgkColorPart[j]); |
faee3b18 | 318 | if(j<(kNParticles-2)){ |
319 | fHistEPDcaZRes[j][i] = new TH1F((const char*)histEPName, (const char*)histTitle, nBins, -maxZBin, maxZBin); | |
69ac0e6f | 320 | fHistEPDcaZRes[j][i]->SetLineColor((int)fgkColorPart[j]); } |
70da6c5a | 321 | } |
faee3b18 | 322 | } // 50 pt bins |
323 | } //12 nparticles | |
324 | } // 2 dca var | |
325 | ||
326 | // TList *fResidualList = 0; | |
327 | residualList->SetOwner(); | |
328 | residualList->SetName("residual"); | |
329 | for(Int_t iPart=0; iPart<kNParticles; iPart++){ | |
330 | for(Int_t iPtBin=0; iPtBin<kNPtBins; iPtBin++){ | |
331 | residualList->Add(fHistDcaXYRes[iPart][iPtBin]); | |
332 | residualList->Add(fHistDcaZRes[iPart][iPtBin]); | |
333 | if(iPart<(kNParticles-2)){ | |
334 | residualList->Add(fHistEPDcaXYRes[iPart][iPtBin]); | |
335 | residualList->Add(fHistEPDcaZRes[iPart][iPtBin]); | |
336 | } | |
337 | } // loop over pt bins | |
338 | } // loop over particles (pos, neg) | |
339 | ||
340 | ||
70da6c5a | 341 | |
342 | } | |
343 | ||
344 | ||
345 | //________________________________________________________________________ | |
346 | void AliHFEdca::CreateHistogramsPull(TList *pullList){ | |
faee3b18 | 347 | // define histogram |
348 | // 2. pull | |
349 | ||
350 | const Int_t nBins = 1000; | |
351 | const Float_t maxXYBin = 20.; | |
352 | const Float_t maxZBin = 20.; | |
353 | ||
354 | ||
355 | // for pull ----------------------------------------------------------------------- | |
69ac0e6f | 356 | // fHistDcaXYPull[kNParticles][kNPtBins]=0x0; |
357 | // fHistDcaZPull[kNParticles][kNPtBins]=0x0; | |
358 | // fHistEPDcaXYPull[kNParticles-2][kNPtBins]=0x0; | |
359 | // fHistEPDcaZPull[kNParticles-2][kNPtBins]=0x0; | |
faee3b18 | 360 | |
361 | ||
362 | for(Int_t k=0; k<kNDcaVar; k++){ | |
363 | TString histTitle((const char*)fgkPullDcaVarTitle[k]); | |
364 | ||
365 | for(Int_t j=0; j<kNParticles; j++){ | |
366 | for(Int_t i=0; i<kNPtBins; i++){ | |
70da6c5a | 367 | |
faee3b18 | 368 | TString histName((const char*)fgkParticles[j]); |
369 | histName += Form("_MCpid_%s_pT-%.2f-%.2f", (const char*)fgkPullDcaVar[k], fgkPtIntv[i], fgkPtIntv[i+1]); | |
370 | TString histEPName((const char*)fgkParticles[j]); | |
371 | histEPName += Form("_ESDpid_%s_pT-%.2f-%.2f", (const char*)fgkPullDcaVar[k], fgkPtIntv[i], fgkPtIntv[i+1]); | |
70da6c5a | 372 | |
373 | if(k==0){ | |
374 | fHistDcaXYPull[j][i] = new TH1F((const char*)histName, (const char*)histTitle, nBins, 1-maxXYBin, 1+maxXYBin); | |
69ac0e6f | 375 | fHistDcaXYPull[j][i]->SetLineColor((int)fgkColorPart[j]); |
faee3b18 | 376 | if(j<(kNParticles-2)) { |
377 | fHistEPDcaXYPull[j][i] = new TH1F((const char*)histEPName, (const char*)histTitle, nBins, 1-maxXYBin, 1+maxXYBin); | |
69ac0e6f | 378 | fHistEPDcaXYPull[j][i]->SetLineColor((int)fgkColorPart[j]);} |
70da6c5a | 379 | } |
380 | if(k==1){ | |
381 | fHistDcaZPull[j][i] = new TH1F((const char*)histName, (const char*)histTitle, nBins, 1-maxZBin, 1+maxZBin); | |
69ac0e6f | 382 | fHistDcaZPull[j][i]->SetLineColor((int)fgkColorPart[j]); |
faee3b18 | 383 | if(j<(kNParticles-2)) { |
384 | fHistEPDcaZPull[j][i] = new TH1F((const char*)histEPName, (const char*)histTitle, nBins, 1-maxZBin, 1+maxZBin); | |
69ac0e6f | 385 | fHistEPDcaZPull[j][i]->SetLineColor((int)fgkColorPart[j]);} |
faee3b18 | 386 | } |
387 | } // 50 pt bins | |
388 | } //6 nparticles | |
389 | } // 2 dca var | |
390 | ||
391 | // TList *fPullList = 0; | |
392 | pullList->SetOwner(); | |
393 | pullList->SetName("pull"); | |
394 | for(Int_t iPart=0; iPart<kNParticles; iPart++){ | |
395 | for(Int_t iPtBin=0; iPtBin<kNPtBins; iPtBin++){ | |
396 | pullList->Add(fHistDcaXYPull[iPart][iPtBin]); | |
397 | pullList->Add(fHistDcaZPull[iPart][iPtBin]); | |
398 | if(iPart<(kNParticles-2)){ | |
399 | pullList->Add(fHistDcaXYPull[iPart][iPtBin]); | |
400 | pullList->Add(fHistDcaZPull[iPart][iPtBin]); } | |
401 | } // loop over pt bins | |
402 | } // loop over particles (pos, neg) | |
403 | ||
404 | } | |
405 | ||
406 | ||
407 | //________________________________________________________________________ | |
408 | void AliHFEdca::CreateHistogramsDca(TList *dcaList){ | |
409 | // | |
410 | // define histograms: MC dca | |
411 | // | |
412 | ||
413 | // statistics | |
414 | fStat = 0x0; | |
415 | fStat = new TH1I("fStatistics", "allStatistics;ID;counts", 7, -3.5, 3.5); | |
416 | fStat->SetMarkerStyle(20); | |
417 | fStat->SetMarkerColor(3); | |
418 | fStat->SetMarkerSize(1); | |
419 | ||
420 | // for dca | |
69ac0e6f | 421 | // fHistDcaXY[kNParticles][kNPtBins]=0x0; |
422 | // fHistDcaZ[kNParticles][kNPtBins]=0x0; | |
423 | // fHistEPDcaXY[kNParticles-2][kNPtBins]=0x0; | |
424 | // fHistEPDcaZ[kNParticles-2][kNPtBins]=0x0; | |
faee3b18 | 425 | |
426 | const Int_t nBins = 1000; | |
427 | const Float_t maxXYBin = 1000.; | |
428 | const Float_t maxZBin = 1000.; | |
429 | ||
430 | ||
431 | for(Int_t k=0; k<kNDcaVar; k++){ | |
432 | TString histTitle((const char*)fgkDcaVarTitle[k]); | |
433 | ||
434 | for(Int_t j=0; j<kNParticles; j++){ | |
435 | for(Int_t i=0; i<kNPtBins; i++){ | |
436 | ||
437 | TString histName((const char*)fgkParticles[j]); | |
438 | histName += Form("_MCpid_%s_pT-%.2f-%.2f", (const char*)fgkDcaVar[k], fgkPtIntv[i], fgkPtIntv[i+1]); | |
439 | ||
440 | TString histNameEP((const char*)fgkParticles[j]); | |
441 | histNameEP += Form("_ESDpid_%s_pT-%.2f-%.2f", (const char*)fgkDcaVar[k], fgkPtIntv[i], fgkPtIntv[i+1]); | |
442 | ||
443 | if(k==0){ | |
444 | fHistDcaXY[j][i] = new TH1F((const char*)histName, (const char*)histTitle, nBins, -maxXYBin, maxXYBin); | |
69ac0e6f | 445 | fHistDcaXY[j][i]->SetLineColor((int)fgkColorPart[j]); |
faee3b18 | 446 | |
447 | if(j<(kNParticles-2)){ | |
448 | fHistEPDcaXY[j][i] = new TH1F((const char*)histNameEP, (const char*)histTitle, nBins, -maxXYBin, maxXYBin); | |
69ac0e6f | 449 | fHistEPDcaXY[j][i]->SetLineColor((int)fgkColorPart[j]);} |
faee3b18 | 450 | } |
451 | if(k==1){ | |
452 | fHistDcaZ[j][i] = new TH1F((const char*)histName, (const char*)histTitle, nBins, -maxZBin, maxZBin); | |
69ac0e6f | 453 | fHistDcaZ[j][i]->SetLineColor((int)fgkColorPart[j]); |
faee3b18 | 454 | if(j<(kNParticles-2)){ |
455 | fHistEPDcaZ[j][i] = new TH1F((const char*)histNameEP, (const char*)histTitle, nBins, -maxZBin, maxZBin); | |
69ac0e6f | 456 | fHistEPDcaZ[j][i]->SetLineColor((int)fgkColorPart[j]);} |
faee3b18 | 457 | } |
458 | } // 50 pt bins | |
459 | } //12 nparticles | |
460 | } // 2 dca var | |
461 | ||
462 | // TList *fDcaList = 0; | |
463 | dcaList->SetOwner(); | |
464 | dcaList->SetName("mcDcaDistr"); | |
465 | for(Int_t iPart=0; iPart<kNParticles; iPart++){ | |
466 | for(Int_t iPtBin=0; iPtBin<kNPtBins; iPtBin++){ | |
467 | dcaList->Add(fHistDcaXY[iPart][iPtBin]); | |
468 | dcaList->Add(fHistDcaZ[iPart][iPtBin]); | |
469 | if(iPart<(kNParticles-2)) { | |
470 | dcaList->Add(fHistEPDcaXY[iPart][iPtBin]); | |
471 | dcaList->Add(fHistEPDcaZ[iPart][iPtBin]); } | |
472 | } // loop over pt bins | |
473 | } // loop over particles (pos, neg) | |
474 | ||
475 | dcaList->Add(fStat); | |
476 | ||
477 | } | |
478 | ||
479 | //________________________________________________________________________ | |
480 | void AliHFEdca::CreateHistogramsKfDca(TList *kfDcaList){ | |
481 | // | |
482 | // define histograms: MC dca | |
483 | // | |
484 | ||
485 | // statistics | |
486 | fStat = 0x0; | |
487 | fStat = new TH1I("fStatistics", "allStatistics;ID;counts", 7, -3.5, 3.5); | |
488 | fStat->SetMarkerStyle(20); | |
489 | fStat->SetMarkerColor(3); | |
490 | fStat->SetMarkerSize(1); | |
491 | ||
492 | // for kf dca | |
69ac0e6f | 493 | // fHistKFDcaXY[kNParticles][kNPtBins]=0x0; |
494 | // fHistKFDcaZ[kNParticles][kNPtBins]=0x0; | |
faee3b18 | 495 | |
496 | const Int_t nBins = 1000; | |
497 | const Float_t maxXYBin = 1000.; | |
498 | const Float_t maxZBin = 1000.; | |
499 | ||
500 | ||
501 | for(Int_t k=0; k<kNDcaVar; k++){ | |
502 | TString histTitle((const char*)fgkDcaVarTitle[k]); | |
503 | ||
504 | for(Int_t j=0; j<kNParticles; j++){ | |
505 | for(Int_t i=0; i<kNPtBins; i++){ | |
506 | TString histNameKF((const char*)fgkParticles[j]); | |
507 | histNameKF += Form("_MCpid_KF%s_pT-%.2f-%.2f", (const char*)fgkDcaVar[k], fgkPtIntv[i], fgkPtIntv[i+1]); | |
508 | ||
509 | if(k==0){ | |
510 | fHistKFDcaXY[j][i] = new TH1F((const char*)histNameKF, (const char*)histTitle, nBins, -maxXYBin, maxXYBin); | |
69ac0e6f | 511 | fHistKFDcaXY[j][i]->SetLineColor((int)fgkColorPart[j]); |
faee3b18 | 512 | } |
513 | if(k==1){ | |
514 | fHistKFDcaZ[j][i] = new TH1F((const char*)histNameKF, (const char*)histTitle, nBins, -maxZBin, maxZBin); | |
69ac0e6f | 515 | fHistKFDcaZ[j][i]->SetLineColor((int)fgkColorPart[j]); |
faee3b18 | 516 | } |
517 | } // 50 pt bins | |
518 | } //12 nparticles | |
519 | } // 2 dca var | |
520 | ||
521 | kfDcaList->SetOwner(); | |
522 | kfDcaList->SetName("mcKfDcaDistr"); | |
523 | for(Int_t iPart=0; iPart<kNParticles; iPart++){ | |
524 | for(Int_t iPtBin=0; iPtBin<kNPtBins; iPtBin++){ | |
525 | kfDcaList->Add(fHistKFDcaXY[iPart][iPtBin]); | |
526 | kfDcaList->Add(fHistKFDcaZ[iPart][iPtBin]); | |
527 | } // loop over pt bins | |
528 | } // loop over particles (pos, neg) | |
529 | ||
530 | kfDcaList->Add(fStat); | |
531 | ||
532 | } | |
533 | ||
534 | ||
535 | //________________________________________________________________________ | |
536 | void AliHFEdca::CreateHistogramsDataDca(TList *dataDcaList){ | |
537 | // | |
538 | // define histograms: real Data | |
539 | // | |
540 | ||
541 | // for dca | |
69ac0e6f | 542 | // fHistDataDcaXY[kNParticles][kNPtBins]=0x0; |
543 | // fHistDataDcaZ[kNParticles][kNPtBins]=0x0; | |
544 | // fHistDataWoDcaXY[kNParticles][kNPtBins]=0x0; | |
545 | // fHistDataWoDcaZ[kNParticles][kNPtBins]=0x0; | |
faee3b18 | 546 | |
547 | const Int_t nBins = 1000; | |
548 | const Float_t maxXYBin = 1000.; | |
549 | const Float_t maxZBin = 1000.; | |
550 | ||
551 | for(Int_t k=0; k<kNDcaVar; k++){ | |
552 | TString histTitle((const char*)fgkDcaVarTitle[k]); | |
553 | for(Int_t j=0; j<kNParticles; j++){ | |
554 | for(Int_t i=0; i<kNPtBins; i++){ | |
555 | ||
556 | TString histName((const char*)fgkParticles[j]); | |
557 | histName += Form("_%s_Data_pT-%.2f-%.2f", (const char*)fgkDcaVar[k], fgkPtIntv[i], fgkPtIntv[i+1]); | |
558 | ||
559 | TString histNameWo((const char*)fgkParticles[j]); | |
560 | histNameWo += Form("_%s_Data_wo_pT-%.2f-%.2f", (const char*)fgkDcaVar[k], fgkPtIntv[i], fgkPtIntv[i+1]); | |
561 | ||
562 | if(k==0){ | |
563 | fHistDataDcaXY[j][i] = new TH1F((const char*)histName, (const char*)histTitle, nBins, -maxXYBin, maxXYBin); | |
69ac0e6f | 564 | fHistDataDcaXY[j][i]->SetLineColor((int)fgkColorPart[j]); |
faee3b18 | 565 | |
566 | fHistDataWoDcaXY[j][i] = new TH1F((const char*)histNameWo, (const char*)histTitle, nBins, -maxXYBin, maxXYBin); | |
69ac0e6f | 567 | fHistDataWoDcaXY[j][i]->SetLineColor((int)fgkColorPart[j]); |
faee3b18 | 568 | } |
569 | if(k==1){ | |
570 | fHistDataDcaZ[j][i] = new TH1F((const char*)histName, (const char*)histTitle, nBins, -maxZBin, maxZBin); | |
69ac0e6f | 571 | fHistDataDcaZ[j][i]->SetLineColor((int)fgkColorPart[j]); |
faee3b18 | 572 | |
573 | fHistDataWoDcaZ[j][i] = new TH1F((const char*)histNameWo, (const char*)histTitle, nBins, -maxZBin, maxZBin); | |
69ac0e6f | 574 | fHistDataWoDcaZ[j][i]->SetLineColor((int)fgkColorPart[j]); |
faee3b18 | 575 | } |
576 | } // 50 pt bins | |
577 | } //12 nparticles | |
578 | } // 2 dca var | |
579 | ||
580 | // TList *fDcaList = 0; | |
581 | dataDcaList->SetOwner(); | |
582 | dataDcaList->SetName("dataDcaDistr"); | |
583 | for(Int_t iPart=0; iPart<kNParticles; iPart++){ | |
584 | for(Int_t iPtBin=0; iPtBin<kNPtBins; iPtBin++){ | |
585 | dataDcaList->Add(fHistDataDcaXY[iPart][iPtBin]); | |
586 | dataDcaList->Add(fHistDataDcaZ[iPart][iPtBin]); | |
587 | ||
588 | dataDcaList->Add(fHistDataWoDcaXY[iPart][iPtBin]); | |
589 | dataDcaList->Add(fHistDataWoDcaZ[iPart][iPtBin]); | |
590 | } // loop over pt bins | |
591 | } // loop over particles (pos, neg) | |
592 | ||
593 | ||
594 | } | |
595 | ||
596 | ||
597 | //________________________________________________________________________ | |
598 | void AliHFEdca::CreateHistogramsDataPull(TList *dataPullList){ | |
599 | // define histogram | |
600 | // 2. pull | |
601 | ||
602 | const Int_t nBins = 1000; | |
603 | const Float_t maxXYBin = 20.; | |
604 | const Float_t maxZBin = 20.; | |
605 | ||
606 | // for pull ----------------------------------------------------------------------- | |
69ac0e6f | 607 | // fHistDataDcaXYPull[kNParticles][kNPtBins]=0x0; |
608 | // fHistDataDcaZPull[kNParticles][kNPtBins]=0x0; | |
faee3b18 | 609 | |
69ac0e6f | 610 | // fHistDataWoDcaXYPull[kNParticles][kNPtBins]=0x0; |
611 | // fHistDataWoDcaZPull[kNParticles][kNPtBins]=0x0; | |
faee3b18 | 612 | |
613 | ||
614 | for(Int_t k=0; k<kNDcaVar; k++){ | |
615 | TString histTitle((const char*)fgkPullDataDcaVarTitle[k]); | |
616 | ||
617 | for(Int_t j=0; j<kNParticles; j++){ | |
618 | for(Int_t i=0; i<kNPtBins; i++){ | |
619 | ||
620 | TString histName((const char*)fgkParticles[j]); | |
621 | histName += Form("_%s_Data_pT-%.2f-%.2f", (const char*)fgkPullDcaVar[k], fgkPtIntv[i], fgkPtIntv[i+1]); | |
622 | ||
623 | TString histNameWo((const char*)fgkParticles[j]); | |
624 | histNameWo += Form("_%s_Data_wo_pT-%.2f-%.2f", (const char*)fgkPullDcaVar[k], fgkPtIntv[i], fgkPtIntv[i+1]); | |
625 | ||
626 | if(k==0){ | |
627 | fHistDataDcaXYPull[j][i] = new TH1F((const char*)histName, (const char*)histTitle, nBins, 1-maxXYBin, 1+maxXYBin); | |
69ac0e6f | 628 | fHistDataDcaXYPull[j][i]->SetLineColor((int)fgkColorPart[j]); |
faee3b18 | 629 | |
630 | fHistDataWoDcaXYPull[j][i] = new TH1F((const char*)histNameWo, (const char*)histTitle, nBins, 1-maxXYBin, 1+maxXYBin); | |
69ac0e6f | 631 | fHistDataWoDcaXYPull[j][i]->SetLineColor((int)fgkColorPart[j]); |
faee3b18 | 632 | } |
633 | if(k==1){ | |
634 | fHistDataDcaZPull[j][i] = new TH1F((const char*)histName, (const char*)histTitle, nBins, 1-maxZBin, 1+maxZBin); | |
69ac0e6f | 635 | fHistDataDcaZPull[j][i]->SetLineColor((int)fgkColorPart[j]); |
faee3b18 | 636 | |
637 | fHistDataWoDcaZPull[j][i] = new TH1F((const char*)histNameWo, (const char*)histTitle, nBins, 1-maxZBin, 1+maxZBin); | |
69ac0e6f | 638 | fHistDataWoDcaZPull[j][i]->SetLineColor((int)fgkColorPart[j]); |
70da6c5a | 639 | } |
faee3b18 | 640 | } // 50 pt bins |
641 | } //6 nparticles | |
642 | } // 2 dca var | |
643 | ||
644 | // TList *fDataPullList = 0; | |
645 | dataPullList->SetOwner(); | |
646 | dataPullList->SetName("dataPull"); | |
647 | for(Int_t iPart=0; iPart<kNParticles; iPart++){ | |
648 | for(Int_t iPtBin=0; iPtBin<kNPtBins; iPtBin++){ | |
649 | dataPullList->Add(fHistDataDcaXYPull[iPart][iPtBin]); | |
650 | dataPullList->Add(fHistDataDcaZPull[iPart][iPtBin]); | |
651 | ||
652 | dataPullList->Add(fHistDataWoDcaXYPull[iPart][iPtBin]); | |
653 | dataPullList->Add(fHistDataWoDcaZPull[iPart][iPtBin]); | |
654 | } // loop over pt bins | |
655 | } // loop over particles (pos, neg) | |
656 | ||
657 | } | |
658 | ||
659 | //________________________________________________________________________ | |
660 | void AliHFEdca::CreateHistogramsVertex(TList *vertexList){ | |
661 | // | |
662 | // define histograms: vertex | |
663 | // | |
664 | // for vertex | |
665 | ||
69ac0e6f | 666 | // fHistMCvertex[kNVertexVar]=0x0; |
667 | // fHistESDvertex[kNVertexVar]=0x0; | |
faee3b18 | 668 | |
669 | const Int_t nBins = 1000; | |
670 | const Float_t minXBin = -0.2e4; | |
671 | const Float_t maxXBin = 0.2e4; | |
672 | const Float_t minYBin = -0.5e4; | |
673 | const Float_t maxYBin = 0.5e4; | |
674 | const Float_t minZBin = -1.5e5; | |
675 | const Float_t maxZBin = 1.5e5; | |
676 | ||
677 | const Float_t minBin[kNVertexVar] = {minXBin, minYBin, minZBin}; | |
678 | const Float_t maxBin[kNVertexVar] = {maxXBin, maxYBin, maxZBin}; | |
679 | ||
680 | for(Int_t k=0; k<kNVertexVar; k++){ | |
681 | TString histTitle((const char*)fgkVertexVarTitle[k]); | |
682 | TString histNameMC((const char*)fgkVertexVar[k]); | |
683 | histNameMC += Form("_MC"); | |
684 | TString histNameESD((const char*)fgkVertexVar[k]); | |
685 | histNameESD += Form("_ESD"); | |
686 | ||
687 | fHistMCvertex[k] = new TH1F((const char*)histNameMC, (const char*)histTitle, nBins, minBin[k], maxBin[k]); | |
688 | fHistMCvertex[k]->SetLineColor(k+2); | |
689 | ||
690 | fHistESDvertex[k] = new TH1F((const char*)histNameESD, (const char*)histTitle, nBins, minBin[k], maxBin[k]); | |
691 | fHistESDvertex[k]->SetLineColor(k+2); | |
692 | } // 3 vertex var | |
693 | ||
694 | vertexList->SetOwner(); | |
695 | vertexList->SetName("vertexDistr"); | |
696 | ||
697 | for(Int_t k=0; k<kNVertexVar; k++){ | |
698 | vertexList->Add(fHistMCvertex[k]); | |
699 | vertexList->Add(fHistESDvertex[k]); | |
700 | } | |
701 | ||
702 | } | |
703 | ||
704 | ||
705 | ||
706 | //________________________________________________________________________ | |
707 | void AliHFEdca::CreateHistogramsDataVertex(TList *dataVertexList){ | |
708 | // | |
709 | // define histograms: vertex | |
710 | // | |
711 | // for data vertex | |
712 | ||
69ac0e6f | 713 | // fHistDatavertex[kNVertexVar]=0x0; |
faee3b18 | 714 | |
715 | const Int_t nBins = 1000; | |
716 | const Float_t minXBin = -0.2e4; | |
717 | const Float_t maxXBin = 0.2e4; | |
718 | const Float_t minYBin = -0.5e4; | |
719 | const Float_t maxYBin = 0.5e4; | |
720 | const Float_t minZBin = -1.5e5; | |
721 | const Float_t maxZBin = 1.5e5; | |
722 | ||
723 | const Float_t minBin[kNVertexVar] = {minXBin, minYBin, minZBin}; | |
724 | const Float_t maxBin[kNVertexVar] = {maxXBin, maxYBin, maxZBin}; | |
725 | ||
726 | for(Int_t k=0; k<kNVertexVar; k++){ | |
727 | TString histTitle((const char*)fgkVertexVarTitle[k]); | |
728 | TString histNameDataESD((const char*)fgkVertexVar[k]); | |
729 | histNameDataESD += Form("_data"); | |
730 | ||
731 | fHistDatavertex[k] = new TH1F((const char*)histNameDataESD, (const char*)histTitle, nBins, minBin[k], maxBin[k]); | |
732 | fHistDatavertex[k]->SetLineColor(k+2); | |
733 | } // 3 vertex var | |
734 | ||
735 | // TList *fVDaraVertexList = 0; | |
736 | dataVertexList->SetOwner(); | |
737 | dataVertexList->SetName("dataVertexDistr"); | |
738 | ||
739 | for(Int_t k=0; k<kNVertexVar; k++){ | |
740 | dataVertexList->Add(fHistDatavertex[k]); | |
741 | } | |
742 | ||
743 | } | |
744 | ||
745 | //_______________________________________________________________________________________________ | |
746 | void AliHFEdca::CreateHistogramsPid(TList *mcPidList){ | |
747 | // | |
748 | // define histograms which fills combined PID | |
749 | // | |
750 | ||
751 | const Char_t *mcOResd[2]={"mcPt", "esdPt"}; | |
752 | ||
753 | for(Int_t iPart=0; iPart<kNParticles; iPart++){ | |
754 | TString histTitleMc((const char*)fgkParticles[iPart]); | |
755 | TString histTitleEsd((const char*)fgkParticles[iPart]); | |
756 | histTitleMc += Form("_McPid_%s;p_{T} [GeV/c];counts", mcOResd[0]); | |
757 | histTitleEsd += Form("_EsdPid_%s;p_{T} [GeV/c];counts", mcOResd[1]); | |
758 | ||
759 | TString histNameMc((const char*)fgkParticles[iPart]); | |
760 | TString histNameEsd((const char*)fgkParticles[iPart]); | |
761 | histNameMc+=Form("_McPid_%s", mcOResd[0]); | |
762 | histNameEsd+=Form("_EsdPid_%s", mcOResd[1]); | |
763 | ||
764 | fHistMcPid[iPart] = new TH1F(histNameMc, histTitleMc, kNPtBins, fgkPtIntv); | |
765 | fHistEsdPid[iPart] = new TH1F(histNameEsd, histTitleEsd, kNPtBins, fgkPtIntv); | |
766 | } | |
767 | ||
768 | ||
769 | mcPidList->SetOwner(); | |
770 | mcPidList->SetName("combinedPid"); | |
771 | ||
772 | for(Int_t iPart=0; iPart<kNParticles; iPart++){ | |
773 | mcPidList->Add(fHistMcPid[iPart]); | |
774 | mcPidList->Add(fHistEsdPid[iPart]); | |
775 | } | |
776 | } | |
777 | ||
778 | ||
779 | //_______________________________________________________________________________________________ | |
780 | void AliHFEdca::CreateHistogramsDataPid(TList *pidList){ | |
781 | // | |
782 | // define histograms which fills combined PID: data | |
783 | // | |
784 | ||
785 | ||
786 | ||
787 | for(Int_t iPart=0; iPart<kNParticles; iPart++){ | |
788 | TString histTitleEsd((const char*)fgkParticles[iPart]); | |
789 | histTitleEsd+=Form("_DataEsdPid_esdPt;p_{T} [GeV/c];counts"); | |
790 | TString histNameEsd((const char*)fgkParticles[iPart]); | |
791 | histNameEsd+=Form("_DataEsdPid"); | |
792 | ||
793 | fHistDataEsdPid[iPart] = new TH1F(histNameEsd, histTitleEsd, kNPtBins, fgkPtIntv); | |
794 | } | |
795 | ||
796 | ||
797 | pidList->SetOwner(); | |
798 | pidList->SetName("dataCombinedPid"); | |
799 | ||
800 | for(Int_t iPart=0; iPart<kNParticles; iPart++) | |
801 | pidList->Add(fHistDataEsdPid[iPart]); | |
802 | ||
803 | ||
804 | } | |
805 | ||
806 | ||
807 | ||
808 | //_______________________________________________________________________________________________ | |
c2690925 | 809 | void AliHFEdca::FillHistogramsDca(const AliESDEvent * const esdEvent, const AliESDtrack * const track, AliMCEvent * const mcEvent) |
faee3b18 | 810 | { |
811 | // the kDca plugin | |
812 | // MC vertex | |
813 | AliMCVertex *mcPrimVtx = (AliMCVertex *)mcEvent->GetPrimaryVertex(); | |
814 | Double_t mcPrimV[3]; | |
815 | mcPrimV[0] = mcPrimVtx->GetX(); | |
816 | mcPrimV[1] = mcPrimVtx->GetY(); | |
817 | mcPrimV[2] = mcPrimVtx->GetZ(); | |
818 | ||
819 | Double_t mcVtxXY = TMath::Abs(mcPrimV[0]*mcPrimV[0] + mcPrimV[1]*mcPrimV[1]); | |
820 | ||
821 | // filling historgams track by track | |
822 | // obtaining reconstructed dca ------------------------------------------------------------------ | |
823 | ||
824 | Float_t esdpx = track->Px(); | |
825 | Float_t esdpy = track->Py(); | |
826 | Float_t esdpt = TMath::Sqrt(esdpx*esdpx+esdpy*esdpy); | |
827 | ||
828 | // obtaining errors of dca ------------------------------------------------------------------ | |
829 | const AliESDVertex *primVtx = esdEvent->GetPrimaryVertex(); | |
830 | Double_t primV[3]; | |
831 | primV[0] = primVtx->GetXv(); | |
832 | primV[1] = primVtx->GetYv(); | |
833 | primV[2] = primVtx->GetZv(); | |
834 | ||
835 | Float_t magneticField = 0; // initialized as 5kG | |
836 | magneticField = esdEvent->GetMagneticField(); // in kG | |
837 | ||
838 | Double_t beampiperadius=3.; | |
839 | Double_t dz[2]; // error of dca in cm | |
840 | Double_t covardz[3]; | |
841 | ||
c2690925 | 842 | AliESDtrack ctrack(*track); |
843 | if(!ctrack.PropagateToDCA(primVtx,magneticField, beampiperadius, dz, covardz)) return; // protection | |
70da6c5a | 844 | |
69ac0e6f | 845 | AliMCParticle *mctrack = dynamic_cast<AliMCParticle *>(mcEvent->GetTrack(TMath::Abs(track->GetLabel()))); |
846 | if(!mctrack) return; | |
faee3b18 | 847 | TParticle *part = mctrack->Particle(); |
848 | ||
849 | Float_t vx = part->Vx(); // in cm | |
850 | Float_t vy = part->Vy(); // in cm | |
851 | Float_t vz = part->Vz(); // in cm | |
852 | ||
853 | Float_t vxy = TMath::Sqrt(vx*vx+vy*vy); | |
854 | ||
855 | Float_t mcpx = part->Px(); | |
856 | Float_t mcpy = part->Py(); | |
857 | Float_t mcpt = TMath::Sqrt(mcpx*mcpx+mcpy*mcpy); | |
858 | ||
859 | Int_t pdg = part->GetPdgCode(); | |
860 | Int_t esdPid = GetCombinedPid(track); | |
861 | ||
862 | Int_t charge = 1; | |
863 | if(pdg==kPDGelectron || pdg==kPDGmuon | |
864 | || pdg==-kPDGpion || pdg==-kPDGkaon || pdg==-kPDGproton) charge = -1; | |
865 | ||
866 | // calculate mcDca ------------------------------------------------------------------ | |
867 | const Float_t conv[2] = {1.783/1.6, 2.99792458}; | |
868 | Float_t radiusMc = mcpt/(TMath::Abs(magneticField)/10.)*conv[0]*conv[1]; // pt in GeV/c, magnetic field in Tesla, radius in meter | |
869 | ||
870 | Float_t nx = esdpx/mcpt; | |
871 | Float_t ny = esdpy/mcpt; | |
872 | ||
873 | Float_t radius; | |
874 | radius = TMath::Abs(radiusMc); | |
875 | ||
876 | Double_t dxy = vxy - mcVtxXY; // in cm | |
877 | Double_t dvx = vx - mcPrimV[0]; // in cm | |
878 | Double_t dvy = vy - mcPrimV[1]; // in cm | |
879 | ||
880 | Float_t mcDcaXY = (radius - TMath::Sqrt(dxy*dxy/100./100. + radius*radius + 2*radius*charge*(dvx*ny-dvy*nx)/100.)) ; // in meters | |
881 | ||
882 | Double_t mcDca[2] = {mcDcaXY*100, vz}; // in cm | |
883 | Double_t residual[2] = {0, 0}; | |
884 | Double_t pull[2] = {0, 0}; | |
885 | Double_t error[2] ={TMath::Sqrt(covardz[0]), TMath::Sqrt(covardz[2])}; | |
886 | for(Int_t i=0; i<2; i++){ | |
887 | residual[i] = dz[i] - mcDca[i]; // in centimeters | |
888 | if(error[i]!=0)pull[i] = residual[i]/error[i]; // unitless | |
889 | } | |
890 | ||
891 | ||
892 | for(Int_t iPart=0; iPart<(kNParticles-2); iPart++){ | |
893 | // identified ones | |
894 | for(Int_t iPtBin=0; iPtBin<kNPtBins; iPtBin++){ | |
895 | if(esdpt>fgkPtIntv[iPtBin] && esdpt<=fgkPtIntv[iPtBin+1]){ | |
896 | if(pdg==fgkPdgParticle[iPart]) { | |
897 | fHistDcaXYRes[iPart][iPtBin]->Fill(residual[0]*1.0e4); | |
898 | fHistDcaZRes[iPart][iPtBin]->Fill(residual[1]*1.0e4); | |
899 | fHistDcaXYPull[iPart][iPtBin]->Fill(pull[0]); | |
900 | fHistDcaZPull[iPart][iPtBin]->Fill(pull[1]); | |
901 | fHistDcaXY[iPart][iPtBin]->Fill(dz[0]*1.0e4); | |
902 | fHistDcaZ[iPart][iPtBin]->Fill(dz[1]*1.0e4); | |
903 | } // mc pdg | |
904 | ||
905 | if(esdPid==fgkPdgParticle[iPart]) { | |
906 | fHistEPDcaXYRes[iPart][iPtBin]->Fill(residual[0]*1.0e4); | |
907 | fHistEPDcaZRes[iPart][iPtBin]->Fill(residual[1]*1.0e4); | |
908 | fHistEPDcaXYPull[iPart][iPtBin]->Fill(pull[0]); | |
909 | fHistEPDcaZPull[iPart][iPtBin]->Fill(pull[1]); | |
910 | fHistEPDcaXY[iPart][iPtBin]->Fill(dz[0]*1.0e4); | |
911 | fHistEPDcaZ[iPart][iPtBin]->Fill(dz[1]*1.0e4); | |
912 | } // esd pid | |
913 | ||
914 | } // pt range | |
915 | ||
916 | else | |
917 | continue; | |
918 | } // pt loop | |
919 | } // particle id loop | |
920 | ||
921 | // for charged particles: no pid | |
922 | for(Int_t iPtBin=0; iPtBin<kNPtBins; iPtBin++){ | |
923 | if(esdpt>fgkPtIntv[iPtBin] && esdpt<=fgkPtIntv[iPtBin+1]){ | |
924 | Int_t iPart = 10; | |
925 | if(charge>0) iPart = 11; | |
926 | fHistDcaXYRes[iPart][iPtBin]->Fill(residual[0]*1e4); | |
927 | fHistDcaZRes[iPart][iPtBin]->Fill(residual[1]*1e4); | |
928 | fHistDcaXYPull[iPart][iPtBin]->Fill(pull[0]); | |
929 | fHistDcaZPull[iPart][iPtBin]->Fill(pull[1]); | |
930 | fHistDcaXY[iPart][iPtBin]->Fill(dz[0]*1e4); | |
931 | fHistDcaZ[iPart][iPtBin]->Fill(dz[1]*1e4); | |
932 | } | |
933 | else | |
934 | continue; | |
935 | } // pt | |
936 | ||
937 | } | |
938 | ||
939 | //_______________________________________________________________________________________________ | |
c2690925 | 940 | void AliHFEdca::FillHistogramsKfDca(const AliESDEvent * const esdEvent, const AliESDtrack * const track, const AliMCEvent * const mcEvent) |
6555e2ad | 941 | { |
faee3b18 | 942 | // the kKfDca plugin |
943 | ||
944 | // filling historgams track by track | |
945 | ||
946 | // obtaining reconstructed dca ------------------------------------------------------------------ | |
947 | const AliESDVertex *primVtx = esdEvent->GetPrimaryVertex(); | |
948 | Float_t magneticField = 0; // initialized as 5kG | |
949 | magneticField = esdEvent->GetMagneticField(); // in kG | |
950 | ||
951 | Float_t esdpx = track->Px(); | |
952 | Float_t esdpy = track->Py(); | |
953 | Float_t esdpt = TMath::Sqrt(esdpx*esdpx+esdpy*esdpy); | |
954 | ||
955 | Int_t charge = track->Charge(); | |
956 | ||
957 | Double_t beampiperadius=3.; | |
958 | Double_t dz[2]; // error of dca in cm | |
959 | Double_t covardz[3]; | |
c2690925 | 960 | AliESDtrack ctrack(*track); |
961 | if(!ctrack.PropagateToDCA(primVtx,magneticField, beampiperadius, dz, covardz)) return; // protection | |
faee3b18 | 962 | |
963 | AliMCParticle *mctrack = dynamic_cast<AliMCParticle *>(mcEvent->GetTrack(TMath::Abs(track->GetLabel()))); | |
69ac0e6f | 964 | if(!mctrack) return; |
faee3b18 | 965 | TParticle *part = mctrack->Particle(); |
966 | Int_t pdg = part->GetPdgCode(); | |
967 | ||
968 | // calculate dca using AliKFParticle class------------------------------------------------------------------ | |
969 | Double_t kfdz[3] = {0, 0, 0}; | |
970 | Double_t kfdzwith[3] = {0, 0, 0}; | |
971 | ||
972 | Int_t trkID = track->GetID(); | |
973 | ||
974 | AliKFParticle::SetField(magneticField); | |
975 | AliKFParticle kfParticle(*track, pdg); | |
976 | ||
977 | // prepare kfprimary vertex | |
978 | AliKFVertex kfESDprimary; | |
979 | // Reconstruct Primary Vertex (with ESD tracks) | |
980 | Int_t n=primVtx->GetNIndices(); | |
981 | if (n>0 && primVtx->GetStatus()){ | |
982 | kfESDprimary = AliKFVertex(*primVtx); | |
983 | ||
984 | Double_t dcaXYWithTrk = kfParticle.GetDistanceFromVertexXY(kfESDprimary); | |
985 | Double_t dcaWithTrk = kfParticle.GetDistanceFromVertex(kfESDprimary); | |
986 | Double_t dcaZWithTrk = 0; | |
987 | if(TMath::Abs(dcaWithTrk)>=TMath::Abs(dcaXYWithTrk)) | |
988 | dcaZWithTrk =TMath::Sqrt(dcaWithTrk*dcaWithTrk-dcaXYWithTrk*dcaXYWithTrk)*((dz[1]*-1<=0)?1:-1); | |
989 | kfdzwith[0] = dcaXYWithTrk; | |
990 | kfdzwith[1] = dcaZWithTrk; | |
991 | kfdzwith[2] = dcaWithTrk; // with current track | |
992 | ||
993 | Double_t dcaXYWoTrk = 0; | |
994 | Double_t dcaZWoTrk = 0; | |
995 | Double_t dcaWoTrk = 0; | |
996 | ||
997 | UShort_t *priIndex = primVtx->GetIndices(); | |
998 | ||
999 | for (Int_t i=0;i<n;i++){ | |
1000 | ||
1001 | Int_t idx = Int_t(priIndex[i]); | |
1002 | if (idx == trkID){ | |
1003 | kfESDprimary -= kfParticle; | |
1004 | dcaXYWoTrk = kfParticle.GetDistanceFromVertexXY(kfESDprimary); | |
1005 | dcaWoTrk = kfParticle.GetDistanceFromVertex(kfESDprimary); | |
1006 | if((dcaWoTrk-dcaXYWoTrk)>=0) | |
1007 | dcaZWoTrk = TMath::Abs(dcaWoTrk*dcaWoTrk - dcaXYWoTrk*dcaXYWoTrk)*((dz[1]*-1<=0)?1:-1); | |
1008 | } // remove current track from this calculation | |
1009 | } // loop over all primary vertex contributors | |
1010 | ||
1011 | ||
1012 | kfdz[0] = dcaXYWoTrk; | |
1013 | kfdz[1] = dcaZWoTrk; | |
1014 | kfdz[2] = dcaWoTrk; | |
1015 | ||
1016 | } // only if n contributor > 0 and primVtx constructed | |
1017 | ||
1018 | fStat->Fill(0); | |
1019 | ||
1020 | if(dz[0]!=0 && dz[0]*kfdzwith[0]>0 && TMath::Abs(kfdzwith[0]/dz[0])>0.9999 && TMath::Abs(kfdzwith[0]/dz[0])<1.0001)fStat->Fill(1);; // same | |
1021 | if(dz[0]!=0 && dz[0]*kfdzwith[0]<0 && TMath::Abs(kfdzwith[0]/dz[0])>0.9999 && TMath::Abs(kfdzwith[0]/dz[0])<1.0001) fStat->Fill(2); // swapped sign | |
1022 | if(kfdzwith[0]==0 && dz[0]!=0) fStat->Fill(3); // 0 from KF particle (with current track) | |
1023 | ||
1024 | if(dz[0]!=0 && dz[0]*kfdz[0]>0 && TMath::Abs(kfdz[0]/dz[0])>0.8 && TMath::Abs(kfdz[0]/dz[0])<1.2) fStat->Fill(-1);; // same | |
1025 | if(dz[0]!=0 && dz[0]*kfdz[0]<0 && TMath::Abs(kfdz[0]/dz[0])>0.8 && TMath::Abs(kfdz[0]/dz[0]) <1.2) fStat->Fill(-2); // swapped sign | |
1026 | if(kfdz[0]==0 && dz[0]!=0) fStat->Fill(-3); // 0 from KF particle (without current track) | |
1027 | ||
1028 | for(Int_t iPart=0; iPart<(kNParticles-2); iPart++){ | |
1029 | // identified ones | |
1030 | for(Int_t iPtBin=0; iPtBin<kNPtBins; iPtBin++){ | |
1031 | if(esdpt>fgkPtIntv[iPtBin] && esdpt<=fgkPtIntv[iPtBin+1]){ | |
1032 | if(pdg==fgkPdgParticle[iPart]) { | |
1033 | fHistKFDcaXY[iPart][iPtBin]->Fill(kfdzwith[0]*1.0e4); | |
1034 | fHistKFDcaZ[iPart][iPtBin]->Fill(kfdzwith[1]*1.0e4); | |
1035 | } // mc pdg | |
1036 | } // pt range | |
1037 | ||
1038 | else | |
1039 | continue; | |
1040 | } // pt loop | |
1041 | } // particle id loop | |
1042 | ||
1043 | // for charged particles: no pid | |
1044 | for(Int_t iPtBin=0; iPtBin<kNPtBins; iPtBin++){ | |
1045 | if(esdpt>fgkPtIntv[iPtBin] && esdpt<=fgkPtIntv[iPtBin+1]){ | |
1046 | Int_t iPart = 10; | |
1047 | if(charge>0) iPart = 11; | |
1048 | fHistKFDcaXY[iPart][iPtBin]->Fill(kfdzwith[0]*1e4); | |
1049 | fHistKFDcaZ[iPart][iPtBin]->Fill(kfdzwith[1]*1e4); | |
1050 | ||
1051 | } | |
1052 | else | |
1053 | continue; | |
1054 | } // pt | |
1055 | ||
1056 | } // KF dca | |
1057 | ||
1058 | ||
1059 | //_______________________________________________________________________________________________ | |
c2690925 | 1060 | void AliHFEdca::FillHistogramsVtx(const AliESDEvent *const esdEvent, const AliMCEvent *const mcEvent) |
faee3b18 | 1061 | { |
1062 | ||
1063 | // MC vertex | |
1064 | AliMCVertex *mcPrimVtx = (AliMCVertex *)mcEvent->GetPrimaryVertex(); | |
1065 | Double_t mcPrimV[3]; | |
1066 | mcPrimV[0] = mcPrimVtx->GetX(); | |
1067 | mcPrimV[1] = mcPrimVtx->GetY(); | |
1068 | mcPrimV[2] = mcPrimVtx->GetZ(); | |
1069 | ||
1070 | // obtaining errors of dca ------------------------------------------------------------------ | |
1071 | const AliESDVertex *primVtx = esdEvent->GetPrimaryVertex(); | |
1072 | Double_t primV[3]; | |
1073 | primV[0] = primVtx->GetXv(); | |
1074 | primV[1] = primVtx->GetYv(); | |
1075 | primV[2] = primVtx->GetZv(); | |
1076 | ||
1077 | for(Int_t i=0; i<kNVertexVar; i++){ | |
1078 | fHistMCvertex[i]->Fill(mcPrimV[i]*1.0e4); | |
1079 | fHistESDvertex[i]->Fill(primV[i]*1.0e4); | |
1080 | } | |
1081 | ||
1082 | } | |
1083 | ||
1084 | //_______________________________________________________________________________________________ | |
c2690925 | 1085 | void AliHFEdca::FillHistogramsPid(const AliESDtrack * const track, const AliMCEvent * const mcEvent) |
faee3b18 | 1086 | { |
1087 | ||
1088 | ||
1089 | // filling historgams track by track | |
1090 | Float_t esdpx = track->Px(); | |
1091 | Float_t esdpy = track->Py(); | |
1092 | Float_t esdpt = TMath::Sqrt(esdpx*esdpx+esdpy*esdpy); | |
1093 | ||
1094 | AliMCParticle *mctrack = dynamic_cast<AliMCParticle *>(mcEvent->GetTrack(TMath::Abs(track->GetLabel()))); | |
69ac0e6f | 1095 | if(!mctrack) return; |
faee3b18 | 1096 | TParticle *part = mctrack->Particle(); |
1097 | ||
1098 | Float_t mcpx = part->Px(); | |
1099 | Float_t mcpy = part->Py(); | |
1100 | Float_t mcpt = TMath::Sqrt(mcpx*mcpx+mcpy*mcpy); | |
1101 | ||
1102 | Int_t pdg = part->GetPdgCode(); | |
1103 | Int_t esdPid = GetCombinedPid(track); | |
1104 | ||
1105 | ||
1106 | Double_t ptMom[2] = {mcpt, esdpt}; | |
1107 | // for combined PID | |
1108 | for(Int_t iPart=0; iPart<(kNParticles-2); iPart++){ | |
1109 | if(pdg==fgkPdgParticle[iPart]) // pid all by MC | |
1110 | fHistMcPid[iPart]->Fill(ptMom[0]); | |
1111 | ||
1112 | if(esdPid==fgkPdgParticle[iPart]) // pid all by combined pid | |
1113 | fHistEsdPid[iPart]->Fill(ptMom[1]); | |
1114 | } // loop over particles | |
1115 | ||
1116 | // for charged | |
1117 | if(pdg==kPDGelectron || pdg==kPDGmuon || pdg==-kPDGpion || pdg==-kPDGkaon || pdg==-kPDGproton) | |
1118 | fHistMcPid[10]->Fill(ptMom[0]); | |
1119 | if(pdg==-kPDGelectron || pdg==-kPDGmuon || pdg==kPDGpion || pdg==kPDGkaon || pdg==kPDGproton) | |
1120 | fHistMcPid[11]->Fill(ptMom[0]); | |
1121 | if(esdPid==kPDGelectron || esdPid==kPDGmuon || esdPid==-kPDGpion || esdPid==-kPDGkaon || esdPid==-kPDGproton) | |
1122 | fHistEsdPid[10]->Fill(ptMom[1]); | |
1123 | if(esdPid==-kPDGelectron || esdPid==-kPDGmuon || esdPid==kPDGpion || esdPid==kPDGkaon || esdPid==kPDGproton) | |
1124 | fHistEsdPid[11]->Fill(ptMom[1]); | |
70da6c5a | 1125 | } |
1126 | ||
1127 | ||
faee3b18 | 1128 | ////_______________________________________________________________________________________________ |
c2690925 | 1129 | void AliHFEdca::FillHistogramsDataDca(const AliESDEvent * const esdEvent, const AliESDtrack * const track, const AliESDVertex * const vtxESDSkip) |
faee3b18 | 1130 | { |
1131 | // filling historgams track by track | |
1132 | // obtaining reconstructed dca -------------------------------------------------------------- | |
1133 | ||
1134 | Float_t esdpx = track->Px(); | |
1135 | Float_t esdpy = track->Py(); | |
1136 | Float_t esdpt = TMath::Sqrt(esdpx*esdpx+esdpy*esdpy); | |
1137 | Int_t charge = track->Charge(); | |
1138 | ||
1139 | // obtaining errors of dca ------------------------------------------------------------------ | |
1140 | const AliESDVertex *primVtx = esdEvent->GetPrimaryVertex(); | |
1141 | Double_t primV[3]; | |
1142 | primV[0] = primVtx->GetXv(); | |
1143 | primV[1] = primVtx->GetYv(); | |
1144 | primV[2] = primVtx->GetZv(); | |
1145 | ||
1146 | ||
1147 | Float_t magneticField = 0; // initialized as 5kG | |
1148 | magneticField = esdEvent->GetMagneticField(); // in kG | |
1149 | ||
1150 | Double_t beampiperadius=3.; | |
1151 | Double_t dz[2]; // error of dca in cm | |
1152 | Double_t covardz[3]; | |
1153 | ||
c2690925 | 1154 | AliESDtrack ctrack(*track); |
1155 | if(!ctrack.PropagateToDCA(primVtx,magneticField, beampiperadius, dz, covardz)) return; // protection | |
6555e2ad | 1156 | |
faee3b18 | 1157 | |
1158 | Double_t pull[2] = {0, 0}; | |
1159 | Double_t error[2] ={TMath::Sqrt(covardz[0]), TMath::Sqrt(covardz[2])}; | |
1160 | for(Int_t i=0; i<2; i++){ | |
1161 | if(error[i]!=0)pull[i] = dz[i]/error[i]; // unitless | |
1162 | } | |
1163 | ||
1164 | // get dca when current track is not included | |
1165 | ||
1166 | Double_t dzwo[2], covardzwo[3]; | |
1167 | Double_t pullwo[2] = {0, 0}; | |
c2690925 | 1168 | if(!ctrack.PropagateToDCA(vtxESDSkip, magneticField, beampiperadius, dzwo, covardzwo)) return; // protection |
faee3b18 | 1169 | |
1170 | Double_t errorwo[2] ={TMath::Sqrt(TMath::Abs(covardzwo[0])), TMath::Sqrt(TMath::Abs(covardzwo[2]))}; | |
1171 | for(Int_t i=0; i<2; i++){ | |
1172 | if(errorwo[i]!=0) pullwo[i] = dzwo[i]/errorwo[i]; // unitless | |
1173 | } | |
1174 | ||
1175 | // do pid | |
1176 | Int_t esdPid = GetCombinedPid(track); | |
1177 | ||
1178 | for(Int_t iPart=0; iPart<(kNParticles-2); iPart++){ | |
1179 | // identified ones | |
1180 | for(Int_t iPtBin=0; iPtBin<kNPtBins; iPtBin++){ | |
1181 | if(esdPid==fgkPdgParticle[iPart] && (esdpt>fgkPtIntv[iPtBin] && esdpt<=fgkPtIntv[iPtBin+1])) { | |
1182 | fHistDataDcaXY[iPart][iPtBin]->Fill(dz[0]*1e4); | |
1183 | fHistDataDcaZ[iPart][iPtBin]->Fill(dz[1]*1e4); | |
1184 | fHistDataDcaXYPull[iPart][iPtBin]->Fill(pull[0]); | |
1185 | fHistDataDcaZPull[iPart][iPtBin]->Fill(pull[1]); | |
1186 | // w/o current track | |
1187 | fHistDataWoDcaXY[iPart][iPtBin]->Fill(dzwo[0]*1e4); | |
1188 | fHistDataWoDcaZ[iPart][iPtBin]->Fill(dzwo[1]*1e4); | |
1189 | fHistDataWoDcaXYPull[iPart][iPtBin]->Fill(pullwo[0]); | |
1190 | fHistDataWoDcaZPull[iPart][iPtBin]->Fill(pullwo[1]); | |
1191 | } | |
1192 | else | |
1193 | continue; | |
1194 | } | |
1195 | } | |
1196 | ||
1197 | // for charged particles | |
1198 | for(Int_t iPtBin=0; iPtBin<kNPtBins; iPtBin++){ | |
1199 | if(esdpt>fgkPtIntv[iPtBin] && esdpt<=fgkPtIntv[iPtBin+1]){ | |
1200 | Int_t iPart = 10; | |
1201 | if(charge>0) iPart = 11; | |
1202 | fHistDataDcaXY[iPart][iPtBin]->Fill(dz[0]*1e4); | |
1203 | fHistDataDcaZ[iPart][iPtBin]->Fill(dz[1]*1e4); | |
1204 | fHistDataDcaXYPull[iPart][iPtBin]->Fill(pull[0]); | |
1205 | fHistDataDcaZPull[iPart][iPtBin]->Fill(pull[1]); | |
1206 | // without current track | |
1207 | fHistDataWoDcaXY[iPart][iPtBin]->Fill(dzwo[0]*1e4); | |
1208 | fHistDataWoDcaZ[iPart][iPtBin]->Fill(dzwo[1]*1e4); | |
1209 | fHistDataWoDcaXYPull[iPart][iPtBin]->Fill(pullwo[0]); | |
1210 | fHistDataWoDcaZPull[iPart][iPtBin]->Fill(pullwo[1]); | |
1211 | ||
1212 | } | |
1213 | else | |
1214 | continue; | |
1215 | } | |
1216 | ||
1217 | } | |
1218 | ||
70da6c5a | 1219 | //_______________________________________________________________________________________________ |
c2690925 | 1220 | void AliHFEdca::FillHistogramsDataVtx(const AliESDEvent * const esdEvent) |
faee3b18 | 1221 | { |
1222 | ||
1223 | ||
1224 | // obtaining errors of dca ------------------------------------------------------------------ | |
1225 | const AliESDVertex *primVtx = esdEvent->GetPrimaryVertex(); | |
1226 | Double_t primV[3]; | |
1227 | primV[0] = primVtx->GetXv(); | |
1228 | primV[1] = primVtx->GetYv(); | |
1229 | primV[2] = primVtx->GetZv(); | |
1230 | ||
1231 | // require events with at least 3 contributors for primary vertex construction | |
1232 | Int_t nEsdPrimVtxCtrb = primVtx->GetNContributors(); | |
1233 | if(nEsdPrimVtxCtrb<1) return; // for pass 1, no diomond constrain, each event has at least 1 contributor to Vtx | |
1234 | for(Int_t i=0; i<kNVertexVar; i++) | |
1235 | fHistDatavertex[i]->Fill(primV[i]*1e4); | |
1236 | ||
1237 | } | |
1238 | ||
1239 | ||
1240 | ////_______________________________________________________________________________________________ | |
c2690925 | 1241 | void AliHFEdca::FillHistogramsDataPid(const AliESDtrack * const track) |
70da6c5a | 1242 | { |
1243 | // filling historgams track by track | |
faee3b18 | 1244 | // obtaining reconstructed dca -------------------------------------------------------------- |
1245 | ||
1246 | Float_t esdpx = track->Px(); | |
1247 | Float_t esdpy = track->Py(); | |
1248 | Float_t esdpt = TMath::Sqrt(esdpx*esdpx+esdpy*esdpy); | |
1249 | Int_t charge = track->Charge(); | |
1250 | ||
1251 | Int_t esdPid = GetCombinedPid(track); | |
1252 | ||
1253 | for(Int_t iPart=0; iPart<kNParticles; iPart++){ | |
1254 | if(iPart<(kNParticles-2)){ | |
1255 | if(esdPid==fgkPdgParticle[iPart]) fHistDataEsdPid[iPart]->Fill(esdpt); | |
1256 | } // for identified | |
70da6c5a | 1257 | |
faee3b18 | 1258 | else { |
1259 | if(charge<0) fHistDataEsdPid[10]->Fill(esdpt); | |
1260 | if(charge>0) fHistDataEsdPid[11]->Fill(esdpt); | |
1261 | } | |
1262 | } | |
1263 | } | |
1264 | ||
1265 | //_________________________________________________________________________________________________ | |
c2690925 | 1266 | void AliHFEdca::ApplyExtraCuts(const AliESDEvent * const esdEvent, Int_t nMinPrimVtxContributor) |
faee3b18 | 1267 | { |
1268 | ||
1269 | // | |
1270 | // only one extra cut, number of contributors to each primary vertex | |
1271 | // | |
1272 | ||
1273 | const AliESDVertex *primVtx = esdEvent->GetPrimaryVertex(); | |
1274 | Int_t nEsdPrimVtxCtrb = primVtx->GetNContributors(); | |
1275 | if(nEsdPrimVtxCtrb<nMinPrimVtxContributor) return; | |
1276 | // for pass 1, no diomond constrain, each event has at least 1 contributor to Vtx | |
1277 | ||
1278 | } | |
1279 | ||
1280 | //_____________________________________________________ | |
c2690925 | 1281 | Int_t AliHFEdca::GetCombinedPid(const AliESDtrack *const track) |
faee3b18 | 1282 | { |
1283 | ||
69ac0e6f | 1284 | // combined detector pid |
faee3b18 | 1285 | Double_t prob[AliPID::kSPECIES]; |
1286 | track->GetESDpid(prob); | |
69ac0e6f | 1287 | |
1288 | // setting priors! | |
1289 | Double_t priors[AliPID::kSPECIESN]; | |
1290 | priors[0] = 0.01; | |
1291 | priors[1] = 0.01; | |
1292 | priors[2] = 0.85; | |
1293 | priors[3] = 0.10; | |
1294 | priors[4] = 0.05; | |
faee3b18 | 1295 | |
1296 | Int_t charge = track->Charge(); | |
1297 | Int_t esdPid = -1; | |
1298 | ||
1299 | AliPID pid; | |
1300 | pid.SetPriors(priors); | |
1301 | pid.SetProbabilities(prob); | |
1302 | ||
1303 | // identify particle as the most probable | |
1304 | ||
1305 | Double_t pelectron = pid.GetProbability(AliPID::kElectron); | |
1306 | if(pelectron > pid.GetProbability(AliPID::kMuon) && | |
1307 | pelectron > pid.GetProbability(AliPID::kPion) && | |
1308 | pelectron > pid.GetProbability(AliPID::kKaon) && | |
1309 | pelectron > pid.GetProbability(AliPID::kProton) ) esdPid = -kPDGelectron; | |
1310 | ||
1311 | Double_t pmuon = pid.GetProbability(AliPID::kMuon); | |
1312 | if(pmuon > pid.GetProbability(AliPID::kElectron) && | |
1313 | pmuon > pid.GetProbability(AliPID::kPion) && | |
1314 | pmuon > pid.GetProbability(AliPID::kKaon) && | |
1315 | pmuon > pid.GetProbability(AliPID::kProton) ) esdPid = -kPDGmuon; | |
1316 | ||
1317 | Double_t ppion = pid.GetProbability(AliPID::kPion); | |
1318 | if(ppion > pid.GetProbability(AliPID::kElectron) && | |
1319 | ppion > pid.GetProbability(AliPID::kMuon) && | |
1320 | ppion > pid.GetProbability(AliPID::kKaon) && | |
1321 | ppion > pid.GetProbability(AliPID::kProton) ) esdPid = kPDGpion; | |
1322 | ||
1323 | Double_t pkaon = pid.GetProbability(AliPID::kKaon); | |
1324 | if(pkaon > pid.GetProbability(AliPID::kElectron) && | |
1325 | pkaon > pid.GetProbability(AliPID::kMuon) && | |
1326 | pkaon > pid.GetProbability(AliPID::kPion) && | |
1327 | pkaon > pid.GetProbability(AliPID::kProton) ) esdPid = kPDGkaon; | |
1328 | ||
1329 | Double_t pproton = pid.GetProbability(AliPID::kProton); | |
1330 | if(pproton > pid.GetProbability(AliPID::kElectron) && | |
1331 | pproton > pid.GetProbability(AliPID::kMuon) && | |
1332 | pproton > pid.GetProbability(AliPID::kPion) && | |
1333 | pproton > pid.GetProbability(AliPID::kKaon) ) esdPid = kPDGproton; | |
1334 | ||
1335 | ||
1336 | return charge*esdPid; | |
1337 | ||
1338 | } | |
1339 | ||
1340 | ||
1341 | // for the HFE pid | |
1342 | ||
1343 | //________________________________________________________________________ | |
1344 | void AliHFEdca::CreateHistogramsHfeDca(TList *hfeDcaList){ | |
1345 | // | |
1346 | // define histograms: hfe pid electrons in MC | |
1347 | // | |
1348 | ||
1349 | const Int_t nBinsDca = 1000; | |
1350 | const Float_t maxXYBinDca = 1000.; | |
1351 | const Float_t maxZBinDca = 1000.; | |
1352 | ||
1353 | const Int_t nBinsPull = 1000; | |
1354 | const Float_t maxXYBinPull = 20.; | |
1355 | const Float_t maxZBinPull = 20.; | |
1356 | ||
1357 | const Char_t *mcOResd[2]={"mcPt", "esdPt"}; | |
1358 | ||
69ac0e6f | 1359 | // fHistHPDcaXY[2][kNPtBins]=0x0; |
1360 | // fHistHPDcaZ[2][kNPtBins]=0x0; | |
1361 | // fHistHPDcaXYRes[2][kNPtBins]=0x0; | |
1362 | // fHistHPDcaZRes[2][kNPtBins]=0x0; | |
1363 | // fHistHPDcaXYPull[2][kNPtBins]=0x0; | |
1364 | // fHistHPDcaZPull[2][kNPtBins]=0x0; | |
faee3b18 | 1365 | |
1366 | ||
1367 | for(Int_t k=0; k<kNDcaVar; k++){ | |
1368 | ||
1369 | TString histTitleDca((const char*)fgkDcaVarTitle[k]); | |
1370 | TString histTitleRes((const char*)fgkPullDcaVarTitle[k]); | |
1371 | TString histTitlePull((const char*)fgkResDcaVarTitle[k]); | |
1372 | ||
1373 | for(Int_t iPart=0; iPart<2; iPart++){ | |
1374 | for(Int_t i=0; i<kNPtBins; i++){ | |
1375 | TString histHPName((const char*)fgkParticles[iPart*5]); | |
1376 | histHPName += Form("_HFEpid_%s_pT-%.2f-%.2f", (const char*)fgkDcaVar[k], fgkPtIntv[i], fgkPtIntv[i+1]); | |
1377 | TString histHPNameRes((const char*)fgkParticles[iPart*5]); | |
1378 | histHPNameRes += Form("_HFEpid_%s_pT-%.2f-%.2f", (const char*)fgkResDcaVar[k], fgkPtIntv[i], fgkPtIntv[i+1]); | |
1379 | TString histHPNamePull((const char*)fgkParticles[iPart*5]); | |
1380 | histHPNamePull += Form("_HFEpid_%s_pT-%.2f-%.2f", (const char*)fgkPullDcaVar[k], fgkPtIntv[i], fgkPtIntv[i+1]); | |
1381 | ||
1382 | if(k==0){ | |
1383 | fHistHPDcaXY[iPart][i] = new TH1F((const char*)histHPName, (const char*)histTitleDca, nBinsDca, 1-maxXYBinDca, 1+maxXYBinDca); | |
69ac0e6f | 1384 | fHistHPDcaXY[iPart][i]->SetLineColor((int)fgkColorPart[iPart*5]); |
faee3b18 | 1385 | fHistHPDcaXYRes[iPart][i] = new TH1F((const char*)histHPNameRes, (const char*)histTitleRes, nBinsDca, 1-maxXYBinDca, 1+maxXYBinDca); |
69ac0e6f | 1386 | fHistHPDcaXYRes[iPart][i]->SetLineColor((int)fgkColorPart[iPart*5]); |
faee3b18 | 1387 | fHistHPDcaXYPull[iPart][i] = new TH1F((const char*)histHPNamePull, (const char*)histTitlePull, nBinsPull, 1-maxXYBinPull, 1+maxXYBinPull); |
69ac0e6f | 1388 | fHistHPDcaXYPull[iPart][i]->SetLineColor((int)fgkColorPart[iPart*5]); |
faee3b18 | 1389 | } |
1390 | ||
1391 | if(k==1){ | |
1392 | fHistHPDcaZ[iPart][i] = new TH1F((const char*)histHPName, (const char*)histTitleDca, nBinsDca, 1-maxZBinDca, 1+maxZBinDca); | |
69ac0e6f | 1393 | fHistHPDcaZ[iPart][i]->SetLineColor((int)fgkColorPart[iPart*5]); |
faee3b18 | 1394 | fHistHPDcaZRes[iPart][i] = new TH1F((const char*)histHPNameRes, (const char*)histTitleRes, nBinsDca, 1-maxZBinDca, 1+maxZBinDca); |
69ac0e6f | 1395 | fHistHPDcaZRes[iPart][i]->SetLineColor((int)fgkColorPart[iPart*5]); |
faee3b18 | 1396 | fHistHPDcaZPull[iPart][i] = new TH1F((const char*)histHPNamePull, (const char*)histTitlePull, nBinsPull, 1-maxZBinPull, 1+maxZBinPull); |
69ac0e6f | 1397 | fHistHPDcaZPull[iPart][i]->SetLineColor((int)fgkColorPart[iPart*5]); |
faee3b18 | 1398 | } |
1399 | } // 50 pt bins | |
1400 | } //2 nparticles | |
1401 | } // 2 dca var | |
1402 | ||
69ac0e6f | 1403 | // fHistHfePid[2][2] = 0x0; //! HFE pid |
faee3b18 | 1404 | for(Int_t id=0; id<2; id++){ |
1405 | for(Int_t iPart=0; iPart<2; iPart++){ | |
1406 | TString histTitleHfe((const char*)fgkParticles[iPart*5]); | |
1407 | histTitleHfe+=Form("_MC_HfePid_esdPt;p_{T} [GeV/c];counts"); | |
1408 | TString histNameHfe((const char*)fgkParticles[iPart*5]); | |
1409 | histNameHfe+=Form("_MC_HfePid_%s", mcOResd[id]); | |
1410 | fHistHfePid[id][iPart] = new TH1F(histNameHfe, histTitleHfe, kNPtBins, fgkPtIntv); | |
1411 | } | |
1412 | } | |
1413 | ||
1414 | hfeDcaList->SetOwner(); | |
1415 | hfeDcaList->SetName("hfeDca"); | |
1416 | for(Int_t iPart=0; iPart<2; iPart++){ | |
1417 | for(Int_t iPtBin=0; iPtBin<kNPtBins; iPtBin++){ | |
1418 | hfeDcaList->Add(fHistHPDcaXY[iPart][iPtBin]); | |
1419 | hfeDcaList->Add(fHistHPDcaZ[iPart][iPtBin]); | |
1420 | hfeDcaList->Add(fHistHPDcaXYRes[iPart][iPtBin]); | |
1421 | hfeDcaList->Add(fHistHPDcaZRes[iPart][iPtBin]); | |
1422 | hfeDcaList->Add(fHistHPDcaXYPull[iPart][iPtBin]); | |
1423 | hfeDcaList->Add(fHistHPDcaZPull[iPart][iPtBin]); | |
1424 | } | |
1425 | for(Int_t id=0; id<2; id++) | |
1426 | hfeDcaList->Add(fHistHfePid[id][iPart]); | |
1427 | } | |
1428 | ||
1429 | } | |
1430 | ||
1431 | ||
1432 | //________________________________________________________________________ | |
1433 | void AliHFEdca::CreateHistogramsHfeDataDca(TList *hfeDataDcaList){ | |
1434 | // | |
1435 | // define histograms: hfe pid electrons in data | |
1436 | // | |
1437 | ||
1438 | const Int_t nBinsDca = 1000; | |
1439 | const Float_t maxXYBinDca = 1000.; | |
1440 | const Float_t maxZBinDca = 1000.; | |
1441 | ||
1442 | const Int_t nBinsPull = 1000; | |
1443 | const Float_t maxXYBinPull = 20.; | |
1444 | const Float_t maxZBinPull = 20.; | |
1445 | ||
1446 | ||
69ac0e6f | 1447 | // fHistHPDataDcaXY[2][kNPtBins]=0x0; |
1448 | // fHistHPDataDcaZ[2][kNPtBins]=0x0; | |
1449 | // fHistHPDataDcaXYPull[2][kNPtBins]=0x0; | |
1450 | // fHistHPDataDcaZPull[2][kNPtBins]=0x0; | |
faee3b18 | 1451 | |
1452 | for(Int_t k=0; k<kNDcaVar; k++){ | |
1453 | TString histTitleDca((const char*)fgkDcaVarTitle[k]); | |
1454 | TString histTitlePull((const char*)fgkPullDcaVarTitle[k]); | |
1455 | for(Int_t iPart=0; iPart<2; iPart++){ | |
1456 | for(Int_t i=0; i<kNPtBins; i++){ | |
1457 | TString histHPName((const char*)fgkParticles[iPart*5]); | |
1458 | histHPName += Form("_HFEpid_%s_pT-%.2f-%.2f", (const char*)fgkDcaVar[k], fgkPtIntv[i], fgkPtIntv[i+1]); | |
1459 | TString histHPNamePull((const char*)fgkParticles[iPart*5]); | |
1460 | histHPNamePull += Form("_HFEpid_%s_pT-%.2f-%.2f", (const char*)fgkPullDcaVar[k], fgkPtIntv[i], fgkPtIntv[i+1]); | |
1461 | ||
1462 | if(k==0){ | |
1463 | fHistHPDataDcaXY[iPart][i] = new TH1F((const char*)histHPName, (const char*)histTitleDca, nBinsDca, 1-maxXYBinDca, 1+maxXYBinDca); | |
69ac0e6f | 1464 | fHistHPDataDcaXY[iPart][i]->SetLineColor((int)fgkColorPart[iPart*5]); |
faee3b18 | 1465 | fHistHPDataDcaXYPull[iPart][i] = new TH1F((const char*)histHPNamePull, (const char*)histTitlePull, nBinsPull, 1-maxXYBinPull, 1+maxXYBinPull); |
69ac0e6f | 1466 | fHistHPDataDcaXYPull[iPart][i]->SetLineColor((int)fgkColorPart[iPart*5]); |
faee3b18 | 1467 | } |
1468 | ||
1469 | if(k==1){ | |
1470 | fHistHPDataDcaZ[iPart][i] = new TH1F((const char*)histHPName, (const char*)histTitleDca, nBinsDca, 1-maxZBinDca, 1+maxZBinDca); | |
69ac0e6f | 1471 | fHistHPDataDcaZ[iPart][i]->SetLineColor((int)fgkColorPart[iPart*5]); |
faee3b18 | 1472 | fHistHPDataDcaZPull[iPart][i] = new TH1F((const char*)histHPNamePull, (const char*)histTitlePull, nBinsDca, 1-maxZBinPull, 1+maxZBinPull); |
69ac0e6f | 1473 | fHistHPDataDcaZPull[iPart][i]->SetLineColor((int)fgkColorPart[iPart*5]); |
faee3b18 | 1474 | } |
1475 | ||
1476 | } // 50 pt bins | |
1477 | } // 2 particle type | |
1478 | } // 2 dca var | |
1479 | ||
6555e2ad | 1480 | //fHistDataHfePid[2] = 0x0; //! HFE pid |
faee3b18 | 1481 | for(Int_t iPart=0; iPart<2; iPart++){ |
1482 | TString histTitleHfe((const char*)fgkParticles[iPart*5]); | |
1483 | histTitleHfe+=Form("_Data_HfePid_esdPt;p_{T} [GeV/c];counts"); | |
1484 | TString histNameHfe((const char*)fgkParticles[iPart*5]); | |
1485 | histNameHfe+=Form("_Data_HfePid"); | |
1486 | fHistDataHfePid[iPart] = new TH1F(histNameHfe, histTitleHfe, kNPtBins, fgkPtIntv); | |
1487 | } | |
1488 | ||
1489 | ||
1490 | hfeDataDcaList->SetOwner(); | |
1491 | hfeDataDcaList->SetName("hfeDataDca"); | |
1492 | for(Int_t iPart=0; iPart<2; iPart++){ | |
1493 | for(Int_t iPtBin=0; iPtBin<kNPtBins; iPtBin++){ | |
1494 | hfeDataDcaList->Add(fHistHPDataDcaXY[iPart][iPtBin]); | |
1495 | hfeDataDcaList->Add(fHistHPDataDcaZ[iPart][iPtBin]); | |
1496 | hfeDataDcaList->Add(fHistHPDataDcaXYPull[iPart][iPtBin]); | |
1497 | hfeDataDcaList->Add(fHistHPDcaZPull[iPart][iPtBin]); | |
1498 | ||
1499 | hfeDataDcaList->Add(fHistDataHfePid[iPart]); | |
1500 | } | |
1501 | } | |
1502 | ||
1503 | ||
1504 | } | |
1505 | ||
1506 | ||
1507 | //_______________________________________________________________________________________________ | |
c2690925 | 1508 | void AliHFEdca::FillHistogramsHfeDca(const AliESDEvent * const esdEvent, const AliESDtrack * const track, const AliMCEvent * const mcEvent) |
faee3b18 | 1509 | { |
1510 | // the kHFEpid plugin | |
1511 | ||
1512 | AliMCVertex *mcPrimVtx = (AliMCVertex *)mcEvent->GetPrimaryVertex(); | |
1513 | Double_t mcPrimV[3]; | |
1514 | mcPrimV[0] = mcPrimVtx->GetX(); | |
1515 | mcPrimV[1] = mcPrimVtx->GetY(); | |
1516 | mcPrimV[2] = mcPrimVtx->GetZ(); | |
1517 | ||
1518 | Double_t mcVtxXY = TMath::Abs(mcPrimV[0]*mcPrimV[0] + mcPrimV[1]*mcPrimV[1]); | |
1519 | ||
1520 | // filling historgams track by track | |
70da6c5a | 1521 | // obtaining reconstructed dca ------------------------------------------------------------------ |
faee3b18 | 1522 | Float_t esdpx = track->Px(); |
1523 | Float_t esdpy = track->Py(); | |
1524 | Float_t esdpt = TMath::Sqrt(esdpx*esdpx+esdpy*esdpy); | |
1525 | ||
1526 | // obtaining errors of dca ------------------------------------------------------------------ | |
1527 | const AliESDVertex *primVtx = esdEvent->GetPrimaryVertex(); | |
1528 | Double_t primV[3]; | |
1529 | primV[0] = primVtx->GetXv(); | |
1530 | primV[1] = primVtx->GetYv(); | |
1531 | primV[2] = primVtx->GetZv(); | |
1532 | ||
1533 | Float_t magneticField = 0; // initialized as 5kG | |
1534 | magneticField = esdEvent->GetMagneticField(); // in kG | |
1535 | ||
1536 | Double_t beampiperadius=3.; | |
1537 | Double_t dz[2]; // error of dca in cm | |
1538 | Double_t covardz[3]; | |
c2690925 | 1539 | AliESDtrack ctrack(*track); |
1540 | if(!ctrack.PropagateToDCA(primVtx,magneticField, beampiperadius, dz, covardz)) return; // protection | |
faee3b18 | 1541 | |
1542 | AliMCParticle *mctrack = dynamic_cast<AliMCParticle *>(mcEvent->GetTrack(TMath::Abs(track->GetLabel()))); | |
69ac0e6f | 1543 | if(!mctrack) return; |
faee3b18 | 1544 | TParticle *part = mctrack->Particle(); |
1545 | ||
1546 | Float_t vx = part->Vx(); // in cm | |
1547 | Float_t vy = part->Vy(); // in cm | |
1548 | Float_t vz = part->Vz(); // in cm | |
1549 | ||
1550 | Float_t vxy = TMath::Sqrt(vx*vx+vy*vy); | |
1551 | ||
1552 | Float_t mcpx = part->Px(); | |
1553 | Float_t mcpy = part->Py(); | |
1554 | Float_t mcpt = TMath::Sqrt(mcpx*mcpx+mcpy*mcpy); | |
1555 | ||
1556 | Int_t pdg = part->GetPdgCode(); | |
1557 | ||
1558 | Int_t charge = 1; | |
1559 | if(pdg==kPDGelectron || pdg==kPDGmuon | |
1560 | || pdg==-kPDGpion || pdg==-kPDGkaon || pdg==-kPDGproton) charge = -1; | |
1561 | ||
1562 | // calculate mcDca ------------------------------------------------------------------ | |
1563 | const Float_t conv[2] = {1.783/1.6, 2.99792458}; | |
1564 | Float_t radiusMc = mcpt/(TMath::Abs(magneticField)/10.)*conv[0]*conv[1]; // pt in GeV/c, magnetic field in Tesla, radius in meter | |
1565 | ||
1566 | Float_t nx = esdpx/mcpt; | |
1567 | Float_t ny = esdpy/mcpt; | |
1568 | ||
1569 | Float_t radius; | |
1570 | radius = TMath::Abs(radiusMc); | |
1571 | ||
1572 | Double_t dxy = vxy - mcVtxXY; // in cm | |
1573 | Double_t dvx = vx - mcPrimV[0]; // in cm | |
1574 | Double_t dvy = vy - mcPrimV[1]; // in cm | |
1575 | ||
1576 | Float_t mcDcaXY = (radius - TMath::Sqrt(dxy*dxy/100./100. + radius*radius + 2*radius*charge*(dvx*ny-dvy*nx)/100.)) ; // in meters | |
1577 | ||
1578 | Double_t mcDca[2] = {mcDcaXY*100, vz}; // in cm | |
1579 | Double_t residual[2] = {0, 0}; | |
1580 | Double_t pull[2] = {0, 0}; | |
1581 | Double_t error[2] ={TMath::Sqrt(covardz[0]), TMath::Sqrt(covardz[2])}; | |
1582 | for(Int_t i=0; i<2; i++){ | |
1583 | residual[i] = dz[i] - mcDca[i]; // in centimeters | |
1584 | if(error[i]!=0)pull[i] = residual[i]/error[i]; // unitless | |
1585 | } | |
1586 | ||
1587 | Int_t iPart = -1; | |
1588 | if(track->Charge()<0) iPart = 0; // electron | |
1589 | if(track->Charge()>0) iPart = 1; // positron | |
1590 | if(track->Charge()==0) { | |
1591 | printf("this is not an electron! Check HFEpid method"); | |
1592 | return; | |
1593 | } | |
1594 | for(Int_t iPtBin=0; iPtBin<kNPtBins; iPtBin++){ | |
1595 | if(esdpt>fgkPtIntv[iPtBin] && esdpt<=fgkPtIntv[iPtBin+1]){ | |
1596 | fHistHPDcaXYRes[iPart][iPtBin]->Fill(residual[0]*1.0e4); | |
1597 | fHistHPDcaZRes[iPart][iPtBin]->Fill(residual[1]*1.0e4); | |
1598 | fHistHPDcaXYPull[iPart][iPtBin]->Fill(pull[0]); | |
1599 | fHistHPDcaZPull[iPart][iPtBin]->Fill(pull[1]); | |
1600 | fHistHPDcaXY[iPart][iPtBin]->Fill(dz[0]*1.0e4); | |
1601 | fHistHPDcaZ[iPart][iPtBin]->Fill(dz[1]*1.0e4); | |
1602 | ||
1603 | } // pt range | |
1604 | ||
1605 | else | |
1606 | continue; | |
1607 | } // pt loop | |
1608 | ||
1609 | fHistHfePid[iPart][0]->Fill(esdpt); | |
1610 | fHistHfePid[iPart][1]->Fill(mcpt); | |
1611 | ||
1612 | } | |
1613 | ||
1614 | ||
1615 | //_______________________________________________________________________________________________ | |
c2690925 | 1616 | void AliHFEdca::FillHistogramsHfeDataDca(const AliESDEvent * const esdEvent, const AliESDtrack * const track, const AliESDVertex * const vtxESDSkip) |
faee3b18 | 1617 | { |
1618 | // filling historgams track by track | |
1619 | // obtaining reconstructed dca -------------------------------------------------------------- | |
1620 | ||
1621 | Float_t esdpx = track->Px(); | |
1622 | Float_t esdpy = track->Py(); | |
1623 | Float_t esdpt = TMath::Sqrt(esdpx*esdpx+esdpy*esdpy); | |
1624 | Int_t charge = track->Charge(); | |
1625 | ||
70da6c5a | 1626 | // obtaining errors of dca ------------------------------------------------------------------ |
6555e2ad | 1627 | const AliESDVertex *primVtx = esdEvent->GetPrimaryVertex(); // UNUSED! |
faee3b18 | 1628 | Double_t primV[3]; |
1629 | primV[0] = primVtx->GetXv(); | |
1630 | primV[1] = primVtx->GetYv(); | |
1631 | primV[2] = primVtx->GetZv(); | |
1632 | ||
1633 | Float_t magneticField = 0; // initialized as 5kG | |
1634 | magneticField = esdEvent->GetMagneticField(); // in kG | |
1635 | Double_t beampiperadius=3.; | |
6555e2ad | 1636 | |
faee3b18 | 1637 | Double_t dz[2]; // error of dca in cm |
1638 | Double_t covardz[3]; | |
6555e2ad | 1639 | |
c2690925 | 1640 | AliESDtrack ctrack(*track); // Propagate on copy track |
1641 | if(!ctrack.PropagateToDCA(vtxESDSkip,magneticField, beampiperadius, dz, covardz)) return; // protection | |
faee3b18 | 1642 | |
1643 | Double_t pull[2] = {0, 0}; | |
1644 | Double_t error[2] ={TMath::Sqrt(covardz[0]), TMath::Sqrt(covardz[2])}; | |
1645 | for(Int_t i=0; i<2; i++){ | |
6555e2ad | 1646 | if(error[i]!=0) pull[i] = dz[i]/error[i]; // unitless |
faee3b18 | 1647 | } |
1648 | ||
1649 | Int_t iPart = -1; | |
1650 | if(charge<0) iPart = 0; // electron | |
1651 | if(charge>0) iPart = 1; // positron | |
1652 | if(charge==0) { | |
1653 | printf("this is not an electron! Check HFEpid method\n"); | |
1654 | return; | |
1655 | } | |
1656 | ||
1657 | for(Int_t iPtBin=0; iPtBin<kNPtBins; iPtBin++){ | |
1658 | if(esdpt>fgkPtIntv[iPtBin] && esdpt<=fgkPtIntv[iPtBin+1]) { | |
1659 | fHistHPDataDcaXY[iPart][iPtBin]->Fill(dz[0]*1e4); | |
1660 | fHistHPDataDcaZ[iPart][iPtBin]->Fill(dz[1]*1e4); | |
1661 | fHistHPDataDcaXYPull[iPart][iPtBin]->Fill(pull[0]); | |
1662 | fHistHPDataDcaZPull[iPart][iPtBin]->Fill(pull[1]); | |
1663 | ||
1664 | } | |
6555e2ad | 1665 | else continue; |
faee3b18 | 1666 | } |
1667 | ||
1668 | fHistDataHfePid[iPart]->Fill(esdpt); | |
1669 | ||
70da6c5a | 1670 | } |
1671 |