PWG3/AliD0toKpi.cxx

   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 //                  for pp and PbPb interactions
  19 // Note: the two decay tracks are labelled: 0 (positive track)
  20 //                                          1 (negative track)
  21 //            Origin: A. Dainese    andrea.dainese@lnl.infn.it
  22 //----------------------------------------------------------------------------
  23
  24 #include <TH1.h>
  25 #include <TH2.h>
  26 #include <TCanvas.h>
  27 #include <TPaveLabel.h>
  28 #include <TVector3.h>
  29
  30 #include "AliD0toKpi.h"
  31
  32 ClassImp(AliD0toKpi)
  33
  34 //----------------------------------------------------------------------------
  35 AliD0toKpi::AliD0toKpi():
  36 fSignal(kFALSE),
  37 fEvent(0),
  38 fV1x(0.),
  39 fV1y(0.),
  40 fV1z(0.),
  41 fV2x(0.),
  42 fV2y(0.),
  43 fV2z(0.),
  44 fDCA(0.),
  45 fWgtAD0(0.),
  46 fWgtAD0bar(0.),
  47 fWgtBD0(0.),
  48 fWgtBD0bar(0.),
  49 fWgtCD0(0.),
  50 fWgtCD0bar(0.),
  51 fWgtDD0(0.),
  52 fWgtDD0bar(0.)
  53 {
  54   // Default constructor
  55
  56   fTrkNum[0] = 0;
  57   fTrkNum[1] = 0;
  58
  59   fPx[0] = 0.;
  60   fPy[0] = 0.;
  61   fPz[0] = 0.;
  62   fPx[1] = 0.;
  63   fPy[1] = 0.;
  64   fPz[1] = 0.;
  65
  66   fd0[0] = 0.;
  67   fd0[1] = 0.;
  68
  69   fPdg[0] = 0;
  70   fPdg[1] = 0;
  71   fMum[0] = 0;
  72   fMum[1] = 0;
  73
  74   fTagPi[0] = 0.;
  75   fTagPi[1] = 0.;
  76   fTagKa[0] = 0.;
  77   fTagKa[1] = 0.;
  78   fTagNid[0] = 0.;
  79   fTagNid[1] = 0.;
  80
  81 }
  82 //----------------------------------------------------------------------------
  83 AliD0toKpi::AliD0toKpi(Int_t ev,Int_t trkNum[2],
  84                        Double_t v1[3],Double_t v2[3],
  85                        Double_t dca,
  86                        Double_t mom[6],Double_t d0[2]):
  87 fSignal(kFALSE),
  88 fEvent(ev),
  89 fV1x(v1[0]),
  90 fV1y(v1[1]),
  91 fV1z(v1[2]),
  92 fV2x(v2[0]),
  93 fV2y(v2[1]),
  94 fV2z(v2[2]),
  95 fDCA(dca),
  96 fWgtAD0(0.),
  97 fWgtAD0bar(0.),
  98 fWgtBD0(0.),
  99 fWgtBD0bar(0.),
 100 fWgtCD0(0.),
 101 fWgtCD0bar(0.),
 102 fWgtDD0(0.),
 103 fWgtDD0bar(0.)
 104 {
 105   // Constructor
 106
 107   fTrkNum[0] = trkNum[0];
 108   fTrkNum[1] = trkNum[1];
 109
 110   fPx[0] = mom[0];
 111   fPy[0] = mom[1];
 112   fPz[0] = mom[2];
 113   fPx[1] = mom[3];
 114   fPy[1] = mom[4];
 115   fPz[1] = mom[5];
 116
 117   fd0[0] = d0[0];
 118   fd0[1] = d0[1];
 119
 120   fPdg[0] = 0;
 121   fPdg[1] = 0;
 122   fMum[0] = 0;
 123   fMum[1] = 0;
 124
 125   fTagPi[0]  = 0.;
 126   fTagPi[1]  = 0.;
 127   fTagKa[0]  = 0.;
 128   fTagKa[1]  = 0.;
 129   fTagNid[0] = 0.;
 130   fTagNid[1] = 0.;
 131 }
 132 //----------------------------------------------------------------------------
 133 AliD0toKpi::~AliD0toKpi() {}
 134 //----------------------------------------------------------------------------
 135 void AliD0toKpi::ApplyPID(TString pidScheme) {
 136   // Applies particle identification
 137
 138   if((!pidScheme.CompareTo("TOFparamPbPb") || !pidScheme.CompareTo("TOFparamPP")) && fPdg[0]==0) {
 139     printf("AliD0toKpi::ApplyPID :\n Warning: TOF parameterized PID can be used only for simulation!\n");
 140     return;
 141   }
 142
 143   if(!pidScheme.CompareTo("TOFparamPbPb")) {
 144     // tagging of the positive track
 145     if(TMath::Abs(fPdg[0])==211 || TMath::Abs(fPdg[0])==13
 146        || TMath::Abs(fPdg[0])==11) { // pion,muon,electron
 147       fTagPi[0]  = LinearInterpolation(PChild(0),kPiBinsPbPb,kPiBinWidthPbPb,kPiTagPiPbPb);
 148       fTagNid[0] = 1.-fTagPi[0];
 149       fTagKa[0]   = 0.;
 150       fTagPr[0]   = 0.;
 151     }
 152     if(TMath::Abs(fPdg[0])==321) { // kaon
 153       fTagKa[0]   = LinearInterpolation(PChild(0),kKBinsPbPb,kKBinWidthPbPb,kKTagKPbPb);
 154       fTagNid[0] = LinearInterpolation(PChild(0),kKBinsPbPb,kKBinWidthPbPb,kKTagNidPbPb);
 155       if((fTagNid[0]+fTagKa[0])>1.) fTagNid[0] = 1.-fTagKa[0];
 156       fTagPi[0] = 1.-fTagNid[0]-fTagKa[0];
 157       fTagPr[0] = 0.;
 158     }
 159     if(TMath::Abs(fPdg[0])==2212) { // proton
 160       fTagPr[0]  = LinearInterpolation(PChild(0),kPBinsPbPb,kPBinWidthPbPb,kPTagPPbPb);
 161       fTagNid[0] = LinearInterpolation(PChild(0),kPBinsPbPb,kPBinWidthPbPb,kPTagNidPbPb);
 162       if((fTagNid[0]+fTagPr[0])>1.) fTagNid[0] = 1.-fTagPr[0];
 163       fTagPi[0] = 1.-fTagNid[0]-fTagPr[0];
 164       fTagKa[0]   = 0.;
 165     }
 166     // tagging of the negative track
 167     if(TMath::Abs(fPdg[1])==211 || TMath::Abs(fPdg[1])==13
 168        || TMath::Abs(fPdg[1])==11) { // pion,muon,electron
 169       fTagPi[1]  = LinearInterpolation(PChild(1),kPiBinsPbPb,kPiBinWidthPbPb,kPiTagPiPbPb);
 170       fTagNid[1] = 1.-fTagPi[1];
 171       fTagKa[1]   = 0.;
 172       fTagPr[1]   = 0.;
 173     }
 174     if(TMath::Abs(fPdg[1])==321) { // kaon
 175       fTagKa[1]   = LinearInterpolation(PChild(1),kKBinsPbPb,kKBinWidthPbPb,kKTagKPbPb);
 176       fTagNid[1] = LinearInterpolation(PChild(1),kKBinsPbPb,kKBinWidthPbPb,kKTagNidPbPb);
 177       if((fTagNid[1]+fTagKa[1])>1.) fTagNid[1] = 1.-fTagKa[1];
 178       fTagPi[1] = 1.-fTagNid[1]-fTagKa[1];
 179       fTagPr[1] = 0.;
 180     }
 181     if(TMath::Abs(fPdg[1])==2212) { // proton
 182       fTagPr[1]  = LinearInterpolation(PChild(1),kPBinsPbPb,kPBinWidthPbPb,kPTagPPbPb);
 183       fTagNid[1] = LinearInterpolation(PChild(1),kPBinsPbPb,kPBinWidthPbPb,kPTagNidPbPb);
 184       if((fTagNid[1]+fTagPr[1])>1.) fTagNid[1] = 1.-fTagPr[1];
 185       fTagPi[1] = 1.-fTagNid[1]-fTagPr[1];
 186       fTagKa[1]   = 0.;
 187     }
 188   }
 189
 190
 191   if(!pidScheme.CompareTo("TOFparamPP")) {
 192     // tagging of the positive track
 193     if(TMath::Abs(fPdg[0])==211 || TMath::Abs(fPdg[0])==13
 194        || TMath::Abs(fPdg[0])==11) { // pion,muon,electron
 195       fTagPi[0]  = LinearInterpolation(PChild(0),kPiBinsPP,kPiBinWidthPP,kPiTagPiPP);
 196       fTagNid[0] = 1.-fTagPi[0];
 197       fTagKa[0]   = 0.;
 198       fTagPr[0]   = 0.;
 199     }
 200     if(TMath::Abs(fPdg[0])==321) { // kaon
 201       fTagKa[0]   = LinearInterpolation(PChild(0),kKBinsPP,kKBinWidthPP,kKTagKPP);
 202       fTagNid[0] = LinearInterpolation(PChild(0),kKBinsPP,kKBinWidthPP,kKTagNidPP);
 203       if((fTagNid[0]+fTagKa[0])>1.) fTagNid[0] = 1.-fTagKa[0];
 204       fTagPi[0] = 1.-fTagNid[0]-fTagKa[0];
 205       fTagPr[0] = 0.;
 206     }
 207     if(TMath::Abs(fPdg[0])==2212) { // proton
 208       fTagPr[0]  = LinearInterpolation(PChild(0),kPBinsPP,kPBinWidthPP,kPTagPPP);
 209       fTagNid[0] = LinearInterpolation(PChild(0),kPBinsPP,kPBinWidthPP,kPTagNidPP);
 210       if((fTagNid[0]+fTagPr[0])>1.) fTagNid[0] = 1.-fTagPr[0];
 211       fTagPi[0] = 1.-fTagNid[0]-fTagPr[0];
 212       fTagKa[0]   = 0.;
 213     }
 214     // tagging of the negative track
 215     if(TMath::Abs(fPdg[1])==211 || TMath::Abs(fPdg[1])==13
 216        || TMath::Abs(fPdg[1])==11) { // pion,muon,electron
 217       fTagPi[1]  = LinearInterpolation(PChild(1),kPiBinsPP,kPiBinWidthPP,kPiTagPiPP);
 218       fTagNid[1] = 1.-fTagPi[1];
 219       fTagKa[1]   = 0.;
 220       fTagPr[1]   = 0.;
 221     }
 222     if(TMath::Abs(fPdg[1])==321) { // kaon
 223       fTagKa[1]   = LinearInterpolation(PChild(1),kKBinsPP,kKBinWidthPP,kKTagKPP);
 224       fTagNid[1] = LinearInterpolation(PChild(1),kKBinsPP,kKBinWidthPP,kKTagNidPP);
 225       if((fTagNid[1]+fTagKa[1])>1.) fTagNid[1] = 1.-fTagKa[1];
 226       fTagPi[1] = 1.-fTagNid[1]-fTagKa[1];
 227       fTagPr[1] = 0.;
 228     }
 229     if(TMath::Abs(fPdg[1])==2212) { // proton
 230       fTagPr[1]  = LinearInterpolation(PChild(1),kPBinsPP,kPBinWidthPP,kPTagPPP);
 231       fTagNid[1] = LinearInterpolation(PChild(1),kPBinsPP,kPBinWidthPP,kPTagNidPP);
 232       if((fTagNid[1]+fTagPr[1])>1.) fTagNid[1] = 1.-fTagPr[1];
 233       fTagPi[1] = 1.-fTagNid[1]-fTagPr[1];
 234       fTagKa[1]   = 0.;
 235     }
 236   }
 237
 238   return;
 239 }
 240 //----------------------------------------------------------------------------
 241 Double_t AliD0toKpi::ChildrenRelAngle() const {
 242   // relative angle between K and pi
 243
 244   TVector3 mom0(fPx[0],fPy[0],fPz[0]);
 245   TVector3 mom1(fPx[1],fPy[1],fPz[1]);
 246
 247   Double_t angle = mom0.Angle(mom1);
 248
 249   return angle;
 250 }
 251 //----------------------------------------------------------------------------
 252 void AliD0toKpi::ComputeWgts() {
 253   // calculate the weights for PID
 254
 255
 256   // assignement of the weights from PID
 257   fWgtAD0    = fTagKa[1]*(fTagPi[0]+fTagNid[0]);
 258   fWgtAD0bar = fTagKa[0]*(fTagPi[1]+fTagNid[1]);
 259   fWgtBD0    = fTagPi[0]*fTagNid[1];
 260   fWgtBD0bar = fTagPi[1]*fTagNid[0];
 261   fWgtCD0    = fTagNid[0]*fTagNid[1];
 262   fWgtCD0bar = fTagNid[0]*fTagNid[1];
 263   fWgtDD0    = 1.-fWgtAD0-fWgtBD0-fWgtCD0;
 264   fWgtDD0bar = 1.-fWgtAD0bar-fWgtBD0bar-fWgtCD0bar;
 265
 266   /*
 267   cerr<<fWgtAD0<<"  "<<fWgtAD0bar<<endl;
 268   cerr<<fWgtBD0<<"  "<<fWgtBD0bar<<endl;
 269   cerr<<fWgtCD0<<"  "<<fWgtCD0bar<<endl;
 270   cerr<<fWgtDD0<<"  "<<fWgtDD0bar<<endl;
 271   */
 272   /*
 273   if(fWgtAD0<0.) cerr<<"AliD0toKpi::ComputeWgts()  Negative weight!!!\n";
 274   if(fWgtAD0bar<0.) cerr<<"AliD0toKpi::ComputeWgts()  Negative weight!!!\n";
 275   if(fWgtBD0<0.) cerr<<"AliD0toKpi::ComputeWgts()  Negative weight!!!\n";
 276   if(fWgtBD0bar<0.) cerr<<"AliD0toKpi::ComputeWgts()  Negative weight!!!\n";
 277   if(fWgtCD0<0.) cerr<<"AliD0toKpi::ComputeWgts()  Negative weight!!!\n";
 278   if(fWgtCD0bar<0.) cerr<<"AliD0toKpi::ComputeWgts()  Negative weight!!!\n";
 279   */
 280
 281   return;
 282 }
 283 //----------------------------------------------------------------------------
 284 void AliD0toKpi::CorrectWgt4BR(Double_t factor) {
 285   // correct weights of background from charm
 286
 287   fWgtAD0    *= factor;
 288   fWgtAD0bar *= factor;
 289   fWgtBD0    *= factor;
 290   fWgtBD0bar *= factor;
 291   fWgtCD0    *= factor;
 292   fWgtCD0bar *= factor;
 293   fWgtDD0    *= factor;
 294   fWgtDD0bar *= factor;
 295
 296   return;
 297 }
 298 //----------------------------------------------------------------------------
 299 Double_t AliD0toKpi::CosPointing() const {
 300   // cosine of pointing angle in space
 301
 302   TVector3 mom(Px(),Py(),Pz());
 303   TVector3 flight(fV2x-fV1x,fV2y-fV1y,fV2z-fV1z);
 304
 305   Double_t pta = mom.Angle(flight);
 306
 307   return TMath::Cos(pta);
 308 }
 309 //----------------------------------------------------------------------------
 310 Double_t AliD0toKpi::CosPointingXY() const {
 311   // cosine of pointing angle in transverse plane
 312
 313   TVector3 momXY(Px(),Py(),0.);
 314   TVector3 flightXY(fV2x-fV1x,fV2y-fV1y,0.);
 315
 316   Double_t ptaXY = momXY.Angle(flightXY);
 317
 318   return TMath::Cos(ptaXY);
 319 }
 320 //----------------------------------------------------------------------------
 321 void AliD0toKpi::CosThetaStar(Double_t &ctsD0,Double_t &ctsD0bar) const {
 322   // cosine of decay angle in the D0 rest frame
 323
 324   Double_t pStar = TMath::Sqrt(TMath::Power(kMD0*kMD0-kMK*kMK-kMPi*kMPi,2.)-4.*kMK*kMK*kMPi*kMPi)/(2.*kMD0);
 325
 326   Double_t beta = P()/Energy();
 327   Double_t gamma = Energy()/kMD0;
 328
 329   ctsD0 = (Ql(1)/gamma-beta*TMath::Sqrt(pStar*pStar+kMK*kMK))/pStar;
 330   // if(ctsD0 > 1.)  { cerr<<"AliD0toKpi::CosThetaStar: > 1 "<<ctsD0<<"!\n"; }
 331   // if(ctsD0 < -1.) { cerr<<"AliD0toKpi::CosThetaStar: < -1 "<<ctsD0<<"!\n"; }
 332
 333   ctsD0bar = (Ql(0)/gamma-beta*TMath::Sqrt(pStar*pStar+kMK*kMK))/pStar;
 334   // if(ctsD0bar > 1.)  { cerr<<"AliD0toKpi::CosThetaStar: > 1 "<<ctsD0bar<<"!\n"; }
 335   // if(ctsD0bar < -1.) { cerr<<"AliD0toKpi::CosThetaStar: < -1 "<<ctsD0bar<<"!\n";}
 336
 337   return;
 338 }
 339 //----------------------------------------------------------------------------
 340 Double_t AliD0toKpi::Eta() const {
 341   // pseudorapidity of the D0
 342
 343   Double_t theta = TMath::Pi()/2.-TMath::ATan2(Pz(),Pt());
 344   Double_t eta = -TMath::Log(TMath::Tan(theta/2.));
 345   return eta;
 346 }
 347 //----------------------------------------------------------------------------
 348 Double_t AliD0toKpi::EtaChild(Int_t child) const {
 349   // pseudorapidity of the decay tracks
 350
 351   Double_t theta = TMath::Pi()/2.-TMath::ATan2(fPz[child],PtChild(child));
 352   Double_t eta = -TMath::Log(TMath::Tan(theta/2.));
 353   return eta;
 354 }
 355 //----------------------------------------------------------------------------
 356 void AliD0toKpi::GetWgts(Double_t &WgtD0,Double_t &WgtD0bar,TString sample)
 357   const {
 358   // returns the weights for pid
 359
 360   if(!sample.CompareTo("A")) { WgtD0 = fWgtAD0;  WgtD0bar = fWgtAD0bar; }
 361   if(!sample.CompareTo("B")) { WgtD0 = fWgtBD0;  WgtD0bar = fWgtBD0bar; }
 362   if(!sample.CompareTo("C")) { WgtD0 = fWgtCD0;  WgtD0bar = fWgtCD0bar; }
 363   if(!sample.CompareTo("D")) { WgtD0 = fWgtDD0;  WgtD0bar = fWgtDD0bar; }
 364   if(!sample.CompareTo("ABCD")) {
 365     WgtD0    = fWgtAD0+fWgtBD0+fWgtCD0+fWgtDD0;
 366     WgtD0bar = fWgtAD0bar+fWgtBD0bar+fWgtCD0bar+fWgtDD0bar;
 367   }
 368   if(!sample.CompareTo("ABC")) {
 369     WgtD0    = fWgtAD0+fWgtBD0+fWgtCD0;
 370     WgtD0bar = fWgtAD0bar+fWgtBD0bar+fWgtCD0bar;
 371   }
 372   if(!sample.CompareTo("BC")) {
 373     WgtD0    = fWgtBD0+fWgtCD0;
 374     WgtD0bar = fWgtBD0bar+fWgtCD0bar;
 375   }
 376
 377   return;
 378 }
 379 //----------------------------------------------------------------------------
 380 Double_t AliD0toKpi::ImpPar() const {
 381   // D0 impact parameter in the bending plane
 382
 383     Double_t k = -(fV2x-fV1x)*Px()-(fV2y-fV1y)*Py();
 384     k /= Pt()*Pt();
 385     Double_t dx = fV2x-fV1x+k*Px();
 386     Double_t dy = fV2y-fV1y+k*Py();
 387     Double_t absDD = TMath::Sqrt(dx*dx+dy*dy);
 388     TVector3 mom(Px(),Py(),Pz());
 389     TVector3 flight(fV2x-fV1x,fV2y-fV1y,fV2z-fV1z);
 390     TVector3 cross = mom.Cross(flight);
 391     return (cross.Z()>0. ? absDD : -absDD);
 392 }
 393 //----------------------------------------------------------------------------
 394 void AliD0toKpi::InvMass(Double_t &mD0,Double_t &mD0bar) const {
 395   // invariant mass as D0 and as D0bar
 396
 397   Double_t energy[2];
 398
 399   // D0 -> K- Pi+
 400   energy[1] = TMath::Sqrt(kMK*kMK+PChild(1)*PChild(1));
 401   energy[0] = TMath::Sqrt(kMPi*kMPi+PChild(0)*PChild(0));
 402
 403   mD0 = TMath::Sqrt((energy[0]+energy[1])*(energy[0]+energy[1])-P()*P());
 404
 405
 406   // D0bar -> K+ Pi-
 407   energy[0] = TMath::Sqrt(kMK*kMK+PChild(0)*PChild(0));
 408   energy[1] = TMath::Sqrt(kMPi*kMPi+PChild(1)*PChild(1));
 409
 410   mD0bar = TMath::Sqrt((energy[0]+energy[1])*(energy[0]+energy[1])-P()*P());
 411
 412   return;
 413
 414 }
 415 //----------------------------------------------------------------------------
 416 Double_t AliD0toKpi::Ql(Int_t child) const {
 417   // longitudinal momentum of decay tracks w.r.t. to D0 momentum
 418
 419   Double_t qL;
 420   TVector3 mom(fPx[child],fPy[child],fPz[child]);
 421   TVector3 momD(Px(),Py(),Pz());
 422
 423   qL = mom.Dot(momD)/momD.Mag();
 424
 425   return qL ;
 426 }
 427 //----------------------------------------------------------------------------
 428 Double_t AliD0toKpi::Qt() const {
 429   // transverse momentum of decay tracks w.r.t. to D0 momentum
 430
 431   TVector3 mom0(fPx[0],fPy[0],fPz[0]);
 432   TVector3 momD(Px(),Py(),Pz());
 433
 434   return mom0.Perp(momD);
 435 }
 436 //----------------------------------------------------------------------------
 437 Bool_t AliD0toKpi::Select(const Double_t* cuts,Int_t& okD0,Int_t& okD0bar)
 438   const {
 439 //
 440 // This function compares the D0 with a set of cuts:
 441 //
 442 // cuts[0] = inv. mass half width [GeV]
 443 // cuts[1] = dca [micron]
 444 // cuts[2] = cosThetaStar
 445 // cuts[3] = pTK [GeV/c]
 446 // cuts[4] = pTPi [GeV/c]
 447 // cuts[5] = d0K [micron]   upper limit!
 448 // cuts[6] = d0Pi [micron]  upper limit!
 449 // cuts[7] = d0d0 [micron^2]
 450 // cuts[8] = cosThetaPoint
 451 //
 452 // If the D0/D0bar doesn't pass the cuts it sets the weights to 0
 453 // If neither D0 nor D0bar pass the cuts return kFALSE
 454 //
 455   Double_t mD0,mD0bar,ctsD0,ctsD0bar;
 456   okD0=1; okD0bar=1;
 457
 458   if(PtChild(1) < cuts[3] || PtChild(0) < cuts[4]) okD0 = 0;
 459   if(PtChild(0) < cuts[3] || PtChild(1) < cuts[4]) okD0bar = 0;
 460   if(!okD0 && !okD0bar) return kFALSE;
 461
 462   if(TMath::Abs(Getd0Child(1)) > cuts[5] ||
 463      TMath::Abs(Getd0Child(0)) > cuts[6]) okD0 = 0;
 464   if(TMath::Abs(Getd0Child(0)) > cuts[6] ||
 465      TMath::Abs(Getd0Child(1)) > cuts[5]) okD0bar = 0;
 466   if(!okD0 && !okD0bar) return kFALSE;
 467
 468   if(GetDCA() > cuts[1]) { okD0 = okD0bar = 0; return kFALSE; }
 469
 470   InvMass(mD0,mD0bar);
 471   if(TMath::Abs(mD0-kMD0)    > cuts[0]) okD0 = 0;
 472   if(TMath::Abs(mD0bar-kMD0) > cuts[0]) okD0bar = 0;
 473   if(!okD0 && !okD0bar) return kFALSE;
 474
 475   CosThetaStar(ctsD0,ctsD0bar);
 476   if(TMath::Abs(ctsD0)    > cuts[2]) okD0 = 0;
 477   if(TMath::Abs(ctsD0bar) > cuts[2]) okD0bar = 0;
 478   if(!okD0 && !okD0bar) return kFALSE;
 479
 480   if(ProdImpParams() > cuts[7]) { okD0 = okD0bar = 0; return kFALSE; }
 481
 482   if(CosPointing()   < cuts[8]) { okD0 = okD0bar = 0; return kFALSE; }
 483
 484   return kTRUE;
 485 }
 486 //-----------------------------------------------------------------------------
 487 void AliD0toKpi::SetPIDresponse(Double_t resp0[5],Double_t resp1[5]) {
 488   // Set combined PID detector response probabilities
 489
 490   fPIDrespEl[0] = resp0[0];
 491   fPIDrespEl[1] = resp1[0];
 492   fPIDrespMu[0] = resp0[1];
 493   fPIDrespMu[1] = resp1[1];
 494   fPIDrespPi[0] = resp0[2];
 495   fPIDrespPi[1] = resp1[2];
 496   fPIDrespKa[0] = resp0[3];
 497   fPIDrespKa[1] = resp1[3];
 498   fPIDrespPr[0] = resp0[4];
 499   fPIDrespPr[1] = resp1[4];
 500
 501   return;
 502 }
 503 //----------------------------------------------------------------------------
 504 Double_t AliD0toKpi::LinearInterpolation(Double_t p,Int_t nBins,Double_t Bin,
 505                                          const Double_t *values) const {
 506   // a linear interpolation method
 507
 508   Double_t value=0;
 509   Double_t slope;
 510
 511   if(p<0.5*Bin) {
 512     value = values[0];
 513   } else if(p>=(nBins-0.5)*Bin) {
 514     slope = (2*values[nBins-1]-values[nBins-2]-values[nBins-3])/Bin/2;
 515     value = values[nBins-2]+slope*(p-Bin*(nBins-1.5));
 516   } else {
 517     for(Int_t i=0; i<nBins; i++) {
 518       if(p<(i+0.5)*Bin) {
 519         slope = (values[i]-values[i-1])/Bin;
 520         value = values[i-1]+slope*(p-Bin*(i-0.5));
 521         break;
 522       }
 523     }
 524   }
 525
 526   if(value<0.) value=0.;
 527   if(value>1.) value=1.;
 528
 529   return value;
 530 }
 531 //----------------------------------------------------------------------------
 532
 533
 534
 535
 536
 537
 538
 539
 540