git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Changes needed on HP and Alpha
[u/mrichter/AliRoot.git] / ANALYSIS / AliD0toKpi.cxx
CommitLineData
3a9a3487 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//----------------------------------------------------------------------------
17// Implementation of the D0toKpi class
18//
19// Note: the two decay tracks are labelled: 0 (positive track)
20// 1 (negative track)
21//
22// Origin: A. Dainese andrea.dainese@pd.infn.it
23//----------------------------------------------------------------------------
24#include <Riostream.h>
25#include <TH1.h>
26#include <TH2.h>
27#include <TCanvas.h>
28#include <TPaveLabel.h>
dab957fb 29#include <TVector3.h>
3a9a3487 30
31#include "AliD0toKpi.h"
32
33ClassImp(AliD0toKpi)
34
35//----------------------------------------------------------------------------
36AliD0toKpi::AliD0toKpi() {
37 // Default constructor
38
39 fSignal = kFALSE;
40
41 fEvent = 0;
42
43 fTrkNum[0] = 0;
44 fTrkNum[1] = 0;
45
46 fV1x = 0.;
47 fV1y = 0.;
48 fV1z = 0.;
49 fV2x = 0.;
50 fV2y = 0.;
51 fV2z = 0.;
52 fDCA = 0.;
53
54 fPx[0] = 0.;
55 fPy[0] = 0.;
56 fPz[0] = 0.;
57 fPx[1] = 0.;
58 fPy[1] = 0.;
59 fPz[1] = 0.;
60
61 fd0[0] = 0.;
62 fd0[1] = 0.;
63
64 fPdg[0] = 0;
65 fPdg[1] = 0;
66 fMum[0] = 0;
67 fMum[1] = 0;
68
69 fTagPi[0] = 0.;
70 fTagPi[1] = 0.;
71 fTagKa[0] = 0.;
72 fTagKa[1] = 0.;
73 fTagNid[0] = 0.;
74 fTagNid[1] = 0.;
75
76 fWgtAD0=fWgtAD0bar=fWgtBD0=fWgtBD0bar=fWgtCD0=fWgtCD0bar=fWgtDD0=fWgtDD0bar=0;
77
78}
79//----------------------------------------------------------------------------
80AliD0toKpi::AliD0toKpi(Int_t ev,Int_t trkNum[2],
81 Double_t v1[3],Double_t v2[3],
82 Double_t dca,
83 Double_t mom[6],Double_t d0[2]) {
84 // Constructor
85
86 fSignal = kFALSE;
87
88 fEvent = ev;
89 fTrkNum[0] = trkNum[0];
90 fTrkNum[1] = trkNum[1];
91
92 fV1x = v1[0];
93 fV1y = v1[1];
94 fV1z = v1[2];
95 fV2x = v2[0];
96 fV2y = v2[1];
97 fV2z = v2[2];
98 fDCA = dca;
99
100 fPx[0] = mom[0];
101 fPy[0] = mom[1];
102 fPz[0] = mom[2];
103 fPx[1] = mom[3];
104 fPy[1] = mom[4];
105 fPz[1] = mom[5];
106
107 fd0[0] = d0[0];
108 fd0[1] = d0[1];
109
110 fPdg[0] = 0;
111 fPdg[1] = 0;
112 fMum[0] = 0;
113 fMum[1] = 0;
114
115 fTagPi[0] = 0.;
116 fTagPi[1] = 0.;
117 fTagKa[0] = 0.;
118 fTagKa[1] = 0.;
119 fTagNid[0] = 0.;
120 fTagNid[1] = 0.;
121
122 fWgtAD0=fWgtAD0bar=fWgtBD0=fWgtBD0bar=fWgtCD0=fWgtCD0bar=fWgtDD0=fWgtDD0bar=0;
123}
124//----------------------------------------------------------------------------
125AliD0toKpi::~AliD0toKpi() {}
126//____________________________________________________________________________
127AliD0toKpi::AliD0toKpi( const AliD0toKpi& d0toKpi):TObject(d0toKpi) {
128 // dummy copy constructor
129}
130//----------------------------------------------------------------------------
131void AliD0toKpi::ApplyPID(TString pidScheme) {
132
133 const char *tofparampbpb = strstr(pidScheme.Data(),"TOFparam_PbPb");
134 const char *tofparampp = strstr(pidScheme.Data(),"TOFparam_pp");
135
136 if((tofparampbpb || tofparampp) && fPdg[0]==0) {
137 printf("AliD0toKpi::ApplyPID :\n Warning: TOF parameterized PID can be used only for simulation!\n");
138 return;
139 }
140
141 if(tofparampbpb) {
142 // tagging of the positive track
143 if(TMath::Abs(fPdg[0])==211 || TMath::Abs(fPdg[0])==13
144 || TMath::Abs(fPdg[0])==11) { // pion,muon,electron
145 fTagPi[0] = LinearInterpolation(PChild(0),kPiBins_PbPb,kPiBinWidth_PbPb,kPiTagPi_PbPb);
146 fTagNid[0] = 1.-fTagPi[0];
147 fTagKa[0] = 0.;
148 fTagPr[0] = 0.;
149 }
150 if(TMath::Abs(fPdg[0])==321) { // kaon
151 fTagKa[0] = LinearInterpolation(PChild(0),kKBins_PbPb,kKBinWidth_PbPb,kKTagK_PbPb);
152 fTagNid[0] = LinearInterpolation(PChild(0),kKBins_PbPb,kKBinWidth_PbPb,kKTagNid_PbPb);
153 if((fTagNid[0]+fTagKa[0])>1.) fTagNid[0] = 1.-fTagKa[0];
154 fTagPi[0] = 1.-fTagNid[0]-fTagKa[0];
155 fTagPr[0] = 0.;
156 }
157 if(TMath::Abs(fPdg[0])==2212) { // proton
158 fTagPr[0] = LinearInterpolation(PChild(0),kPBins_PbPb,kPBinWidth_PbPb,kPTagP_PbPb);
159 fTagNid[0] = LinearInterpolation(PChild(0),kPBins_PbPb,kPBinWidth_PbPb,kPTagNid_PbPb);
160 if((fTagNid[0]+fTagPr[0])>1.) fTagNid[0] = 1.-fTagPr[0];
161 fTagPi[0] = 1.-fTagNid[0]-fTagPr[0];
162 fTagKa[0] = 0.;
163 }
164 // tagging of the negative track
165 if(TMath::Abs(fPdg[1])==211 || TMath::Abs(fPdg[1])==13
166 || TMath::Abs(fPdg[1])==11) { // pion,muon,electron
167 fTagPi[1] = LinearInterpolation(PChild(1),kPiBins_PbPb,kPiBinWidth_PbPb,kPiTagPi_PbPb);
168 fTagNid[1] = 1.-fTagPi[1];
169 fTagKa[1] = 0.;
170 fTagPr[1] = 0.;
171 }
172 if(TMath::Abs(fPdg[1])==321) { // kaon
173 fTagKa[1] = LinearInterpolation(PChild(1),kKBins_PbPb,kKBinWidth_PbPb,kKTagK_PbPb);
174 fTagNid[1] = LinearInterpolation(PChild(1),kKBins_PbPb,kKBinWidth_PbPb,kKTagNid_PbPb);
175 if((fTagNid[1]+fTagKa[1])>1.) fTagNid[1] = 1.-fTagKa[1];
176 fTagPi[1] = 1.-fTagNid[1]-fTagKa[1];
177 fTagPr[1] = 0.;
178 }
179 if(TMath::Abs(fPdg[1])==2212) { // proton
180 fTagPr[1] = LinearInterpolation(PChild(1),kPBins_PbPb,kPBinWidth_PbPb,kPTagP_PbPb);
181 fTagNid[1] = LinearInterpolation(PChild(1),kPBins_PbPb,kPBinWidth_PbPb,kPTagNid_PbPb);
182 if((fTagNid[1]+fTagPr[1])>1.) fTagNid[1] = 1.-fTagPr[1];
183 fTagPi[1] = 1.-fTagNid[1]-fTagPr[1];
184 fTagKa[1] = 0.;
185 }
186 }
187
188
189 if(tofparampp) {
190 // tagging of the positive track
191 if(TMath::Abs(fPdg[0])==211 || TMath::Abs(fPdg[0])==13
192 || TMath::Abs(fPdg[0])==11) { // pion,muon,electron
193 fTagPi[0] = LinearInterpolation(PChild(0),kPiBins_pp,kPiBinWidth_pp,kPiTagPi_pp);
194 fTagNid[0] = 1.-fTagPi[0];
195 fTagKa[0] = 0.;
196 fTagPr[0] = 0.;
197 }
198 if(TMath::Abs(fPdg[0])==321) { // kaon
199 fTagKa[0] = LinearInterpolation(PChild(0),kKBins_pp,kKBinWidth_pp,kKTagK_pp);
200 fTagNid[0] = LinearInterpolation(PChild(0),kKBins_pp,kKBinWidth_pp,kKTagNid_pp);
201 if((fTagNid[0]+fTagKa[0])>1.) fTagNid[0] = 1.-fTagKa[0];
202 fTagPi[0] = 1.-fTagNid[0]-fTagKa[0];
203 fTagPr[0] = 0.;
204 }
205 if(TMath::Abs(fPdg[0])==2212) { // proton
206 fTagPr[0] = LinearInterpolation(PChild(0),kPBins_pp,kPBinWidth_pp,kPTagP_pp);
207 fTagNid[0] = LinearInterpolation(PChild(0),kPBins_pp,kPBinWidth_pp,kPTagNid_pp);
208 if((fTagNid[0]+fTagPr[0])>1.) fTagNid[0] = 1.-fTagPr[0];
209 fTagPi[0] = 1.-fTagNid[0]-fTagPr[0];
210 fTagKa[0] = 0.;
211 }
212 // tagging of the negative track
213 if(TMath::Abs(fPdg[1])==211 || TMath::Abs(fPdg[1])==13
214 || TMath::Abs(fPdg[1])==11) { // pion,muon,electron
215 fTagPi[1] = LinearInterpolation(PChild(1),kPiBins_pp,kPiBinWidth_pp,kPiTagPi_pp);
216 fTagNid[1] = 1.-fTagPi[1];
217 fTagKa[1] = 0.;
218 fTagPr[1] = 0.;
219 }
220 if(TMath::Abs(fPdg[1])==321) { // kaon
221 fTagKa[1] = LinearInterpolation(PChild(1),kKBins_pp,kKBinWidth_pp,kKTagK_pp);
222 fTagNid[1] = LinearInterpolation(PChild(1),kKBins_pp,kKBinWidth_pp,kKTagNid_pp);
223 if((fTagNid[1]+fTagKa[1])>1.) fTagNid[1] = 1.-fTagKa[1];
224 fTagPi[1] = 1.-fTagNid[1]-fTagKa[1];
225 fTagPr[1] = 0.;
226 }
227 if(TMath::Abs(fPdg[1])==2212) { // proton
228 fTagPr[1] = LinearInterpolation(PChild(1),kPBins_pp,kPBinWidth_pp,kPTagP_pp);
229 fTagNid[1] = LinearInterpolation(PChild(1),kPBins_pp,kPBinWidth_pp,kPTagNid_pp);
230 if((fTagNid[1]+fTagPr[1])>1.) fTagNid[1] = 1.-fTagPr[1];
231 fTagPi[1] = 1.-fTagNid[1]-fTagPr[1];
232 fTagKa[1] = 0.;
233 }
234 }
235
236 return;
237}
238//----------------------------------------------------------------------------
239Double_t AliD0toKpi::ChildrenRelAngle() const {
240 // relative angle between K and pi
241
242 TVector3 mom0(fPx[0],fPy[0],fPz[0]);
243 TVector3 mom1(fPx[1],fPy[1],fPz[1]);
244
245 Double_t angle = mom0.Angle(mom1);
246
247 return angle;
248}
249//----------------------------------------------------------------------------
250void AliD0toKpi::ComputeWgts() {
251 // calculate the weights for PID
252
253
254 // assignement of the weights from PID
255 fWgtAD0 = fTagKa[1]*(fTagPi[0]+fTagNid[0]);
256 fWgtAD0bar = fTagKa[0]*(fTagPi[1]+fTagNid[1]);
257 fWgtBD0 = fTagPi[0]*fTagNid[1];
258 fWgtBD0bar = fTagPi[1]*fTagNid[0];
259 fWgtCD0 = fTagNid[0]*fTagNid[1];
260 fWgtCD0bar = fTagNid[0]*fTagNid[1];
261 fWgtDD0 = 1.-fWgtAD0-fWgtBD0-fWgtCD0;
262 fWgtDD0bar = 1.-fWgtAD0bar-fWgtBD0bar-fWgtCD0bar;
263
264 /*
265 cerr<<fWgtAD0<<" "<<fWgtAD0bar<<endl;
266 cerr<<fWgtBD0<<" "<<fWgtBD0bar<<endl;
267 cerr<<fWgtCD0<<" "<<fWgtCD0bar<<endl;
268
269 if(fWgtAD0<0.) cerr<<"AliD0toKpi::ComputeWgts() Negative weight!!!\n";
270 if(fWgtAD0bar<0.) cerr<<"AliD0toKpi::ComputeWgts() Negative weight!!!\n";
271 if(fWgtBD0<0.) cerr<<"AliD0toKpi::ComputeWgts() Negative weight!!!\n";
272 if(fWgtBD0bar<0.) cerr<<"AliD0toKpi::ComputeWgts() Negative weight!!!\n";
273 if(fWgtCD0<0.) cerr<<"AliD0toKpi::ComputeWgts() Negative weight!!!\n";
274 if(fWgtCD0bar<0.) cerr<<"AliD0toKpi::ComputeWgts() Negative weight!!!\n";
275 */
276
277 return;
278}
279//----------------------------------------------------------------------------
280void AliD0toKpi::CorrectWgt4BR(Double_t factor) {
281 // correct weights of background from charm
282
283 fWgtAD0 *= factor;
284 fWgtAD0bar *= factor;
285 fWgtBD0 *= factor;
286 fWgtBD0bar *= factor;
287 fWgtCD0 *= factor;
288 fWgtCD0bar *= factor;
289 fWgtDD0 *= factor;
290 fWgtDD0bar *= factor;
291
292 return;
293}
294//----------------------------------------------------------------------------
295Double_t AliD0toKpi::CosPointing() const {
296 // cosine of pointing angle in space
297
298 TVector3 mom(Px(),Py(),Pz());
299 TVector3 flight(fV2x-fV1x,fV2y-fV1y,fV2z-fV1z);
300
301 Double_t pta = mom.Angle(flight);
302
303 return TMath::Cos(pta);
304}
305//----------------------------------------------------------------------------
306Double_t AliD0toKpi::CosPointingXY() const {
307 // cosine of pointing angle in transverse plane
308
309 TVector3 momXY(Px(),Py(),0.);
310 TVector3 flightXY(fV2x-fV1x,fV2y-fV1y,0.);
311
312 Double_t ptaXY = momXY.Angle(flightXY);
313
314 return TMath::Cos(ptaXY);
315}
316//----------------------------------------------------------------------------
317void AliD0toKpi::CosThetaStar(Double_t &ctsD0,Double_t &ctsD0bar) const {
318 // cosine of decay angle in the D0 rest frame
319
320 Double_t pStar = TMath::Sqrt(TMath::Power(kMD0*kMD0-kMK*kMK-kMPi*kMPi,2.)-4.*kMK*kMK*kMPi*kMPi)/(2.*kMD0);
321
322 Double_t beta = P()/Energy();
323 Double_t gamma = Energy()/kMD0;
324
325 ctsD0 = (Ql(1)/gamma-beta*TMath::Sqrt(pStar*pStar+kMK*kMK))/pStar;
326 // if(ctsD0 > 1.) { cerr<<"AliD0toKpi::CosThetaStar: > 1 "<<ctsD0<<"!\n"; }
327 // if(ctsD0 < -1.) { cerr<<"AliD0toKpi::CosThetaStar: < -1 "<<ctsD0<<"!\n"; }
328
329 ctsD0bar = (Ql(0)/gamma-beta*TMath::Sqrt(pStar*pStar+kMK*kMK))/pStar;
330 // if(ctsD0bar > 1.) { cerr<<"AliD0toKpi::CosThetaStar: > 1 "<<ctsD0bar<<"!\n"; }
331 // if(ctsD0bar < -1.) { cerr<<"AliD0toKpi::CosThetaStar: < -1 "<<ctsD0bar<<"!\n";}
332
333 return;
334}
335//----------------------------------------------------------------------------
336Double_t AliD0toKpi::Eta() const {
337 // pseudorapidity of the D0
338
339 Double_t theta = TMath::Pi()/2.-TMath::ATan2(Pz(),Pt());
340 Double_t eta = -TMath::Log(TMath::Tan(theta/2.));
341 return eta;
342}
343//----------------------------------------------------------------------------
344Double_t AliD0toKpi::EtaChild(Int_t child) const {
345 // pseudorapidity of the decay tracks
346
347 Double_t theta = TMath::Pi()/2.-TMath::ATan2(fPz[child],PtChild(child));
348 Double_t eta = -TMath::Log(TMath::Tan(theta/2.));
349 return eta;
350}
351//----------------------------------------------------------------------------
352void AliD0toKpi::GetWgts(Double_t &WgtD0,Double_t &WgtD0bar,TString sample)
353 const {
354 // returns the weights for pid
355
356 const char *sampleA = strstr(sample.Data(),"A");
357 const char *sampleB = strstr(sample.Data(),"B");
358 const char *sampleC = strstr(sample.Data(),"C");
359 const char *sampleD = strstr(sample.Data(),"D");
360 const char *sampleABCD = strstr(sample.Data(),"ABCD");
361 const char *sampleABC = strstr(sample.Data(),"ABC");
362 const char *sampleBC = strstr(sample.Data(),"BC");
363
364 if(sampleA) { WgtD0 = fWgtAD0; WgtD0bar = fWgtAD0bar; }
365 if(sampleB) { WgtD0 = fWgtBD0; WgtD0bar = fWgtBD0bar; }
366 if(sampleC) { WgtD0 = fWgtCD0; WgtD0bar = fWgtCD0bar; }
367 if(sampleD) { WgtD0 = fWgtDD0; WgtD0bar = fWgtDD0bar; }
368 if(sampleABCD) {
369 WgtD0 = fWgtAD0+fWgtBD0+fWgtCD0+fWgtDD0;
370 WgtD0bar = fWgtAD0bar+fWgtBD0bar+fWgtCD0bar+fWgtDD0bar;
371 }
372 if(sampleABC) {
373 WgtD0 = fWgtAD0+fWgtBD0+fWgtCD0;
374 WgtD0bar = fWgtAD0bar+fWgtBD0bar+fWgtCD0bar;
375 }
376 if(sampleBC) {
377 WgtD0 = fWgtBD0+fWgtCD0;
378 WgtD0bar = fWgtBD0bar+fWgtCD0bar;
379 }
380
381
382 if(fSignal) {
383 if(fMum[0]==421) WgtD0bar = 0.;
384 if(fMum[0]==-421) WgtD0 = 0.;
385 }
386
387 return;
388}
389//----------------------------------------------------------------------------
390void AliD0toKpi::InvMass(Double_t &mD0,Double_t &mD0bar) const {
391 // invariant mass as D0 and as D0bar
392
393 Double_t energy[2];
394
395 // D0 -> K- Pi+
396 energy[1] = TMath::Sqrt(kMK*kMK+PChild(1)*PChild(1));
397 energy[0] = TMath::Sqrt(kMPi*kMPi+PChild(0)*PChild(0));
398
399 mD0 = TMath::Sqrt((energy[0]+energy[1])*(energy[0]+energy[1])-P()*P());
400
401
402 // D0bar -> K+ Pi-
403 energy[0] = TMath::Sqrt(kMK*kMK+PChild(0)*PChild(0));
404 energy[1] = TMath::Sqrt(kMPi*kMPi+PChild(1)*PChild(1));
405
406 mD0bar = TMath::Sqrt((energy[0]+energy[1])*(energy[0]+energy[1])-P()*P());
407
408 return;
409
410}
411//----------------------------------------------------------------------------
412Double_t AliD0toKpi::Ql(Int_t child) const {
413 // longitudinal momentum of decay tracks w.r.t. to D0 momentum
414
415 Double_t qL;
416 TVector3 mom(fPx[child],fPy[child],fPz[child]);
417 TVector3 momD(Px(),Py(),Pz());
418
419 qL = mom.Dot(momD)/momD.Mag();
420
421 return qL ;
422}
423//----------------------------------------------------------------------------
424Double_t AliD0toKpi::Qt() const {
425 // transverse momentum of decay tracks w.r.t. to D0 momentum
426
427 TVector3 mom0(fPx[0],fPy[0],fPz[0]);
428 TVector3 momD(Px(),Py(),Pz());
429
430 return mom0.Perp(momD);
431}
432//----------------------------------------------------------------------------
433Bool_t AliD0toKpi::Select(const Double_t* cuts,Int_t& okD0,Int_t& okD0bar)
434 const {
435//
436// This function compares the D0 with a set of cuts:
437//
438// cuts[0] = inv. mass half width [GeV]
439// cuts[1] = dca [micron]
440// cuts[2] = cosThetaStar
441// cuts[3] = pTK [GeV/c]
442// cuts[4] = pTPi [GeV/c]
443// cuts[5] = d0K [micron] upper limit!
444// cuts[6] = d0Pi [micron] upper limit!
445// cuts[7] = d0d0 [micron^2]
446// cuts[8] = cosThetaPoint
447//
448// If the D0/D0bar doesn't pass the cuts it sets the weights to 0
449// If neither D0 nor D0bar pass the cuts return kFALSE
450//
451 Double_t mD0,mD0bar,ctsD0,ctsD0bar;
452 okD0=1; okD0bar=1;
453
454 if(PtChild(1) < cuts[3] || PtChild(0) < cuts[4]) okD0 = 0;
455 if(PtChild(0) < cuts[3] || PtChild(1) < cuts[4]) okD0bar = 0;
456 if(!okD0 && !okD0bar) return kFALSE;
457
458 if(TMath::Abs(Getd0Child(1)) > cuts[5] ||
459 TMath::Abs(Getd0Child(0)) > cuts[6]) okD0 = 0;
460 if(TMath::Abs(Getd0Child(0)) > cuts[6] ||
461 TMath::Abs(Getd0Child(1)) > cuts[5]) okD0bar = 0;
462 if(!okD0 && !okD0bar) return kFALSE;
463
464 if(GetDCA() > cuts[1]) { okD0 = okD0bar = 0; return kFALSE; }
465
466 InvMass(mD0,mD0bar);
467 if(TMath::Abs(mD0-kMD0) > cuts[0]) okD0 = 0;
468 if(TMath::Abs(mD0bar-kMD0) > cuts[0]) okD0bar = 0;
469 if(!okD0 && !okD0bar) return kFALSE;
470
471 CosThetaStar(ctsD0,ctsD0bar);
472 if(TMath::Abs(ctsD0) > cuts[2]) okD0 = 0;
473 if(TMath::Abs(ctsD0bar) > cuts[2]) okD0bar = 0;
474 if(!okD0 && !okD0bar) return kFALSE;
475
476 if(ProdImpParams() > cuts[7]) { okD0 = okD0bar = 0; return kFALSE; }
477
478 if(CosPointing() < cuts[8]) { okD0 = okD0bar = 0; return kFALSE; }
479
480 return kTRUE;
481}
482//-----------------------------------------------------------------------------
483void AliD0toKpi::SetPIDresponse(Double_t resp0[5],Double_t resp1[5]) {
484 // Set combined PID detector response probabilities
485
486 fPIDrespEl[0] = resp0[0];
487 fPIDrespEl[1] = resp1[0];
488 fPIDrespMu[0] = resp0[1];
489 fPIDrespMu[1] = resp1[1];
490 fPIDrespPi[0] = resp0[2];
491 fPIDrespPi[1] = resp1[2];
492 fPIDrespKa[0] = resp0[3];
493 fPIDrespKa[1] = resp1[3];
494 fPIDrespPr[0] = resp0[4];
495 fPIDrespPr[1] = resp1[4];
496
497 return;
498}
499//-----------------------------------------------------------------------------
500void AliD0toKpi::DrawPIDinTOF(TString pidScheme) const {
501 // Draw parameterized PID probabilities in TOF
502
503 const char *tofparampbpb = strstr(pidScheme.Data(),"TOFparam_PbPb");
504 const char *tofparampp = strstr(pidScheme.Data(),"TOFparam_pp");
505
506 TH2F* framePi = new TH2F("framePi","Tag probabilities for PIONS",2,0,2.5,2,0,1);
507 framePi->SetXTitle("p [GeV/c]");
508 framePi->SetStats(0);
509 TH2F* frameK = new TH2F("frameK","Tag probabilities for KAONS",2,0,2.5,2,0,1);
510 frameK->SetXTitle("p [GeV/c]");
511 frameK->SetStats(0);
512 TH2F* frameP = new TH2F("frameP","Tag probabilities for PROTONS",2,0,4.5,2,0,1);
513 frameP->SetXTitle("p [GeV/c]");
514 frameP->SetStats(0);
515
516 TH1F* hPiPi = new TH1F("hPiPi","Tag probabilities for PIONS",kPiBins_PbPb,0,2.5);
517 TH1F* hPiNid = new TH1F("hPiNid","Tag probabilities for PIONS",kPiBins_PbPb,0,2.5);
518
519 TH1F* hKK = new TH1F("hKK","Tag probabilities for KAONS",kKBins_PbPb,0,2.5);
520 TH1F* hKNid = new TH1F("hKNid","Tag probabilities for KAONS",kKBins_PbPb,0,2.5);
521 TH1F* hKPi = new TH1F("hKPi","Tag probabilities for KAONS",kKBins_PbPb,0,2.5);
522
523 TH1F* hPP = new TH1F("hPP","Tag probabilities for PROTONS",kPBins_PbPb,0,4.5);
524 TH1F* hPNid = new TH1F("hPNid","Tag probabilities for PROTONS",kPBins_PbPb,0,4.5);
525 TH1F* hPPi = new TH1F("hPPi","Tag probabilities for PROTONS",kPBins_PbPb,0,4.5);
526
527
528 if(tofparampbpb) {
529
530 for(Int_t i=1; i<=kPiBins_PbPb; i++) {
531 hPiPi->SetBinContent(i,kPiTagPi_PbPb[i-1]);
532 hPiNid->SetBinContent(i,kPiTagPi_PbPb[i-1]+kPiTagNid_PbPb[i-1]);
533
534 hKK->SetBinContent(i,kKTagK_PbPb[i-1]);
535 hKPi->SetBinContent(i,kKTagK_PbPb[i-1]+kKTagPi_PbPb[i-1]);
536 hKNid->SetBinContent(i,kKTagK_PbPb[i-1]+kKTagPi_PbPb[i-1]+kKTagNid_PbPb[i-1]);
537 }
538 for(Int_t i=1; i<=kPBins_PbPb; i++) {
539 hPP->SetBinContent(i,kPTagP_PbPb[i-1]);
540 hPPi->SetBinContent(i,kPTagP_PbPb[i-1]+kPTagPi_PbPb[i-1]);
541 hPNid->SetBinContent(i,kPTagP_PbPb[i-1]+kPTagPi_PbPb[i-1]+kPTagNid_PbPb[i-1]);
542 }
543
544 } else if(tofparampp) {
545
546 for(Int_t i=1; i<=kPiBins_pp; i++) {
547 hPiPi->SetBinContent(i,kPiTagPi_pp[i-1]);
548 hPiNid->SetBinContent(i,kPiTagPi_pp[i-1]+kPiTagNid_pp[i-1]);
549
550 hKK->SetBinContent(i,kKTagK_pp[i-1]);
551 hKPi->SetBinContent(i,kKTagK_pp[i-1]+kKTagPi_pp[i-1]);
552 hKNid->SetBinContent(i,kKTagK_pp[i-1]+kKTagPi_pp[i-1]+kKTagNid_pp[i-1]);
553 }
554 for(Int_t i=1; i<=kPBins_pp; i++) {
555 hPP->SetBinContent(i,kPTagP_pp[i-1]);
556 hPPi->SetBinContent(i,kPTagP_pp[i-1]+kPTagPi_pp[i-1]);
557 hPNid->SetBinContent(i,kPTagP_pp[i-1]+kPTagPi_pp[i-1]+kPTagNid_pp[i-1]);
558 }
559
560 }
561
562
563 TCanvas* c = new TCanvas("c","Parameterized PID in TOF",0,0,1000,400);
564 c->Divide(3,1);
565 c->cd(1);
566 framePi->Draw();
567 hPiNid->SetFillColor(18); hPiNid->Draw("same");
568 hPiPi->SetFillColor(4); hPiPi->Draw("same");
569 TPaveLabel* pav1 = new TPaveLabel(1,.2,1.4,.3,"#pi");
570 pav1->SetBorderSize(0);
571 pav1->Draw("same");
572 TPaveLabel* pav2 = new TPaveLabel(1,.8,1.8,.9,"non-id");
573 pav2->SetBorderSize(0);
574 pav2->Draw("same");
575
576 c->cd(2);
577 frameK->Draw();
578 hKNid->SetFillColor(18); hKNid->Draw("same");
579 hKPi->SetFillColor(4); hKPi->Draw("same");
580 hKK->SetFillColor(7); hKK->Draw("same");
581 TPaveLabel* pav3 = new TPaveLabel(1,.2,1.5,.3,"K");
582 pav3->SetBorderSize(0);
583 pav3->Draw("same");
584 TPaveLabel* pav4 = new TPaveLabel(1,.8,1.8,.9,"non-id");
585 pav4->SetBorderSize(0);
586 pav4->Draw("same");
587 TPaveLabel* pav5 = new TPaveLabel(.4,.5,.8,.6,"#pi");
588 pav5->SetBorderSize(0);
589 pav5->Draw("same");
590
591 c->cd(3);
592 frameP->Draw();
593 hPNid->SetFillColor(18); hPNid->Draw("same");
594 hPPi->SetFillColor(4); hPPi->Draw("same");
595 hPP->SetFillColor(3); hPP->Draw("same");
596 TPaveLabel* pav6 = new TPaveLabel(1,.2,1.5,.3,"p");
597 pav6->SetBorderSize(0);
598 pav6->Draw("same");
599 TPaveLabel* pav7 = new TPaveLabel(1,.8,2.6,.9,"non-id");
600 pav7->SetBorderSize(0);
601 pav7->Draw("same");
602 TPaveLabel* pav8 = new TPaveLabel(.2,.5,1,.6,"#pi");
603 pav8->SetBorderSize(0);
604 pav8->Draw("same");
605
606
607 return;
608}
609//----------------------------------------------------------------------------
610Double_t AliD0toKpi::LinearInterpolation(Double_t p,Int_t nBins,Double_t Bin,
611 const Double_t *values) const {
612 // a linear interpolation method
613
614 Double_t value=0;
615 Double_t slope;
616
617 if(p<0.5*Bin) {
618 value = values[0];
619 } else if(p>=(nBins-0.5)*Bin) {
620 slope = (2*values[nBins-1]-values[nBins-2]-values[nBins-3])/Bin/2;
621 value = values[nBins-2]+slope*(p-Bin*(nBins-1.5));
622 } else {
623 for(Int_t i=0; i<nBins; i++) {
624 if(p<(i+0.5)*Bin) {
625 slope = (values[i]-values[i-1])/Bin;
626 value = values[i-1]+slope*(p-Bin*(i-0.5));
627 break;
628 }
629 }
630 }
631
632 if(value<0.) value=0.;
633 if(value>1.) value=1.;
634
635 return value;
636}
637//----------------------------------------------------------------------------
638
639
640
641
642
643
644/*
645//____________________________________________________________________________
646void AliD0toKpi::SetPtWgts4pp() {
647 // Correct pt distribution in order to reproduce MNR pt slope
648 // (for pp generated with PYTHIA min. bias)
649
650 if(TMath::Abs(fMum[0]) != 421 && TMath::Abs(fMum[1]) != 421 &&
651 TMath::Abs(fMum[0]) != 411 && TMath::Abs(fMum[1]) != 411) return;
652
653 Double_t ptWgt = 1.;
654 ptWgt = 2.05-0.47*Pt()+0.02*Pt()*Pt();
655 if(Pt() >= 5.) ptWgt = 0.56*TMath::Exp(-0.12*Pt());
656
657 fWgtAD0 *= ptWgt;
658 fWgtAD0bar *= ptWgt;
659 fWgtBD0 *= ptWgt;
660 fWgtBD0bar *= ptWgt;
661 fWgtCD0 *= ptWgt;
662 fWgtCD0bar *= ptWgt;
663 fWgtDD0 *= ptWgt;
664 fWgtDD0bar *= ptWgt;
665
666 return;
667}
668//____________________________________________________________________________
669*/
670
671
672
673
674
675
676
677
678