PWG3/vertexingHF/AliAODRecoCascadeHF.cxx

   1 /**************************************************************************
   2  * Copyright(c) 1998-2008, 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 AOD reconstructed heavy-flavour cascades
  19 //
  20 // Author: X-M. Zhang, zhangxm@ccnu.iop.edu.cn
  21 /////////////////////////////////////////////////////////////
  22
  23 #include <TVector3.h>
  24 #include <TDatabasePDG.h>
  25 #include "AliAODRecoDecay.h"
  26 #include "AliAODVertex.h"
  27 #include "AliAODRecoDecayHF2Prong.h"
  28 #include "AliAODRecoCascadeHF.h"
  29
  30 ClassImp(AliAODRecoCascadeHF)
  31 //-----------------------------------------------------------------------------
  32
  33 AliAODRecoCascadeHF::AliAODRecoCascadeHF() :
  34   AliAODRecoDecayHF2Prong(),
  35   f2Prong()
  36 {
  37   //
  38   // Default Constructor
  39   //
  40 }
  41 //-----------------------------------------------------------------------------
  42 AliAODRecoCascadeHF::AliAODRecoCascadeHF(AliAODVertex *vtx2, Short_t charge,
  43                                          Double_t *px, Double_t *py, Double_t *pz,
  44                                          Double_t *d0, Double_t *d0err, Double_t dca) :
  45   AliAODRecoDecayHF2Prong(vtx2, px, py, pz, d0, d0err, dca),
  46   f2Prong()
  47 {
  48   //
  49   //  Constructor with AliAODVertex for decay vertex
  50   //
  51   SetCharge(charge);
  52 }
  53 //-----------------------------------------------------------------------------
  54 AliAODRecoCascadeHF::AliAODRecoCascadeHF(AliAODVertex *vtx2, Short_t charge,
  55                                          Double_t *d0, Double_t *d0err, Double_t dca) :
  56   AliAODRecoDecayHF2Prong(vtx2, d0, d0err, dca),
  57   f2Prong()
  58 {
  59   //
  60   //  Constructor with decay vertex and without prongs momenta
  61   //
  62   SetCharge(charge);
  63 }
  64 //-----------------------------------------------------------------------------
  65 AliAODRecoCascadeHF::AliAODRecoCascadeHF(const AliAODRecoCascadeHF &source) :
  66   AliAODRecoDecayHF2Prong(source),
  67   f2Prong()
  68 {
  69   //
  70   // Copy constructor
  71   //
  72 }
  73 //-----------------------------------------------------------------------------
  74 AliAODRecoCascadeHF &AliAODRecoCascadeHF::operator=(const AliAODRecoCascadeHF &source)
  75 {
  76   //
  77   // assignment operator
  78   //
  79   if(&source == this) return *this;
  80
  81   AliAODRecoDecayHF2Prong::operator=(source);
  82
  83   f2Prong = source.f2Prong;
  84
  85   return *this;
  86 }
  87 //-----------------------------------------------------------------------------
  88 AliAODRecoCascadeHF::~AliAODRecoCascadeHF()
  89 {
  90   //
  91   // Default Destructor
  92   //
  93 }
  94 //-----------------------------------------------------------------------------
  95 Double_t AliAODRecoCascadeHF::InvMassDstarKpipi() const
  96 {
  97   //
  98   // 3 prong invariant mass of the D0 daughters and the soft pion
  99   //
 100
 101   Double_t px[3],py[3],pz[3];
 102   UInt_t pdg[3]={321,211,211};
 103   pdg[0] = (Charge()>0 ? 211 : 321); // positive daughter of D0
 104   px[0] = Get2Prong()->PxProng(0);
 105   py[0] = Get2Prong()->PyProng(0);
 106   pz[0] = Get2Prong()->PzProng(0);
 107   pdg[1] = (Charge()>0 ? 321 : 211); // negative daughter of D0
 108   px[1] = Get2Prong()->PxProng(1);
 109   py[1] = Get2Prong()->PyProng(1);
 110   pz[1] = Get2Prong()->PzProng(1);
 111   pdg[2] = 211; // soft pion
 112   px[2] = PxProng(0);
 113   py[2] = PyProng(0);
 114   pz[2] = PzProng(0);
 115   Short_t dummycharge=0;
 116   Double_t dummyd0[3]={0,0,0};
 117   AliAODRecoDecay *rd = new AliAODRecoDecay(0x0,3,dummycharge,px,py,pz,dummyd0);
 118
 119   Double_t minv = rd->InvMass(3,pdg);
 120
 121   delete rd; rd=NULL;
 122
 123   return minv;
 124 }
 125 //-----------------------------------------------------------------------------
 126 Bool_t AliAODRecoCascadeHF::SelectDstar(const Double_t *cutsDstar,
 127                                         const Double_t *cutsD0,
 128                                         Bool_t testD0) const
 129 {
 130   //
 131   // cutsDstar[0] = inv. mass half width of D* [GeV]
 132   // cutsDstar[1] = half width of (M_Kpipi-M_D0) [GeV]
 133   // cutsDstar[2] = PtMin of pi_s [GeV/c]
 134   // cutsDstar[3] = PtMax of pi_s [GeV/c]
 135   // cutsDstar[4] = theta, angle between the pi_s and decay plane of the D0 [rad]
 136   //
 137   // cutsD0[0] = inv. mass half width [GeV]
 138   // cutsD0[1] = dca [cm]
 139   // cutsD0[2] = cosThetaStar
 140   // cutsD0[3] = pTK [GeV/c]
 141   // cutsD0[4] = pTPi [GeV/c]
 142   // cutsD0[5] = d0K [cm]   upper limit!
 143   // cutsD0[6] = d0Pi [cm]  upper limit!
 144   // cutsD0[7] = d0d0 [cm^2]
 145   // cutsD0[8] = cosThetaPoint
 146
 147
 148   // check that the D0 passes the cuts
 149   // (if we have a D*+, it has to pass as D0,
 150   //  if we have a D*-, it has to pass as D0bar)
 151
 152   if(testD0) {
 153     Int_t okD0=0,okD0bar=0;
 154     Get2Prong()->SelectD0(cutsD0,okD0,okD0bar);
 155     if((Charge()==+1 && !okD0) || (Charge()==-1 && !okD0bar)) return kFALSE;
 156   }
 157
 158   if( (PtProng(0)<cutsDstar[2]) || (PtProng(0)>cutsDstar[3]) ) return kFALSE;
 159
 160   Double_t mDstar = TDatabasePDG::Instance()->GetParticle(413)->Mass();
 161   Double_t invmDstar = InvMassDstarKpipi();
 162   if(TMath::Abs(mDstar-invmDstar)>cutsDstar[0]) return kFALSE;
 163
 164   Double_t mD0 = TDatabasePDG::Instance()->GetParticle(421)->Mass();
 165   if(TMath::Abs((mDstar-mD0)-DeltaInvMass())>cutsDstar[1]) return kFALSE;
 166
 167   TVector3 p3Trk0(Get2Prong()->PxProng(0),Get2Prong()->PyProng(0),Get2Prong()->PzProng(0)); // from D0
 168   TVector3 p3Trk1(Get2Prong()->PxProng(1),Get2Prong()->PyProng(1),Get2Prong()->PzProng(1)); // from D0
 169   TVector3 p3Trk2(PxProng(0),PyProng(0),PzProng(0)); // pi_s
 170
 171   TVector3 perp = p3Trk0.Cross(p3Trk1);
 172   Double_t theta = p3Trk2.Angle(perp);
 173   if(theta>(TMath::Pi()-theta)) theta = TMath::Pi() - theta;
 174   theta = TMath::Pi()/2. - theta;
 175
 176   if(theta>cutsDstar[4]) return kFALSE;
 177
 178   Double_t alpha = p3Trk0.Angle(p3Trk2);
 179   Double_t belta = p3Trk1.Angle(p3Trk2);
 180
 181   Double_t cosphi01 = TMath::Cos(alpha) / TMath::Cos(theta);
 182   Double_t cosphi02 = TMath::Cos(belta) / TMath::Cos(theta);
 183
 184   Double_t phi01 = TMath::ACos(cosphi01);
 185   Double_t phi02 = TMath::ACos(cosphi02);
 186   Double_t phi00 = p3Trk0.Angle(p3Trk1);
 187
 188   if((phi01>phi00) || (phi02>phi00)) return kFALSE;
 189
 190   return kTRUE;
 191 }
 192 //-----------------------------------------------------------------------------