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