PWG3/vertexingHF/AliRDHFCutsDStartoKpipi.cxx

   1 /**************************************************************************
   2  * Copyright(c) 1998-2010, 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 /* $Id$ */
  17
  18 /////////////////////////////////////////////////////////////
  19 //
  20 // Class for cuts on AOD reconstructed DStar->Kpipi
  21 //
  22 // Author: A.Grelli, alessandro.grelli@uu.nl
  23 //
  24 // PID method implemented by   Y.Wang, yifei@physi.uni-heidelberg.de
  25 //
  26 /////////////////////////////////////////////////////////////
  27
  28 #include <TDatabasePDG.h>
  29 #include <Riostream.h>
  30
  31 #include "AliAODRecoDecayHF2Prong.h"
  32 #include "AliAODRecoCascadeHF.h"
  33 #include "AliRDHFCutsD0toKpi.h"
  34 #include "AliRDHFCutsDStartoKpipi.h"
  35 #include "AliAODTrack.h"
  36 #include "AliESDtrack.h"
  37 #include "AliAODPid.h"
  38 #include "AliTPCPIDResponse.h"
  39 #include "AliAODVertex.h"
  40 #include "AliESDVertex.h"
  41
  42 ClassImp(AliRDHFCutsDStartoKpipi)
  43
  44 //--------------------------------------------------------------------------
  45 AliRDHFCutsDStartoKpipi::AliRDHFCutsDStartoKpipi(const char* name) :
  46   AliRDHFCuts(name),
  47   fTrackCutsSoftPi(0)
  48 {
  49   //
  50   // Default Constructor
  51   //
  52
  53   Int_t nvars=14;
  54   SetNVars(nvars);
  55   TString varNames[14]={
  56     "inv. mass [GeV]",
  57     "dca [cm]",
  58     "cosThetaStar",
  59     "pTK [GeV/c]",
  60     "pTPi [GeV/c]",
  61     "d0K [cm]",
  62     "d0Pi [cm]",
  63     "d0d0 [cm^2]",
  64     "cosThetaPoint",
  65     "inv. mass half width of D* [GeV]",
  66     "half width of (M_Kpipi-M_D0) [GeV]",
  67     "PtMin of pi_s [GeV/c]",
  68     "PtMax of pi_s [GeV/c]",
  69     "theta, angle between the pi_s and decay plane of the D0 [rad]"};
  70   Bool_t isUpperCut[14]={
  71     kTRUE,
  72     kTRUE,
  73     kTRUE,
  74     kFALSE,
  75     kFALSE,
  76     kTRUE,
  77     kTRUE,
  78     kTRUE,
  79     kFALSE,
  80     kTRUE,
  81     kTRUE,
  82     kTRUE,
  83     kTRUE,
  84     kFALSE};
  85   SetVarNames(nvars,varNames,isUpperCut);
  86   Bool_t forOpt[14]={
  87     kFALSE,
  88     kTRUE,
  89     kTRUE,
  90     kFALSE,
  91     kFALSE,
  92     kFALSE,
  93     kFALSE,
  94     kTRUE,
  95     kTRUE,
  96     kFALSE,
  97     kTRUE,
  98     kFALSE,
  99     kFALSE,
 100     kFALSE};
 101   SetVarsForOpt(5,forOpt);
 102   Float_t limits[2]={0,999999999.};
 103   SetPtBins(2,limits);
 104 }
 105 //--------------------------------------------------------------------------
 106 AliRDHFCutsDStartoKpipi::AliRDHFCutsDStartoKpipi(const AliRDHFCutsDStartoKpipi &source) :
 107   AliRDHFCuts(source),
 108   fTrackCutsSoftPi(0)
 109 {
 110   //
 111   // Copy constructor
 112   //
 113
 114   if(source.GetTrackCutsSoftPi()) AddTrackCutsSoftPi(source.GetTrackCutsSoftPi());
 115
 116 }
 117 //--------------------------------------------------------------------------
 118 AliRDHFCutsDStartoKpipi &AliRDHFCutsDStartoKpipi::operator=(const AliRDHFCutsDStartoKpipi &source)
 119 {
 120   //
 121   // assignment operator
 122   //
 123   if(&source == this) return *this;
 124
 125   AliRDHFCuts::operator=(source);
 126   if(source.GetTrackCutsSoftPi()) AddTrackCutsSoftPi(source.GetTrackCutsSoftPi());
 127
 128   return *this;
 129 }
 130
 131
 132 //---------------------------------------------------------------------------
 133 void AliRDHFCutsDStartoKpipi::GetCutVarsForOpt(AliAODRecoDecayHF *d,Float_t *vars,Int_t nvars,Int_t *pdgdaughters) {
 134   //
 135   // Fills in vars the values of the variables
 136   //
 137   if(nvars!=fnVarsForOpt) {
 138     printf("AliRDHFCutsDStartoKpipi::GetCutsVarsForOpt: wrong number of variables\n");
 139     return;
 140   }
 141
 142
 143   AliAODRecoCascadeHF* dstarD0pi = (AliAODRecoCascadeHF*)d;
 144
 145   AliAODTrack *softPi = (AliAODTrack*)dstarD0pi->GetBachelor();
 146
 147   AliAODRecoDecayHF2Prong* dd = (AliAODRecoDecayHF2Prong*)dstarD0pi->Get2Prong();
 148
 149    Int_t iter=-1;
 150   if(fVarsForOpt[0]){
 151     iter++;
 152     if(TMath::Abs(pdgdaughters[0])==211) {
 153       vars[iter]=dd->InvMassD0();
 154     } else {
 155       vars[iter]=dd->InvMassD0bar();
 156     }
 157   }
 158   if(fVarsForOpt[1]){
 159     iter++;
 160     vars[iter]=dd->GetDCA();
 161   }
 162   if(fVarsForOpt[2]){
 163     iter++;
 164     if(TMath::Abs(pdgdaughters[0])==211) {
 165       vars[iter] = dd->CosThetaStarD0();
 166     } else {
 167       vars[iter] = dd->CosThetaStarD0bar();
 168     }
 169   }
 170   if(fVarsForOpt[3]){
 171     iter++;
 172    if(TMath::Abs(pdgdaughters[0])==321) {
 173      vars[iter]=dd->PtProng(0);
 174    }
 175    else{
 176      vars[iter]=dd->PtProng(1);
 177    }
 178   }
 179   if(fVarsForOpt[4]){
 180     iter++;
 181    if(TMath::Abs(pdgdaughters[0])==211) {
 182      vars[iter]=dd->PtProng(0);
 183    }
 184    else{
 185      vars[iter]=dd->PtProng(1);
 186    }
 187   }
 188   if(fVarsForOpt[5]){
 189     iter++;
 190     if(TMath::Abs(pdgdaughters[0])==321) {
 191      vars[iter]=dd->Getd0Prong(0);
 192    }
 193    else{
 194      vars[iter]=dd->Getd0Prong(1);
 195    }
 196   }
 197   if(fVarsForOpt[6]){
 198     iter++;
 199      if(TMath::Abs(pdgdaughters[0])==211) {
 200      vars[iter]=dd->Getd0Prong(0);
 201    }
 202    else{
 203      vars[iter]=dd->Getd0Prong(1);
 204    }
 205   }
 206   if(fVarsForOpt[7]){
 207     iter++;
 208     vars[iter]= dd->Prodd0d0();
 209   }
 210   if(fVarsForOpt[8]){
 211     iter++;
 212     vars[iter]=dd->CosPointingAngle();
 213   }
 214   if(fVarsForOpt[9]){
 215     iter++;
 216     vars[iter]=dstarD0pi->InvMassDstarKpipi();
 217   }
 218   if(fVarsForOpt[10]){
 219     iter++;
 220     vars[iter]=dstarD0pi->DeltaInvMass();
 221   }
 222   if(fVarsForOpt[11]){
 223     iter++;
 224     vars[iter] = softPi->Pt();
 225   }
 226   if(fVarsForOpt[12]){
 227     iter++;
 228     vars[iter] = softPi->Pt();
 229   }
 230   if(fVarsForOpt[13]){
 231     iter++;
 232     vars[iter] =dstarD0pi->AngleD0dkpPisoft();
 233   }
 234
 235   return;
 236 }
 237 //---------------------------------------------------------------------------
 238 Int_t AliRDHFCutsDStartoKpipi::IsSelected(TObject* obj,Int_t selectionLevel) {
 239   //
 240   // Apply selection for D*.
 241   //
 242   if(!fCutsRD){
 243     cout<<"Cut matrice not inizialized. Exit..."<<endl;
 244     return 0;
 245   }
 246
 247   AliAODRecoCascadeHF* d = (AliAODRecoCascadeHF*)obj;
 248   if(!d){
 249     cout<<"AliAODRecoCascadeHF null"<<endl;
 250     return 0;
 251   }
 252
 253   Double_t ptD=d->Pt();
 254   if(ptD<fMinPtCand) return 0;
 255   if(ptD>fMaxPtCand) return 0;
 256
 257   AliAODRecoDecayHF2Prong* dd = (AliAODRecoDecayHF2Prong*)d->Get2Prong();
 258   if(!dd){
 259     cout<<"AliAODRecoDecayHF2Prong null"<<endl;
 260     return 0;
 261   }
 262
 263   AliAODTrack *b = (AliAODTrack*)d->GetBachelor();
 264
 265
 266   Int_t returnvalue=1;
 267   Int_t returnvaluePID=3;
 268
 269
 270   // selection on candidate
 271   if(selectionLevel==AliRDHFCuts::kAll ||
 272      selectionLevel==AliRDHFCuts::kCandidate) {
 273
 274     Double_t pt=d->Pt();
 275     Int_t ptbin=PtBin(pt);
 276
 277     // DStarMass and D0mass
 278     Double_t mDSPDG = TDatabasePDG::Instance()->GetParticle(413)->Mass();
 279     Double_t mD0PDG = TDatabasePDG::Instance()->GetParticle(421)->Mass();
 280     // delta mass PDG
 281     Double_t deltaPDG = mDSPDG-mD0PDG;
 282
 283     // Half width DStar mass
 284     if(TMath::Abs(mDSPDG - (d->InvMassDstarKpipi()))>fCutsRD[GetGlobalIndex(9,ptbin)]) return 0;
 285     // Half width Delta mass
 286
 287     if(TMath::Abs(deltaPDG-(d->DeltaInvMass())) > fCutsRD[GetGlobalIndex(10,ptbin)]) return 0;
 288
 289     // cut on soft pion pt
 290     if(b->Pt() < fCutsRD[GetGlobalIndex(11,ptbin)] || b->Pt() > fCutsRD[GetGlobalIndex(12,ptbin)]) return 0;
 291     // cut on the angle between D0 decay plane and soft pion
 292     if(d->AngleD0dkpPisoft() > fCutsRD[GetGlobalIndex(13,ptbin)]) return 0;
 293
 294     // select D0 that passes D* cuts
 295     returnvalue = IsD0FromDStarSelected(pt,dd,selectionLevel);
 296     if((b->Charge()==+1 && returnvalue==2) || (b->Charge()==-1 && returnvalue==1)) return 0;
 297
 298   }
 299
 300   // selection on PID
 301   if(selectionLevel==AliRDHFCuts::kAll ||
 302      selectionLevel==AliRDHFCuts::kCandidate ||
 303      selectionLevel==AliRDHFCuts::kPID) {
 304     returnvaluePID = IsSelectedPID(d);
 305   }
 306   if(returnvaluePID!=3) returnvalue =0;
 307
 308
 309   // selection on daughter tracks
 310   if(selectionLevel==AliRDHFCuts::kAll ||
 311      selectionLevel==AliRDHFCuts::kTracks) {
 312     if(!AreDaughtersSelected(dd)) return 0;
 313     if(fTrackCutsSoftPi) {
 314       AliAODVertex *vAOD = d->GetPrimaryVtx();
 315       Double_t pos[3],cov[6];
 316       vAOD->GetXYZ(pos);
 317       vAOD->GetCovarianceMatrix(cov);
 318       const AliESDVertex vESD(pos,cov,100.,100);
 319       if(!IsDaughterSelected(b,&vESD,fTrackCutsSoftPi)) return 0;
 320     }
 321   }
 322
 323   return returnvalue;
 324
 325 }
 326 //_________________________________________________________________________________________________
 327 Int_t AliRDHFCutsDStartoKpipi::IsD0FromDStarSelected(Double_t pt, TObject* obj,Int_t selectionLevel) const {
 328   //
 329   // Apply selection for D0 from D*. The selection in on D0 prongs
 330   //
 331
 332   if(!fCutsRD){
 333     cout<<"Cut matrice not inizialized. Exit..."<<endl;
 334     return 0;
 335   }
 336
 337   AliAODRecoDecayHF2Prong* dd = (AliAODRecoDecayHF2Prong*)obj;
 338
 339   if(!dd){
 340     cout<<"AliAODRecoDecayHF2Prong null"<<endl;
 341     return 0;
 342   }
 343
 344   // selection on daughter tracks is done in IsSelected()
 345
 346   Int_t returnvalue=1;
 347
 348   // selection on candidate
 349   if(selectionLevel==AliRDHFCuts::kAll ||
 350      selectionLevel==AliRDHFCuts::kCandidate) {
 351
 352     // D0 mass
 353     Double_t mD0PDG = TDatabasePDG::Instance()->GetParticle(421)->Mass();
 354     // delta mass PDG
 355
 356     Int_t ptbin=PtBin(pt);
 357
 358     Double_t mD0,mD0bar,ctsD0,ctsD0bar;
 359
 360     Int_t okD0     =0;
 361     Int_t okD0bar  =0;
 362     okD0=1; okD0bar=1;
 363
 364     if(dd->PtProng(1) < fCutsRD[GetGlobalIndex(3,ptbin)] || dd->PtProng(0) < fCutsRD[GetGlobalIndex(4,ptbin)]) okD0 = 0;
 365     if(dd->PtProng(0) < fCutsRD[GetGlobalIndex(3,ptbin)] || dd->PtProng(1) < fCutsRD[GetGlobalIndex(4,ptbin)]) okD0bar = 0;
 366
 367     if(!okD0 && !okD0bar) return 0;
 368
 369     if(TMath::Abs(dd->Getd0Prong(1)) > fCutsRD[GetGlobalIndex(5,ptbin)] ||
 370        TMath::Abs(dd->Getd0Prong(0)) > fCutsRD[GetGlobalIndex(6,ptbin)]) okD0 = 0;
 371     if(TMath::Abs(dd->Getd0Prong(0)) > fCutsRD[GetGlobalIndex(6,ptbin)] ||
 372        TMath::Abs(dd->Getd0Prong(1)) > fCutsRD[GetGlobalIndex(5,ptbin)]) okD0bar = 0;
 373     if(!okD0 && !okD0bar) return 0;
 374
 375     if(dd->GetDCA() > fCutsRD[GetGlobalIndex(1,ptbin)]) return 0;
 376
 377     dd->InvMassD0(mD0,mD0bar);
 378     if(TMath::Abs(mD0-mD0PDG) > fCutsRD[GetGlobalIndex(0,ptbin)]) okD0 = 0;
 379     if(TMath::Abs(mD0bar-mD0PDG) > fCutsRD[GetGlobalIndex(0,ptbin)]) okD0bar = 0;
 380     if(!okD0 && !okD0bar) return 0;
 381
 382     dd->CosThetaStarD0(ctsD0,ctsD0bar);
 383     if(TMath::Abs(ctsD0) > fCutsRD[GetGlobalIndex(2,ptbin)]) okD0 = 0;
 384     if(TMath::Abs(ctsD0bar) > fCutsRD[GetGlobalIndex(2,ptbin)]) okD0bar = 0;
 385     if(!okD0 && !okD0bar) return 0;
 386
 387     if(dd->Prodd0d0() > fCutsRD[GetGlobalIndex(7,ptbin)]) return 0;
 388
 389     if(dd->CosPointingAngle() < fCutsRD[GetGlobalIndex(8,ptbin)]) return 0;
 390
 391     if (okD0) returnvalue=1; //cuts passed as D0
 392     if (okD0bar) returnvalue=2; //cuts passed as D0bar
 393     if (okD0 && okD0bar) returnvalue=3; //both
 394   }
 395
 396   return returnvalue;
 397 }
 398 //----------------------------------------------------------------------------------
 399 Bool_t AliRDHFCutsDStartoKpipi::IsInFiducialAcceptance(Double_t pt, Double_t y) const
 400 {
 401   //
 402   // D* fiducial acceptance region
 403   //
 404
 405   if(pt > 5.) {
 406     // applying cut for pt > 5 GeV
 407     AliDebug(4,Form("pt of D* = %f (> 5), cutting at |y| < 0.8\n",pt));
 408     if (TMath::Abs(y) > 0.8){
 409       return kFALSE;
 410     }
 411   } else {
 412     // appliying smooth cut for pt < 5 GeV
 413     Double_t maxFiducialY = -0.2/15*pt*pt+1.9/15*pt+0.5;
 414     Double_t minFiducialY = 0.2/15*pt*pt-1.9/15*pt-0.5;
 415     AliDebug(2,Form("pt of D* = %f (< 5), cutting  according to the fiducial zone [%f, %f]\n",pt,minFiducialY,maxFiducialY));
 416     if (y < minFiducialY || y > maxFiducialY){
 417       return kFALSE;
 418     }
 419   }
 420
 421   return kTRUE;
 422 }
 423
 424 //_______________________________________________________________________________-
 425 Int_t AliRDHFCutsDStartoKpipi::IsSelectedPID(AliAODRecoDecayHF* obj)
 426 {
 427   //
 428   // PID method, n signa approach default
 429   //
 430
 431   if(!fUsePID) return 3;
 432
 433   AliAODRecoCascadeHF* dstar = (AliAODRecoCascadeHF*)obj;
 434   if(!dstar){
 435     cout<<"AliAODRecoCascadeHF null"<<endl;
 436     return 0;
 437   }
 438   AliAODRecoDecayHF2Prong* d0 = (AliAODRecoDecayHF2Prong*)dstar->Get2Prong();
 439   if(!d0){
 440     cout<<"AliAODRecoDecayHF2Prong null"<<endl;
 441     return 0;
 442   }
 443
 444   //  here the PID
 445   AliAODTrack *pos = (AliAODTrack*)dstar->Get2Prong()->GetDaughter(0);
 446   AliAODTrack *neg = (AliAODTrack*)dstar->Get2Prong()->GetDaughter(1);
 447
 448   if (dstar->Charge()>0){
 449     if(!SelectPID(pos,2)) return 0;//pion+
 450     if(!SelectPID(neg,3)) return 0;//kaon-
 451   }else{
 452     if(!SelectPID(pos,3)) return 0;//kaon+
 453     if(!SelectPID(neg,2)) return 0;//pion-
 454   }
 455
 456   return 3;
 457 }
 458
 459 //_______________________________________________________________________________-
 460 Int_t AliRDHFCutsDStartoKpipi::SelectPID(AliAODTrack *track, Int_t type)
 461 {
 462   //
 463   //  here the PID
 464
 465   Bool_t isParticle=kTRUE;
 466
 467   if(fPidHF->GetMatch()==1){//n-sigma
 468     Bool_t TPCon=TMath::Abs(2)>1e-4?kTRUE:kFALSE;
 469     Bool_t TOFon=TMath::Abs(3)>1e-4?kTRUE:kFALSE;
 470
 471     Bool_t isTPC=kTRUE;
 472     Bool_t isTOF=kTRUE;
 473
 474     if (TPCon){//TPC
 475       if(fPidHF->CheckStatus(track,"TPC")){
 476         if(type==2) isTPC=fPidHF->IsPionRaw(track,"TPC");
 477         if(type==3) isTPC=fPidHF->IsKaonRaw(track,"TPC");
 478       }
 479     }
 480     if (TOFon){//TOF
 481       if(fPidHF->CheckStatus(track,"TOF")){
 482         if(type==2) isTOF=fPidHF->IsPionRaw(track,"TOF");
 483         if(type==3) isTOF=fPidHF->IsKaonRaw(track,"TOF");
 484       }
 485     }
 486     isParticle = isTPC&&isTOF;
 487   }
 488
 489   if(fPidHF->GetMatch()==2){//bayesian
 490     //Double_t priors[5]={0.01,0.001,0.3,0.3,0.3};
 491     Double_t prob[5]={1.,1.,1.,1.,1.};
 492
 493     //fPidHF->SetPriors(priors);
 494     fPidHF->BayesianProbability(track,prob);
 495
 496     Double_t max=0.;
 497     Int_t k=-1;
 498     for (Int_t i=0; i<5; i++) {
 499       if (prob[i]>max) {k=i; max=prob[i];}
 500     }
 501     isParticle = Bool_t(k==type);
 502   }
 503
 504   return isParticle;
 505
 506 }
 507 //__________________________________________________________________________________-
 508 void  AliRDHFCutsDStartoKpipi::SetStandardCutsPP2010() {
 509   //
 510   //STANDARD CUTS USED FOR 2010 pp analysis
 511   //
 512   // Need to be updated for the final cut version
 513   //
 514
 515   SetName("DStartoD0piCutsStandard");
 516   SetTitle("Standard Cuts for D* analysis");
 517
 518   // PILE UP REJECTION
 519   SetOptPileup(AliRDHFCuts::kRejectPileupEvent);
 520
 521   // EVENT CUTS
 522   SetMinVtxContr(1);
 523
 524   // CUTS ON SINGLE TRACKS
 525   AliESDtrackCuts *esdTrackCuts = new AliESDtrackCuts("AliESDtrackCuts","default");
 526   esdTrackCuts->SetRequireSigmaToVertex(kFALSE);
 527   esdTrackCuts->SetRequireTPCRefit(kTRUE);
 528   esdTrackCuts->SetRequireITSRefit(kTRUE);
 529   esdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,AliESDtrackCuts::kAny);
 530   esdTrackCuts->SetMinDCAToVertexXY(0.);
 531   esdTrackCuts->SetEtaRange(-0.8,0.8);
 532   esdTrackCuts->SetPtRange(0.3,1.e10);
 533
 534   // CUTS on SOFT PION
 535   AliESDtrackCuts* esdSoftPicuts=new AliESDtrackCuts();
 536   esdSoftPicuts->SetRequireSigmaToVertex(kFALSE);
 537   esdSoftPicuts->SetRequireTPCRefit(kFALSE);
 538   esdSoftPicuts->SetRequireITSRefit(kFALSE);
 539   esdSoftPicuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,
 540                                           AliESDtrackCuts::kAny);
 541   esdSoftPicuts->SetPtRange(0.0,1.e10);
 542
 543   AddTrackCuts(esdTrackCuts);
 544   AddTrackCutsSoftPi(esdSoftPicuts);
 545
 546   const Int_t nptbins =13;
 547   const Double_t ptmax = 9999.;
 548   const Int_t nvars=14;
 549   Float_t ptbins[nptbins+1];
 550   ptbins[0]=0.;
 551   ptbins[1]=0.5;
 552   ptbins[2]=1.;
 553   ptbins[3]=2.;
 554   ptbins[4]=3.;
 555   ptbins[5]=4.;
 556   ptbins[6]=5.;
 557   ptbins[7]=6.;
 558   ptbins[8]=8.;
 559   ptbins[9]=12.;
 560   ptbins[10]=16.;
 561   ptbins[11]=20.;
 562   ptbins[12]=24.;
 563   ptbins[13]=ptmax;
 564
 565   SetGlobalIndex(nvars,nptbins);
 566   SetPtBins(nptbins+1,ptbins);
 567
 568   Float_t cutsMatrixD0toKpiStand[nptbins][nvars]={{0.7,220.*1E-4,0.7,0.21,0.21,500.*1E-4,500.*1E-4,-2000.*1E-8,0.85,0.3,0.1,0.05,100,0.5},/* pt<0.5*/
 569                                                   {0.7,220.*1E-4,0.7,0.21,0.21,500.*1E-4,500.*1E-4,-16000.*1E-8,0.85,0.3,0.1,0.05,100,0.5},/* 0.5<pt<1*/
 570                                                   {0.7,400.*1E-4,0.8,0.7,0.7,400.*1E-4,400.*1E-4,-36000.*1E-8,0.82,0.3,0.1,0.05,100,0.5},/* 1<pt<2 */
 571                                                   {0.7,200.*1E-4,0.8,0.7,0.7,800.*1E-4,800.*1E-4,-16000.*1E-8,0.9,0.3,0.1,0.05,100,0.5},/* 2<pt<3 */
 572                                                   {0.7,500.*1E-4,0.8,1.0,1.0,420.*1E-4,560.*1E-4,-6500.*1E-8,0.9,0.3,0.1,0.05,100,0.5},/* 3<pt<4 */
 573                                                   {0.7,800.*1E-4,0.9,1.2,1.2,700.*1E-4,700.*1E-4,1000.*1E-8,0.9,0.3,0.1,0.05,100,0.5},/* 4<pt<5 */
 574                                                   {0.7,1000.*1E-4,1.0,1.0,1.0,800.*1E-4,800.*1E-4,50000.*1E-8,0.8,0.3,0.1,0.05,100,0.5},/* 5<pt<6 */
 575                                                   {0.7,1000.*1E-4,1.0,1.0,1.0,1000.*1E-4,1000.*1E-4,100000.*1E-8,0.7,0.3,0.1,0.05,100,0.5},/* 6<pt<8 */
 576                                                   {0.7,1000.*1E-4,1.0,1.0,1.0,1000.*1E-4,1000.*1E-4,600000.*1E-8,0.7,0.3,0.1,0.05,100,0.5},/* 8<pt<12 */
 577                                                   {0.7,1000.*1E-4,1.0,1.0,1.0,1500.*1E-4,1500.*1E-4,1000000.*1E-8,0.7,0.3,0.1,0.05,100,0.5},/* 12<pt<16 */
 578                                                   {0.7,1000.*1E-4,1.0,1.0,1.0,1500.*1E-4,1500.*1E-4,1000000.*1E-8,0.7,0.3,0.1,0.05,100,0.5},/* 16<pt<20 */
 579                                                   {0.7,1000.*1E-4,1.0,1.0,1.0,1500.*1E-4,1500.*1E-4,1000000.*1E-8,0.7,0.3,0.1,0.05,100,0.5},/* 20<pt<24 */
 580                                                   {0.7,1000.*1E-4,1.0,1.0,1.0,1500.*1E-4,1500.*1E-4,1000000.*1E-8,0.7,0.3,0.1,0.05,100,0.5}};/* pt>24 */
 581
 582
 583   //CREATE TRANSPOSE MATRIX...REVERSE INDICES as required by AliRDHFCuts
 584   Float_t **cutsMatrixTransposeStand=new Float_t*[nvars];
 585   for(Int_t iv=0;iv<nvars;iv++)cutsMatrixTransposeStand[iv]=new Float_t[nptbins];
 586
 587   for (Int_t ibin=0;ibin<nptbins;ibin++){
 588     for (Int_t ivar = 0; ivar<nvars; ivar++){
 589       cutsMatrixTransposeStand[ivar][ibin]=cutsMatrixD0toKpiStand[ibin][ivar];
 590     }
 591   }
 592
 593   SetCuts(nvars,nptbins,cutsMatrixTransposeStand);
 594
 595   for(Int_t iv=0;iv<nvars;iv++) delete [] cutsMatrixTransposeStand[iv];
 596   delete [] cutsMatrixTransposeStand;
 597   cutsMatrixTransposeStand=NULL;
 598
 599   // PID SETTINGS FOR D* analysis
 600   AliAODPidHF* pidObj=new AliAODPidHF();
 601   //pidObj->SetName("pid4DSatr");
 602   Int_t mode=1;
 603   Double_t priors[5]={0.01,0.001,0.3,0.3,0.3};
 604   pidObj->SetPriors(priors);
 605   pidObj->SetMatch(mode);
 606   pidObj->SetSigma(0,2); // TPC
 607   pidObj->SetSigma(3,3); // TOF
 608   pidObj->SetTPC(kTRUE);
 609   pidObj->SetTOF(kTRUE);
 610
 611   SetPidHF(pidObj);
 612   SetUsePID(kTRUE);
 613
 614   PrintAll();
 615
 616   delete pidObj;
 617   pidObj=NULL;
 618
 619   return;
 620 }
 621 //_____________________________________________________________________________-
 622 void  AliRDHFCutsDStartoKpipi::SetStandardCutsPbPb2010(){
 623   //
 624   // TEMPORARY, WORK IN PROGRESS ... BUT WORKING!
 625   //
 626   //  Lead Lead
 627   //
 628
 629   SetName("DStartoD0piCutsStandard");
 630   SetTitle("Standard Cuts for D* analysis in PbPb 2010");
 631
 632   // EVENT CUTS
 633   SetMinVtxContr(1);
 634
 635   // CUTS ON SINGLE TRACKS
 636   AliESDtrackCuts *esdTrackCuts = new AliESDtrackCuts("AliESDtrackCuts","default");
 637   esdTrackCuts->SetRequireSigmaToVertex(kFALSE);
 638   esdTrackCuts->SetRequireTPCRefit(kTRUE);
 639   esdTrackCuts->SetRequireITSRefit(kTRUE);
 640   esdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,AliESDtrackCuts::kAny);
 641   esdTrackCuts->SetMinDCAToVertexXY(0.);
 642   esdTrackCuts->SetEtaRange(-0.8,0.8);
 643   esdTrackCuts->SetPtRange(0.3,1.e10);
 644
 645   // CUTS on SOFT PION
 646   AliESDtrackCuts* esdSoftPicuts=new AliESDtrackCuts();
 647   esdSoftPicuts->SetRequireSigmaToVertex(kFALSE);
 648   esdSoftPicuts->SetRequireTPCRefit(kTRUE);
 649   esdSoftPicuts->SetRequireITSRefit(kTRUE);
 650   esdSoftPicuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,
 651                                           AliESDtrackCuts::kAny); //test d0 asimmetry
 652   esdSoftPicuts->SetPtRange(0.25,5);
 653
 654   AddTrackCuts(esdTrackCuts);
 655   AddTrackCutsSoftPi(esdSoftPicuts);
 656
 657   const Int_t nptbins =13;
 658   const Double_t ptmax = 9999.;
 659   const Int_t nvars=14;
 660   Float_t ptbins[nptbins+1];
 661   ptbins[0]=0.;
 662   ptbins[1]=0.5;
 663   ptbins[2]=1.;
 664   ptbins[3]=2.;
 665   ptbins[4]=3.;
 666   ptbins[5]=4.;
 667   ptbins[6]=5.;
 668   ptbins[7]=6.;
 669   ptbins[8]=8.;
 670   ptbins[9]=12.;
 671   ptbins[10]=16.;
 672   ptbins[11]=20.;
 673   ptbins[12]=24.;
 674   ptbins[13]=ptmax;
 675
 676   SetGlobalIndex(nvars,nptbins);
 677   SetPtBins(nptbins+1,ptbins);
 678
 679   Float_t cutsMatrixD0toKpiStand[nptbins][nvars]={{0.7,220.*1E-4,0.7,0.21,0.21,500.*1E-4,500.*1E-4,-2000.*1E-8,0.85,0.3,0.1,0.05,100,0.5},/* pt<0.5*/
 680                                                   {0.7,220.*1E-4,0.7,0.21,0.21,500.*1E-4,500.*1E-4,-16000.*1E-8,0.85,0.3,0.1,0.05,100,0.5},/* 0.5<pt<1*/
 681                                                   {0.7,400.*1E-4,0.8,0.7,0.7,800.*1E-4,800.*1E-4,-36000.*1E-8,0.82,0.3,0.1,0.05,100,0.5},/* 1<pt<2 */
 682                                                   {0.7,200.*1E-4,0.8,0.7,0.7,800.*1E-4,800.*1E-4,-16000.*1E-8,0.9,0.3,0.1,0.05,100,0.5},/* 2<pt<3 */
 683                                                   {0.7,500.*1E-4,0.8,1.0,1.0,420.*1E-4,560.*1E-4,-6500.*1E-8,0.9,0.3,0.1,0.05,100,0.5},/* 3<pt<4 */
 684                                                   {0.7,800.*1E-4,0.9,1.2,1.2,700.*1E-4,700.*1E-4,1000.*1E-8,0.9,0.3,0.1,0.05,100,0.5},/* 4<pt<5 */
 685                                                   {0.7,1000.*1E-4,1.0,1.0,1.0,800.*1E-4,800.*1E-4,50000.*1E-8,0.8,0.3,0.1,0.05,100,0.5},/* 5<pt<6 */
 686                                                   {0.7,1000.*1E-4,1.0,1.0,1.0,1000.*1E-4,1000.*1E-4,100000.*1E-8,0.7,0.3,0.1,0.05,100,0.5},/* 6<pt<8 */
 687                                                   {0.7,1000.*1E-4,1.0,1.0,1.0,1000.*1E-4,1000.*1E-4,600000.*1E-8,0.7,0.3,0.1,0.05,100,0.5},/* 8<pt<12 */
 688                                                   {0.7,1000.*1E-4,1.0,1.0,1.0,1500.*1E-4,1500.*1E-4,1000000.*1E-8,0.7,0.3,0.1,0.05,100,0.5},/* 12<pt<16 */
 689                                                   {0.7,1000.*1E-4,1.0,1.0,1.0,1500.*1E-4,1500.*1E-4,1000000.*1E-8,0.7,0.3,0.1,0.05,100,0.5},/* 16<pt<20 */
 690                                                   {0.7,1000.*1E-4,1.0,1.0,1.0,1500.*1E-4,1500.*1E-4,1000000.*1E-8,0.7,0.3,0.1,0.05,100,0.5},/* 20<pt<24 */
 691                                                   {0.7,1000.*1E-4,1.0,1.0,1.0,1500.*1E-4,1500.*1E-4,1000000.*1E-8,0.7,0.3,0.1,0.05,100,0.5}};/* pt>24 */
 692
 693
 694   //CREATE TRANSPOSE MATRIX...REVERSE INDICES as required by AliRDHFCuts
 695   Float_t **cutsMatrixTransposeStand=new Float_t*[nvars];
 696   for(Int_t iv=0;iv<nvars;iv++)cutsMatrixTransposeStand[iv]=new Float_t[nptbins];
 697
 698   for (Int_t ibin=0;ibin<nptbins;ibin++){
 699     for (Int_t ivar = 0; ivar<nvars; ivar++){
 700       cutsMatrixTransposeStand[ivar][ibin]=cutsMatrixD0toKpiStand[ibin][ivar];
 701     }
 702   }
 703
 704   SetCuts(nvars,nptbins,cutsMatrixTransposeStand);
 705
 706   for(Int_t iv=0;iv<nvars;iv++) delete [] cutsMatrixTransposeStand[iv];
 707   delete [] cutsMatrixTransposeStand;
 708   cutsMatrixTransposeStand=NULL;
 709
 710   // PID SETTINGS
 711   AliAODPidHF* pidObj=new AliAODPidHF();
 712   // pidObj->SetName("pid4DSatr");
 713   Int_t mode=1;
 714   Double_t priors[5]={0.01,0.001,0.3,0.3,0.3};
 715   pidObj->SetPriors(priors);
 716   pidObj->SetMatch(mode);
 717   pidObj->SetSigma(0,2); // TPC
 718   pidObj->SetSigma(3,3); // TOF
 719   pidObj->SetTPC(kTRUE);
 720   pidObj->SetTOF(kTRUE);
 721
 722   SetPidHF(pidObj);
 723   SetUsePID(kTRUE);
 724
 725   PrintAll();
 726
 727   delete pidObj;
 728   pidObj=NULL;
 729
 730   return;
 731
 732 }