PWG3/vertexingHF/AliCFVertexingHF2Prong.cxx

   1 /**************************************************************************
   2  * Copyright(c) 1998-2009, 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 // Class for HF corrections as a function of many variables and step
  18 // Author : C. Zampolli, CERN
  19 // D. Caffarri, Univ & INFN Padova caffarri@pd.infn.it
  20 // Base class for HF Unfolding - agrelli@uu.nl
  21 //-----------------------------------------------------------------------
  22
  23 #include "AliAODRecoDecayHF2Prong.h"
  24 #include "AliAODMCParticle.h"
  25 #include "AliAODEvent.h"
  26 #include "TClonesArray.h"
  27 #include "AliCFVertexingHF.h"
  28 #include "AliESDtrack.h"
  29 #include "TDatabasePDG.h"
  30
  31 #include "AliCFVertexingHF2Prong.h"
  32 #include "AliCFContainer.h"
  33
  34 ClassImp(AliCFVertexingHF2Prong)
  35
  36
  37 //_________________________________________
  38   AliCFVertexingHF2Prong::AliCFVertexingHF2Prong(TClonesArray *mcArray, UShort_t originDselection):
  39           AliCFVertexingHF(mcArray, originDselection)
  40 {
  41         //
  42         // constructor
  43         //
  44
  45         SetNProngs(2);
  46 }
  47
  48
  49 //_____________________________________
  50 AliCFVertexingHF2Prong& AliCFVertexingHF2Prong::operator=(const AliCFVertexingHF2Prong& c)
  51 {
  52         //
  53         // copy constructor
  54         //
  55
  56         if  (this != &c) {
  57                 AliCFVertexingHF::operator=(c);
  58         }
  59         return *this;
  60 }
  61
  62 //__________________________________________
  63 Bool_t AliCFVertexingHF2Prong::SetRecoCandidateParam(AliAODRecoDecayHF *recoCand)
  64 {
  65         //
  66         // setting the recontructed candidate
  67         //
  68
  69         Bool_t bSignAssoc = kFALSE;
  70         fRecoCandidate = recoCand;
  71         if (!fRecoCandidate) {
  72                 AliError("fRecoCandidate not found, problem in assignement\n");
  73                 return bSignAssoc;
  74         }
  75
  76         if (fRecoCandidate->GetPrimaryVtx()) AliDebug(3,"fReco Candidate has a pointer to PrimVtx\n");
  77         if (recoCand->GetPrimaryVtx()) AliDebug(3,"Reco Cand has a pointer to PrimVtx\n");
  78
  79         Int_t pdgCand = 421;
  80         Int_t pdgDgD0toKpi[2]={321,211};
  81         Int_t nentries = fmcArray->GetEntriesFast();
  82
  83         AliDebug(3,Form("nentries = %d\n", nentries));
  84
  85         Int_t mcLabel = fRecoCandidate->MatchToMC(pdgCand,fmcArray,2,pdgDgD0toKpi);
  86         if (mcLabel == -1) return bSignAssoc;
  87         SetMCLabel(mcLabel);
  88         fmcPartCandidate = dynamic_cast<AliAODMCParticle*>(fmcArray->At(fmcLabel));
  89         if (!fmcPartCandidate){
  90                 AliDebug(3,"No part candidate");
  91                 return bSignAssoc;
  92         }
  93         bSignAssoc = kTRUE;
  94         return bSignAssoc;
  95 }
  96
  97 //______________________________________________
  98 Bool_t AliCFVertexingHF2Prong::GetGeneratedValuesFromMCParticle(Double_t* vectorMC)
  99 {
 100         //
 101         // collecting all the necessary info (pt, y, cosThetaStar, ptPi, ptKa, cT) from MC particle
 102         //
 103
 104         Bool_t bGenValues = kFALSE;
 105         Double_t vtx1[3] = {0,0,0};   // primary vertex
 106         Double_t vtx2daughter0[3] = {0,0,0};   // secondary vertex from daughter 0
 107         Double_t vtx2daughter1[3] = {0,0,0};   // secondary vertex from daughter 1
 108         fmcPartCandidate->XvYvZv(vtx1);  // cm
 109
 110         Int_t daughter0 = fmcPartCandidate->GetDaughter(0);
 111         Int_t daughter1 = fmcPartCandidate->GetDaughter(1);
 112         AliAODMCParticle* mcPartDaughter0 = dynamic_cast<AliAODMCParticle*>(fmcArray->At(daughter0));
 113         AliAODMCParticle* mcPartDaughter1 = dynamic_cast<AliAODMCParticle*>(fmcArray->At(daughter1));
 114
 115         // getting vertex from daughters
 116         mcPartDaughter0->XvYvZv(vtx2daughter0);  // cm
 117         mcPartDaughter1->XvYvZv(vtx2daughter1);  //cm
 118         if (TMath::Abs(vtx2daughter0[0] - vtx2daughter1[0]) > 1E-5 || TMath::Abs(vtx2daughter0[1] - vtx2daughter1[1]) > 1E-5 || TMath::Abs(vtx2daughter0[2] - vtx2daughter1[2])>1E-5) {
 119                 AliError("Daughters have different secondary vertex, skipping the track");
 120                 return bGenValues;
 121         }
 122
 123         Int_t nprongs = 2;
 124         Short_t charge = 0;
 125         // always instantiate the AliAODRecoDecay with the positive daughter first, the negative second
 126         AliAODMCParticle* positiveDaugh = mcPartDaughter0;
 127         AliAODMCParticle* negativeDaugh = mcPartDaughter1;
 128         if (mcPartDaughter0->GetPdgCode()<0 && mcPartDaughter1->GetPdgCode()>0){
 129                 // inverting in case the positive daughter is the second one
 130                 positiveDaugh = mcPartDaughter1;
 131                 negativeDaugh = mcPartDaughter0;
 132         }
 133         // getting the momentum from the daughters
 134         Double_t px[2] = {positiveDaugh->Px(), negativeDaugh->Px()};
 135         Double_t py[2] = {positiveDaugh->Py(), negativeDaugh->Py()};
 136         Double_t pz[2] = {positiveDaugh->Pz(), negativeDaugh->Pz()};
 137
 138         Double_t d0[2] = {0.,0.};
 139
 140         AliAODRecoDecayHF* decay = new AliAODRecoDecayHF(vtx1,vtx2daughter0,nprongs,charge,px,py,pz,d0);
 141
 142         Double_t cosThetaStar = 0.;
 143         Double_t cosThetaStarD0 = 0.;
 144         Double_t cosThetaStarD0bar = 0.;
 145         cosThetaStarD0 = decay->CosThetaStar(1,421,211,321);
 146         cosThetaStarD0bar = decay->CosThetaStar(0,421,321,211);
 147         if (fmcPartCandidate->GetPdgCode() == 421){  // D0
 148                 AliDebug(3, Form("D0, with pdgprong0 = %d, pdgprong1 = %d",mcPartDaughter0->GetPdgCode(),mcPartDaughter1->GetPdgCode()));
 149                 cosThetaStar = cosThetaStarD0;
 150         }
 151         else if (fmcPartCandidate->GetPdgCode() == -421){  // D0bar{
 152                 AliDebug(3, Form("D0bar, with pdgprong0 = %d, pdgprong1 = %d",mcPartDaughter0->GetPdgCode(),mcPartDaughter1->GetPdgCode()));
 153                 cosThetaStar = cosThetaStarD0bar;
 154         }
 155         else{
 156                 AliWarning("There are problems!! particle was expected to be either a D0 or a D0bar, check...");
 157                 return vectorMC;
 158         }
 159         if (cosThetaStar < -1 || cosThetaStar > 1) {
 160                 AliWarning("Invalid value for cosine Theta star, returning");
 161                 return bGenValues;
 162         }
 163
 164         //ct
 165         Double_t cT = decay->Ct(421);
 166         // get the pT of the daughters
 167         Double_t pTpi = 0.;
 168         Double_t pTK = 0.;
 169
 170         if (TMath::Abs(mcPartDaughter0->GetPdgCode()) == 211) {
 171                 pTpi = mcPartDaughter0->Pt();
 172                 pTK = mcPartDaughter1->Pt();
 173         }
 174         else {
 175                 pTpi = mcPartDaughter1->Pt();
 176                 pTK = mcPartDaughter0->Pt();
 177         }
 178
 179         vectorMC[0] = fmcPartCandidate->Pt();
 180         vectorMC[1] = fmcPartCandidate->Y() ;
 181         vectorMC[2] = cosThetaStar ;
 182         vectorMC[3] = pTpi ;
 183         vectorMC[4] = pTK ;
 184         vectorMC[5] = cT*1.E4 ;  // in micron
 185         vectorMC[6] = 0.;   // dummy value, meaningless in MC
 186         vectorMC[7] = 0.;   // dummy value, meaningless in MC, in micron
 187         vectorMC[8] = 0.;   // dummy value, meaningless in MC, in micron
 188         vectorMC[9] = -100000.; // dummy value, meaningless in MC, in micron^2
 189         vectorMC[10] = 1.01;    // dummy value, meaningless in MC
 190         vectorMC[11] = fmcPartCandidate->Phi();
 191         vectorMC[12] = fzMCVertex;    // z of reconstructed of primary vertex
 192         delete decay;
 193         bGenValues = kTRUE;
 194         return bGenValues;
 195 }
 196 //____________________________________________
 197 Bool_t AliCFVertexingHF2Prong::GetRecoValuesFromCandidate(Double_t *vectorReco) const
 198 {
 199         //
 200         // Getting the reconstructed values from the candidate
 201         //
 202
 203         Bool_t bFillRecoValues=kFALSE;
 204
 205         AliAODRecoDecayHF2Prong *d0toKpi = (AliAODRecoDecayHF2Prong*)fRecoCandidate;
 206
 207         if (d0toKpi->GetPrimaryVtx())AliDebug(3,"d0toKpi has primary vtx\n");
 208         if (fRecoCandidate->GetPrimaryVtx())AliDebug(3,"fRecoCandidate has primary vtx\n");
 209
 210         Double_t pt = d0toKpi->Pt();
 211         Double_t rapidity = d0toKpi->YD0();
 212         Double_t invMass=0.;
 213         Double_t cosThetaStar = 9999.;
 214         Double_t pTpi = 0.;
 215         Double_t pTK = 0.;
 216         Double_t dca = d0toKpi->GetDCA();
 217         Double_t d0pi = 0.;
 218         Double_t d0K = 0.;
 219         Double_t d0xd0 = d0toKpi->Prodd0d0();
 220         Double_t cosPointingAngle = d0toKpi->CosPointingAngle();
 221         Double_t phi = d0toKpi->Phi();
 222         Int_t pdgCode = fmcPartCandidate->GetPdgCode();
 223
 224         if (pdgCode > 0){
 225                 cosThetaStar = d0toKpi->CosThetaStarD0();
 226                 pTpi = d0toKpi->PtProng(0);
 227                 pTK = d0toKpi->PtProng(1);
 228                 d0pi = d0toKpi->Getd0Prong(0);
 229                 d0K = d0toKpi->Getd0Prong(1);
 230                 invMass=d0toKpi->InvMassD0();
 231         }
 232         else {
 233                 cosThetaStar = d0toKpi->CosThetaStarD0bar();
 234                 pTpi = d0toKpi->PtProng(1);
 235                 pTK = d0toKpi->PtProng(0);
 236                 d0pi = d0toKpi->Getd0Prong(1);
 237                 d0K = d0toKpi->Getd0Prong(0);
 238                 invMass= d0toKpi->InvMassD0bar();
 239         }
 240
 241         Double_t cT = d0toKpi->CtD0();
 242
 243         vectorReco[0] = pt;
 244         vectorReco[1] = rapidity;
 245         vectorReco[2] = cosThetaStar;
 246         vectorReco[3] = pTpi;
 247         vectorReco[4] = pTK;
 248         vectorReco[5] = cT*1.E4;  // in micron
 249         vectorReco[6] = dca*1.E4;  // in micron
 250         vectorReco[7] = d0pi*1.E4;  // in micron
 251         vectorReco[8] = d0K*1.E4;  // in micron
 252         vectorReco[9] = d0xd0*1.E8;  // in micron^2
 253         vectorReco[10] = cosPointingAngle;  // in micron
 254         vectorReco[11] = phi;
 255         vectorReco[12] = fzPrimVertex;    // z of reconstructed of primary vertex
 256         bFillRecoValues = kTRUE;
 257
 258         return bFillRecoValues;
 259 }
 260
 261 //_____________________________________________________________
 262 Bool_t AliCFVertexingHF2Prong::CheckMCChannelDecay() const
 263 {
 264         //
 265         // checking the MC decay channel
 266         //
 267         Bool_t checkCD = kFALSE;
 268         Int_t daughter0 = fmcPartCandidate->GetDaughter(0);
 269         Int_t daughter1 = fmcPartCandidate->GetDaughter(1);
 270         AliAODMCParticle* mcPartDaughter0 = dynamic_cast<AliAODMCParticle*>(fmcArray->At(daughter0));
 271         AliAODMCParticle* mcPartDaughter1 = dynamic_cast<AliAODMCParticle*>(fmcArray->At(daughter1));
 272
 273         if (!mcPartDaughter0 || !mcPartDaughter1) {
 274                 AliDebug (2,"Problems in the MC Daughters\n");
 275                 return checkCD;
 276         }
 277
 278         if (!(TMath::Abs(mcPartDaughter0->GetPdgCode())==321 &&
 279               TMath::Abs(mcPartDaughter1->GetPdgCode())==211) &&
 280             !(TMath::Abs(mcPartDaughter0->GetPdgCode())==211 &&
 281               TMath::Abs(mcPartDaughter1->GetPdgCode())==321)) {
 282                 AliDebug(2, "The D0 MC doesn't come from a Kpi decay, skipping!!");
 283                 return checkCD;
 284         }
 285
 286         checkCD = kTRUE;
 287         return checkCD;
 288
 289 }
 290