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 /////////////////////////////////////////////////////////////
  17 //
  18 // Class for cuts on AOD reconstructed DStar->Kpipi
  19 //
  20 // Author: A.Grelli, alessandro.grelli@uu.nl
  21 /////////////////////////////////////////////////////////////
  22
  23 #include <TDatabasePDG.h>
  24 #include <Riostream.h>
  25
  26 #include "AliAODRecoDecayHF2Prong.h"
  27 #include "AliAODRecoCascadeHF.h"
  28 #include "AliRDHFCutsD0toKpi.h"
  29 #include "AliRDHFCutsDStartoKpipi.h"
  30 #include "AliAODTrack.h"
  31 #include "AliESDtrack.h"
  32 #include "AliAODPid.h"
  33 #include "AliTPCPIDResponse.h"
  34 #include "AliAODVertex.h"
  35 #include "AliESDVertex.h"
  36
  37
  38
  39 ClassImp(AliRDHFCutsDStartoKpipi)
  40
  41 //--------------------------------------------------------------------------
  42 AliRDHFCutsDStartoKpipi::AliRDHFCutsDStartoKpipi(const char* name) :
  43   AliRDHFCuts(name),
  44   fTrackCutsSoftPi(0)
  45 {
  46   //
  47   // Default Constructor
  48   //
  49
  50   Int_t nvars=14;
  51   SetNVars(nvars);
  52   TString varNames[14]={
  53     "inv. mass [GeV]",
  54     "dca [cm]",
  55     "cosThetaStar",
  56     "pTK [GeV/c]",
  57     "pTPi [GeV/c]",
  58     "d0K [cm]",
  59     "d0Pi [cm]",
  60     "d0d0 [cm^2]",
  61     "cosThetaPoint",
  62     "inv. mass half width of D* [GeV]",
  63     "half width of (M_Kpipi-M_D0) [GeV]",
  64     "PtMin of pi_s [GeV/c]",
  65     "PtMax of pi_s [GeV/c]",
  66     "theta, angle between the pi_s and decay plane of the D0 [rad]"};
  67   Bool_t isUpperCut[14]={
  68     kTRUE,
  69     kTRUE,
  70     kTRUE,
  71     kFALSE,
  72     kFALSE,
  73     kTRUE,
  74     kTRUE,
  75     kTRUE,
  76     kFALSE,
  77     kTRUE,
  78     kTRUE,
  79     kTRUE,
  80     kTRUE,
  81     kFALSE};
  82   SetVarNames(nvars,varNames,isUpperCut);
  83   Bool_t forOpt[14]={
  84     kFALSE,
  85     kTRUE,
  86     kTRUE,
  87     kFALSE,
  88     kFALSE,
  89     kFALSE,
  90     kFALSE,
  91     kTRUE,
  92     kTRUE,
  93     kFALSE,
  94     kTRUE,
  95     kFALSE,
  96     kFALSE,
  97     kFALSE};
  98   SetVarsForOpt(5,forOpt);
  99   Float_t limits[2]={0,999999999.};
 100   SetPtBins(2,limits);
 101 }
 102 //--------------------------------------------------------------------------
 103 AliRDHFCutsDStartoKpipi::AliRDHFCutsDStartoKpipi(const AliRDHFCutsDStartoKpipi &source) :
 104   AliRDHFCuts(source),
 105   fTrackCutsSoftPi(0)
 106 {
 107   //
 108   // Copy constructor
 109   //
 110
 111   if(source.GetTrackCutsSoftPi()) AddTrackCutsSoftPi(source.GetTrackCutsSoftPi());
 112
 113 }
 114 //--------------------------------------------------------------------------
 115 AliRDHFCutsDStartoKpipi &AliRDHFCutsDStartoKpipi::operator=(const AliRDHFCutsDStartoKpipi &source)
 116 {
 117   //
 118   // assignment operator
 119   //
 120   if(&source == this) return *this;
 121
 122   AliRDHFCuts::operator=(source);
 123   if(source.GetTrackCutsSoftPi()) AddTrackCutsSoftPi(source.GetTrackCutsSoftPi());
 124
 125   return *this;
 126 }
 127
 128
 129 //---------------------------------------------------------------------------
 130 void AliRDHFCutsDStartoKpipi::GetCutVarsForOpt(AliAODRecoDecayHF *d,Float_t *vars,Int_t nvars,Int_t *pdgdaughters) {
 131   //
 132   // Fills in vars the values of the variables
 133   //
 134   if(nvars!=fnVarsForOpt) {
 135     printf("AliRDHFCutsDStartoKpipi::GetCutsVarsForOpt: wrong number of variables\n");
 136     return;
 137   }
 138
 139
 140   AliAODRecoCascadeHF* dstarD0pi = (AliAODRecoCascadeHF*)d;
 141
 142   AliAODTrack *softPi = (AliAODTrack*)dstarD0pi->GetBachelor();
 143
 144   AliAODRecoDecayHF2Prong* dd = (AliAODRecoDecayHF2Prong*)dstarD0pi->Get2Prong();
 145
 146    Int_t iter=-1;
 147   if(fVarsForOpt[0]){
 148     iter++;
 149     if(TMath::Abs(pdgdaughters[0])==211) {
 150       vars[iter]=dd->InvMassD0();
 151     } else {
 152       vars[iter]=dd->InvMassD0bar();
 153     }
 154   }
 155   if(fVarsForOpt[1]){
 156     iter++;
 157     vars[iter]=dd->GetDCA();
 158   }
 159   if(fVarsForOpt[2]){
 160     iter++;
 161     if(TMath::Abs(pdgdaughters[0])==211) {
 162       vars[iter] = dd->CosThetaStarD0();
 163     } else {
 164       vars[iter] = dd->CosThetaStarD0bar();
 165     }
 166   }
 167   if(fVarsForOpt[3]){
 168     iter++;
 169    if(TMath::Abs(pdgdaughters[0])==321) {
 170      vars[iter]=dd->PtProng(0);
 171    }
 172    else{
 173      vars[iter]=dd->PtProng(1);
 174    }
 175   }
 176   if(fVarsForOpt[4]){
 177     iter++;
 178    if(TMath::Abs(pdgdaughters[0])==211) {
 179      vars[iter]=dd->PtProng(0);
 180    }
 181    else{
 182      vars[iter]=dd->PtProng(1);
 183    }
 184   }
 185   if(fVarsForOpt[5]){
 186     iter++;
 187     if(TMath::Abs(pdgdaughters[0])==321) {
 188      vars[iter]=dd->Getd0Prong(0);
 189    }
 190    else{
 191      vars[iter]=dd->Getd0Prong(1);
 192    }
 193   }
 194   if(fVarsForOpt[6]){
 195     iter++;
 196      if(TMath::Abs(pdgdaughters[0])==211) {
 197      vars[iter]=dd->Getd0Prong(0);
 198    }
 199    else{
 200      vars[iter]=dd->Getd0Prong(1);
 201    }
 202   }
 203   if(fVarsForOpt[7]){
 204     iter++;
 205     vars[iter]= dd->Prodd0d0();
 206   }
 207   if(fVarsForOpt[8]){
 208     iter++;
 209     vars[iter]=dd->CosPointingAngle();
 210   }
 211   if(fVarsForOpt[9]){
 212     iter++;
 213     vars[iter]=dstarD0pi->InvMassDstarKpipi();
 214   }
 215   if(fVarsForOpt[10]){
 216     iter++;
 217     vars[iter]=dstarD0pi->DeltaInvMass();
 218   }
 219   if(fVarsForOpt[11]){
 220     iter++;
 221     vars[iter] = softPi->Pt();
 222   }
 223   if(fVarsForOpt[12]){
 224     iter++;
 225     vars[iter] = softPi->Pt();
 226   }
 227   if(fVarsForOpt[13]){
 228     iter++;
 229     vars[iter] =dstarD0pi->AngleD0dkpPisoft();
 230   }
 231
 232   return;
 233 }
 234 //---------------------------------------------------------------------------
 235 Int_t AliRDHFCutsDStartoKpipi::IsSelected(TObject* obj,Int_t selectionLevel) {
 236   //
 237   // Apply selection for D*.
 238   //
 239   if(!fCutsRD){
 240     cout<<"Cut matrice not inizialized. Exit..."<<endl;
 241     return 0;
 242   }
 243
 244   AliAODRecoCascadeHF* d = (AliAODRecoCascadeHF*)obj;
 245   if(!d){
 246     cout<<"AliAODRecoCascadeHF null"<<endl;
 247     return 0;
 248   }
 249
 250   AliAODRecoDecayHF2Prong* dd = (AliAODRecoDecayHF2Prong*)d->Get2Prong();
 251   if(!dd){
 252     cout<<"AliAODRecoDecayHF2Prong null"<<endl;
 253     return 0;
 254   }
 255
 256   AliAODTrack *b = (AliAODTrack*)d->GetBachelor();
 257
 258   // selection on daughter tracks
 259   if(selectionLevel==AliRDHFCuts::kAll ||
 260      selectionLevel==AliRDHFCuts::kTracks) {
 261     if(!AreDaughtersSelected(dd)) return 0;
 262     if(fTrackCutsSoftPi) {
 263       AliAODVertex *vAOD = d->GetPrimaryVtx();
 264       Double_t pos[3],cov[6];
 265       vAOD->GetXYZ(pos);
 266       vAOD->GetCovarianceMatrix(cov);
 267       const AliESDVertex vESD(pos,cov,100.,100);
 268       if(!IsDaughterSelected(b,&vESD,fTrackCutsSoftPi)) return 0;
 269     }
 270   }
 271
 272   Int_t returnvalue=1;
 273   Int_t returnvaluePID=3;
 274
 275   // selection on PID
 276   if(selectionLevel==AliRDHFCuts::kAll ||
 277      selectionLevel==AliRDHFCuts::kCandidate ||
 278      selectionLevel==AliRDHFCuts::kPID) {
 279     returnvaluePID = IsSelectedPID(dd);
 280   }
 281
 282
 283   // selection on candidate
 284   if(selectionLevel==AliRDHFCuts::kAll ||
 285      selectionLevel==AliRDHFCuts::kCandidate) {
 286
 287     Double_t pt=d->Pt();
 288     Int_t ptbin=PtBin(pt);
 289
 290     // select D0 that passes D* cuts
 291     returnvalue = IsD0FromDStarSelected(pt,dd,selectionLevel);
 292
 293     if((b->Charge()==+1 && returnvalue==2) || (b->Charge()==-1 && returnvalue==1)) return 0;
 294
 295     // DStarMass and D0mass
 296     Double_t mDSPDG = TDatabasePDG::Instance()->GetParticle(413)->Mass();
 297     Double_t mD0PDG = TDatabasePDG::Instance()->GetParticle(421)->Mass();
 298     // delta mass PDG
 299     Double_t deltaPDG = mDSPDG-mD0PDG;
 300
 301     // Half width DStar mass
 302     if(TMath::Abs(mDSPDG - (d->InvMassDstarKpipi()))>fCutsRD[GetGlobalIndex(9,ptbin)]) return 0;
 303     // Half width Delta mass
 304
 305     if(TMath::Abs(deltaPDG-(d->DeltaInvMass())) > fCutsRD[GetGlobalIndex(10,ptbin)]) return 0;
 306
 307     // cut on soft pion pt
 308     if(b->Pt() < fCutsRD[GetGlobalIndex(11,ptbin)] || b->Pt() > fCutsRD[GetGlobalIndex(12,ptbin)]) return 0;
 309     // cut on the angle between D0 decay plane and soft pion
 310     if(d->AngleD0dkpPisoft() > fCutsRD[GetGlobalIndex(13,ptbin)]) return 0;
 311
 312   }
 313
 314   return returnvalue;
 315 }
 316 //_________________________________________________________________________________________________
 317 Int_t AliRDHFCutsDStartoKpipi::IsD0FromDStarSelected(Double_t pt, TObject* obj,Int_t selectionLevel) const {
 318   //
 319   // Apply selection for D0 from D*. The selection in on D0 prongs
 320   //
 321
 322   if(!fCutsRD){
 323     cout<<"Cut matrice not inizialized. Exit..."<<endl;
 324     return 0;
 325   }
 326
 327   AliAODRecoDecayHF2Prong* dd = (AliAODRecoDecayHF2Prong*)obj;
 328
 329   if(!dd){
 330     cout<<"AliAODRecoDecayHF2Prong null"<<endl;
 331     return 0;
 332   }
 333
 334   // selection on daughter tracks is done in IsSelected()
 335
 336   Int_t returnvalue=1;
 337
 338   // selection on candidate
 339   if(selectionLevel==AliRDHFCuts::kAll ||
 340      selectionLevel==AliRDHFCuts::kCandidate) {
 341
 342     // D0 mass
 343     Double_t mD0PDG = TDatabasePDG::Instance()->GetParticle(421)->Mass();
 344     // delta mass PDG
 345
 346     Int_t ptbin=PtBin(pt);
 347
 348     Double_t mD0,mD0bar,ctsD0,ctsD0bar;
 349
 350     Int_t okD0     =0;
 351     Int_t okD0bar  =0;
 352     okD0=1; okD0bar=1;
 353
 354     if(dd->PtProng(1) < fCutsRD[GetGlobalIndex(3,ptbin)] || dd->PtProng(0) < fCutsRD[GetGlobalIndex(4,ptbin)]) okD0 = 0;
 355     if(dd->PtProng(0) < fCutsRD[GetGlobalIndex(3,ptbin)] || dd->PtProng(1) < fCutsRD[GetGlobalIndex(4,ptbin)]) okD0bar = 0;
 356
 357     if(!okD0 && !okD0bar) return 0;
 358
 359     if(TMath::Abs(dd->Getd0Prong(1)) > fCutsRD[GetGlobalIndex(5,ptbin)] ||
 360        TMath::Abs(dd->Getd0Prong(0)) > fCutsRD[GetGlobalIndex(6,ptbin)]) okD0 = 0;
 361     if(TMath::Abs(dd->Getd0Prong(0)) > fCutsRD[GetGlobalIndex(6,ptbin)] ||
 362        TMath::Abs(dd->Getd0Prong(1)) > fCutsRD[GetGlobalIndex(5,ptbin)]) okD0bar = 0;
 363     if(!okD0 && !okD0bar) return 0;
 364
 365     if(dd->GetDCA() > fCutsRD[GetGlobalIndex(1,ptbin)]) return 0;
 366
 367     dd->InvMassD0(mD0,mD0bar);
 368     if(TMath::Abs(mD0-mD0PDG) > fCutsRD[GetGlobalIndex(0,ptbin)]) okD0 = 0;
 369     if(TMath::Abs(mD0bar-mD0PDG) > fCutsRD[GetGlobalIndex(0,ptbin)]) okD0bar = 0;
 370     if(!okD0 && !okD0bar) return 0;
 371
 372     dd->CosThetaStarD0(ctsD0,ctsD0bar);
 373     if(TMath::Abs(ctsD0) > fCutsRD[GetGlobalIndex(2,ptbin)]) okD0 = 0;
 374     if(TMath::Abs(ctsD0bar) > fCutsRD[GetGlobalIndex(2,ptbin)]) okD0bar = 0;
 375     if(!okD0 && !okD0bar) return 0;
 376
 377     if(dd->Prodd0d0() > fCutsRD[GetGlobalIndex(7,ptbin)]) return 0;
 378
 379     if(dd->CosPointingAngle() < fCutsRD[GetGlobalIndex(8,ptbin)]) return 0;
 380
 381     if (okD0) returnvalue=1; //cuts passed as D0
 382     if (okD0bar) returnvalue=2; //cuts passed as D0bar
 383     if (okD0 && okD0bar) returnvalue=3; //both
 384   }
 385
 386   return returnvalue;
 387 }
 388 //----------------------------------------------------------------------------------
 389 Bool_t AliRDHFCutsDStartoKpipi::IsInFiducialAcceptance(Double_t pt, Double_t y) const
 390 {
 391   //
 392   // D* fiducial acceptance region
 393   //
 394
 395   if(pt > 5.) {
 396     // applying cut for pt > 5 GeV
 397     AliDebug(4,Form("pt of D* = %f (> 5), cutting at |y| < 0.8\n",pt));
 398     if (TMath::Abs(y) > 0.8){
 399       return kFALSE;
 400     }
 401   } else {
 402
 403     // appliying smooth cut for pt < 5 GeV
 404     Double_t maxFiducialY = -0.2/15*pt*pt+1.9/15*pt+0.5;
 405     Double_t minFiducialY = 0.2/15*pt*pt-1.9/15*pt-0.5;
 406     AliDebug(2,Form("pt of D* = %f (< 5), cutting  according to the fiducial zone [%f, %f]\n",pt,minFiducialY,maxFiducialY));
 407     if (y < minFiducialY || y > maxFiducialY){
 408       return kFALSE;
 409     }
 410   }
 411
 412
 413   return kTRUE;
 414 }
 415 //_______________________________________________________________________________-
 416 Int_t AliRDHFCutsDStartoKpipi::IsSelectedPID(AliAODRecoDecayHF* dd)
 417 {
 418   //
 419   //  here the PID
 420
 421   // Double_t cw = dd->Pt();
 422
 423   if(dd) return 1;
 424   return 1;
 425 }