PWGHF/vertexingHF/AliCFVertexingHFLctoV0bachelor.cxx

   1 /**************************************************************************
   2  * Copyright(c) 1998-2011, 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 // Class for HF corrections as a function of many variables and steps
  20 // For Lc->V0+bachelor
  21 //
  22 // - Based on AliCFVertexingHFCascade -
  23 //
  24 // Contact : A.De Caro - decaro@sa.infn.it
  25 //           Centro 'E.Fermi' - Rome (Italy)
  26 //           INFN and University of Salerno (Italy)
  27 //
  28 //-----------------------------------------------------------------------
  29
  30 #include "TDatabasePDG.h"
  31 #include "TClonesArray.h"
  32 #include "AliAODv0.h"
  33 #include "AliAODMCParticle.h"
  34 #include "AliAODRecoDecayHF.h"
  35 #include "AliAODRecoCascadeHF.h"
  36 #include "AliCFTaskVertexingHF.h"
  37 #include "AliCFContainer.h"
  38 #include "AliCFVertexingHF.h"
  39 #include "AliCFVertexingHFLctoV0bachelor.h"
  40
  41 #include <Riostream.h>
  42
  43 using std::cout;
  44 using std::endl;
  45
  46 ClassImp(AliCFVertexingHFLctoV0bachelor)
  47
  48 //_________________________________________
  49 AliCFVertexingHFLctoV0bachelor::AliCFVertexingHFLctoV0bachelor():
  50 fGenLcOption(0)
  51 {
  52   // standard constructor
  53 }
  54
  55 //_____________________________________
  56 AliCFVertexingHFLctoV0bachelor::AliCFVertexingHFLctoV0bachelor(TClonesArray *mcArray, UShort_t originDselection, Int_t lcDecay):
  57 AliCFVertexingHF(mcArray, originDselection),
  58   fGenLcOption(lcDecay)
  59 {
  60   // standard constructor
  61
  62   SetNProngs(3);
  63   fPtAccCut=new Float_t[fProngs];
  64   fEtaAccCut=new Float_t[fProngs];
  65   for(Int_t iP=0; iP<fProngs; iP++){
  66     fPtAccCut[iP]=0.1;
  67     fEtaAccCut[iP]=0.9;
  68   }
  69
  70 }
  71
  72
  73 //_____________________________________
  74 AliCFVertexingHFLctoV0bachelor& AliCFVertexingHFLctoV0bachelor::operator=(const AliCFVertexingHFLctoV0bachelor& c)
  75 {
  76   // operator =
  77
  78   if  (this != &c) {
  79     AliCFVertexingHF::operator=(c);
  80   }
  81
  82   return *this;
  83
  84 }
  85
  86 //__________________________________________
  87 Bool_t AliCFVertexingHFLctoV0bachelor::SetRecoCandidateParam(AliAODRecoDecayHF *recoCand)
  88 {
  89   // set the AliAODRecoDecay candidate
  90
  91   Bool_t bSignAssoc = kFALSE;
  92
  93   fRecoCandidate = recoCand;
  94   if (!fRecoCandidate) {
  95     AliError("fRecoCandidate not found, problem in assignement\n");
  96     return bSignAssoc;
  97   }
  98
  99   if (fRecoCandidate->GetPrimaryVtx()) AliDebug(4,"fReco Candidate has a pointer to PrimVtx\n");
 100
 101   AliAODRecoCascadeHF* lcV0bachelor = dynamic_cast<AliAODRecoCascadeHF*>(fRecoCandidate);
 102   if (!lcV0bachelor) {
 103     AliDebug(2,"cascade doesn't exist, skipped!");
 104     return bSignAssoc;
 105   }
 106
 107   if (!lcV0bachelor->GetSecondaryVtx()) {
 108     AliDebug(2,"No secondary vertex for cascade");
 109     return bSignAssoc;
 110   }
 111
 112   if (lcV0bachelor->GetNDaughters()!=2) {
 113     AliDebug(2,Form("No 2 daughters for current cascade (nDaughters=%d)",lcV0bachelor->GetNDaughters()));
 114     return bSignAssoc;
 115   }
 116
 117   AliVTrack *cascTrk0 = dynamic_cast<AliVTrack*>(lcV0bachelor->GetDaughter(0));
 118   AliVTrack *cascTrk1 = dynamic_cast<AliVTrack*>(lcV0bachelor->GetDaughter(1));
 119   if (!cascTrk0 || !cascTrk1) {
 120     AliDebug(2,"At least one of V0daughters doesn't exist");
 121     return bSignAssoc;
 122   }
 123
 124   AliAODv0 * v0part = dynamic_cast<AliAODv0*>(lcV0bachelor->Getv0());
 125   AliAODTrack * bachPart = dynamic_cast<AliAODTrack*>(lcV0bachelor->GetBachelor());
 126   if (!v0part || !bachPart) {
 127     AliDebug(2,"No V0 or no bachelor for current cascade");
 128     return bSignAssoc;
 129   }
 130
 131   if (bachPart->GetID()<0) {
 132     AliDebug(2,Form("Bachelor has negative ID %d",bachPart->GetID()));
 133     return bSignAssoc;
 134   }
 135
 136   if (!v0part->GetSecondaryVtx()) {
 137     AliDebug(2,"No secondary vertex for V0 by cascade");
 138     return bSignAssoc;
 139   }
 140
 141   if (v0part->GetNDaughters()!=2) {
 142     AliDebug(2,Form("No 2 daughters for V0 of current cascade (onTheFly=%d, nDaughters=%d)",v0part->GetOnFlyStatus(),v0part->GetNDaughters()));
 143     return bSignAssoc;
 144   }
 145
 146   AliVTrack *trk0 = dynamic_cast<AliVTrack*>(v0part->GetDaughter(0));
 147   AliVTrack *trk1 = dynamic_cast<AliVTrack*>(v0part->GetDaughter(1));
 148   if (!trk0 || !trk1) {
 149     AliDebug(2,"At least one of V0daughters doesn't exist");
 150     return bSignAssoc;
 151   }
 152
 153   //if (trk0->GetLabel()<0 || trk1->GetLabel()<0) {
 154   //AliDebug(2,Form("At least one of V0daughters has label negative (%d %d)",trk0->GetLabel(),trk1->GetLabel()));
 155   //return bSignAssoc;
 156   //}
 157
 158  if (trk0->GetID()<0 || trk1->GetID()<0) {
 159     AliDebug(2,Form("At least one of V0 daughters has negative ID %d %d",trk0->GetID(),trk1->GetID()));
 160     return bSignAssoc;
 161   }
 162
 163   Int_t pdgCand = 4122;
 164   Int_t mcLabel = -1;
 165   Int_t mcLabelK0S = -1;
 166   Int_t mcLabelLambda = -1;
 167
 168   // Lc->K0S+p and cc
 169   Int_t pdgDgLctoV0bachelor[2]={2212,310}; // first bachelor, second V0
 170   Int_t pdgDgV0toDaughters[2]={211,211};
 171   mcLabelK0S = lcV0bachelor->MatchToMC(pdgCand,pdgDgLctoV0bachelor[1],pdgDgLctoV0bachelor,pdgDgV0toDaughters,fmcArray,kTRUE);
 172
 173   // Lc->Lambda+pi and cc
 174   pdgDgLctoV0bachelor[0]=211, pdgDgLctoV0bachelor[1]=3122; // first bachelor, second V0
 175   pdgDgV0toDaughters[0]=2212,  pdgDgV0toDaughters[1]=211;
 176   mcLabelLambda = lcV0bachelor->MatchToMC(pdgCand,pdgDgLctoV0bachelor[1],pdgDgLctoV0bachelor,pdgDgV0toDaughters,fmcArray,kTRUE);
 177
 178   if (mcLabelK0S!=-1 && mcLabelLambda!=-1)
 179     AliDebug(2,"Strange: current Lc->V0+bachelor candidate has two MC different labels!");
 180
 181   if (fGenLcOption==kCountK0Sp) {
 182     if (mcLabelK0S!=-1) mcLabel=mcLabelK0S;
 183     if (mcLabelLambda!=-1) {
 184       mcLabel=-1;
 185       fFake = 0;    // fake candidate
 186       if (fFakeSelection==1) return bSignAssoc;
 187     }
 188   }
 189   else if (fGenLcOption==kCountLambdapi) {
 190     if (mcLabelLambda!=-1) mcLabel=mcLabelLambda;
 191     if (mcLabelK0S!=-1) {
 192       mcLabel=-1;
 193       fFake = 0;    // fake candidate
 194       if (fFakeSelection==1) return bSignAssoc;
 195     }
 196   }
 197
 198   if (mcLabel==-1) {
 199     AliDebug(4,"No mcLabel found for current candidate");
 200     return bSignAssoc;
 201   }
 202   AliDebug(1,Form("Found mcLabel (%d) for current candidate (onTheFly=%1d)",mcLabel,v0part->GetOnFlyStatus()));
 203
 204   if (fRecoCandidate->NumberOfFakeDaughters()>0){
 205     fFake = 0;    // fake candidate
 206     if (fFakeSelection==1) return bSignAssoc;
 207   }
 208   if (fRecoCandidate->NumberOfFakeDaughters()==0){
 209     fFake = 2;    // non-fake candidate
 210     if (fFakeSelection==2) return bSignAssoc;
 211   }
 212
 213   SetMCLabel(mcLabel); // fmcLabel=mcLabel
 214   fmcPartCandidate = dynamic_cast<AliAODMCParticle*>(fmcArray->At(fmcLabel));
 215   if (!fmcPartCandidate){
 216     AliDebug(3,"No MC object for current candidate");
 217     return bSignAssoc;
 218   }
 219
 220   bSignAssoc = kTRUE;
 221   return bSignAssoc;
 222
 223 }
 224
 225 //______________________________________________
 226 Bool_t AliCFVertexingHFLctoV0bachelor::GetGeneratedValuesFromMCParticle(Double_t* vectorMC)
 227 {
 228   //
 229   // collecting all the necessary info (pt, y, invMassV0, cosPAwrtPVxV0, onTheFlyStatusV0) from MC particle
 230   // (additional infos: pTbachelor, pTV0pos, pTV0neg, phi, dcaV0, cTV0, cT, cosPA)
 231   //
 232
 233   Bool_t bGenValues = kFALSE;
 234
 235   if (fmcPartCandidate->GetNDaughters()!=2) {
 236     AliDebug(2,"Lc MC particle doesn't decay in 2 daughters");
 237     return bGenValues;
 238   }
 239
 240   Int_t daughter0lc = fmcPartCandidate->GetDaughter(0);
 241   Int_t daughter1lc = fmcPartCandidate->GetDaughter(1);
 242   if (daughter0lc<=0 || daughter1lc<=0) {
 243     AliDebug(2,"Lc daughters are not in MC array");
 244     return bGenValues;
 245   }
 246
 247   AliAODMCParticle* mcPartDaughter0 = dynamic_cast<AliAODMCParticle*>(fmcArray->At(daughter0lc));
 248   AliAODMCParticle* mcPartDaughter1 = dynamic_cast<AliAODMCParticle*>(fmcArray->At(daughter1lc));
 249   if (!mcPartDaughter0 || !mcPartDaughter1) {
 250     AliDebug(2,"Problems in the MC Daughters\n");
 251     return bGenValues;
 252   }
 253
 254   if ( fGenLcOption==kCountLambdapi &&
 255        !(TMath::Abs(mcPartDaughter0->GetPdgCode())==3122 &&
 256          TMath::Abs(mcPartDaughter1->GetPdgCode())==211) &&
 257        !(TMath::Abs(mcPartDaughter1->GetPdgCode())==3122 &&
 258          TMath::Abs(mcPartDaughter0->GetPdgCode())==211) ) return bGenValues;
 259   if ( fGenLcOption==kCountK0Sp &&
 260        !(TMath::Abs(mcPartDaughter0->GetPdgCode())==2212 &&
 261          TMath::Abs(mcPartDaughter1->GetPdgCode())==311) &&
 262        !(TMath::Abs(mcPartDaughter1->GetPdgCode())==2212 &&
 263          TMath::Abs(mcPartDaughter0->GetPdgCode())==311) ) return bGenValues;
 264
 265   if ( (TMath::Abs(mcPartDaughter0->GetPdgCode())==311   &&
 266         TMath::Abs(mcPartDaughter1->GetPdgCode())==2212) ||
 267        (TMath::Abs(mcPartDaughter0->GetPdgCode())==3122  &&
 268         TMath::Abs(mcPartDaughter1->GetPdgCode())==211) )
 269     bGenValues = FillVectorFromMCarray(mcPartDaughter1,mcPartDaughter0,vectorMC);
 270   else if ( (TMath::Abs(mcPartDaughter1->GetPdgCode())==311   &&
 271              TMath::Abs(mcPartDaughter0->GetPdgCode())==2212) ||
 272             (TMath::Abs(mcPartDaughter1->GetPdgCode())==3122  &&
 273              TMath::Abs(mcPartDaughter0->GetPdgCode())==211) )
 274     bGenValues = FillVectorFromMCarray(mcPartDaughter0,mcPartDaughter1,vectorMC);
 275
 276   if (!bGenValues)
 277     AliDebug(2,"There is something wrong in filling MC vector");
 278
 279   return bGenValues;
 280
 281 }
 282
 283 //____________________________________________
 284 Bool_t AliCFVertexingHFLctoV0bachelor::GetRecoValuesFromCandidate(Double_t *vectorReco) const
 285 {
 286   // read the variables for the container
 287
 288   Bool_t bFillRecoValues = kFALSE;
 289
 290   //Get the Lc and the V0 from Lc
 291   AliAODRecoCascadeHF* lcV0bachelor = dynamic_cast<AliAODRecoCascadeHF*>(fRecoCandidate);
 292
 293   if (!lcV0bachelor) {
 294     AliDebug(2,"Current cascade doesn't exist, skipped");
 295     return bFillRecoValues;
 296   }
 297
 298   if (!lcV0bachelor->GetSecondaryVtx()) {
 299     AliDebug(2,"No secondary vertex for cascade");
 300     return bFillRecoValues;
 301   }
 302
 303   if (lcV0bachelor->GetNDaughters()!=2) {
 304     AliDebug(2,Form("No 2 daughters for current cascade (nDaughters=%d)",lcV0bachelor->GetNDaughters()));
 305     return bFillRecoValues;
 306   }
 307
 308   AliAODTrack* bachelor = dynamic_cast<AliAODTrack*>(lcV0bachelor->GetBachelor());
 309   AliAODv0* v0toDaughters = dynamic_cast<AliAODv0*>(lcV0bachelor->Getv0());
 310   if (!bachelor || !v0toDaughters) {
 311     AliDebug(2,"No V0 or bachelor in this reco candidate, skipping!");
 312     return bFillRecoValues;
 313   }
 314
 315   if (!v0toDaughters->GetSecondaryVtx()) {
 316     AliDebug(2,"No secondary vertex for V0 by cascade");
 317     return bFillRecoValues;
 318   }
 319
 320   if (v0toDaughters->GetNDaughters()!=2) {
 321     AliDebug(2,Form("current V0 has not 2 daughters (onTheFly=%d, nDaughters=%d)",v0toDaughters->GetOnFlyStatus(),v0toDaughters->GetNDaughters()));
 322     return bFillRecoValues;
 323   }
 324
 325   Bool_t onTheFlyStatus = v0toDaughters->GetOnFlyStatus();
 326   AliAODTrack* v0positiveTrack;
 327   AliAODTrack* v0negativeTrack;
 328   if (onTheFlyStatus) {
 329     v0positiveTrack = dynamic_cast<AliAODTrack*>(lcV0bachelor->Getv0NegativeTrack());
 330     v0negativeTrack = dynamic_cast<AliAODTrack*>(lcV0bachelor->Getv0PositiveTrack());
 331   } else {
 332     v0positiveTrack = dynamic_cast<AliAODTrack*>(lcV0bachelor->Getv0PositiveTrack());
 333     v0negativeTrack = dynamic_cast<AliAODTrack*>(lcV0bachelor->Getv0NegativeTrack());
 334   }
 335
 336   if (!v0positiveTrack || !v0negativeTrack) {
 337     AliDebug(2,"No V0daughters in this reco candidate, skipping!");
 338     return bFillRecoValues;
 339   }
 340
 341   Double_t pt = lcV0bachelor->Pt();
 342   Double_t rapidity = lcV0bachelor->Y(4122);
 343
 344   Double_t cosPAwrtPrimVtxV0 = lcV0bachelor->CosV0PointingAngle();
 345
 346   //Double_t pTbachelor = bachelor->Pt();
 347   Double_t pbachelor = bachelor->P();
 348   //Double_t pTV0pos = v0positiveTrack->Pt();
 349   //Double_t pTV0neg = v0negativeTrack->Pt();
 350   Double_t phi = lcV0bachelor->Phi();
 351   Double_t dcaV0 = v0toDaughters->GetDCA();
 352   //Double_t cTLc = lcV0bachelor->Ct(4122); // wrt PrimVtx
 353   //Double_t dcaLc = lcV0bachelor->GetDCA();
 354   Double_t cosPointingAngleLc = lcV0bachelor->CosPointingAngle();
 355
 356   /*
 357   Double_t cTV0 = 0.;
 358   AliAODVertex *vtx0 = (AliAODVertex*)lcV0bachelor->GetPrimaryVtx();
 359   if (!vtx0) {
 360     AliDebug(2,"Candidate has not primary vtx");
 361   } else {
 362     Double_t primVtxPos[3] = {0.,0.,0.}; vtx0->GetXYZ(primVtxPos);
 363     if (fGenLcOption==kCountK0Sp) {
 364       cTV0 = v0toDaughters->Ct(310,primVtxPos);
 365     } else if (fGenLcOption==kCountLambdapi) {
 366       cTV0 = v0toDaughters->Ct(3122,primVtxPos);
 367     }
 368   }
 369   */
 370
 371   Double_t invMassV0 = 0.;
 372   if (fGenLcOption==kCountLambdapi) {
 373
 374     Short_t bachelorCharge = bachelor->Charge();
 375     if (bachelorCharge==1) {
 376       invMassV0 = v0toDaughters->MassLambda();
 377     } else if (bachelorCharge==-1) {
 378       invMassV0 = v0toDaughters->MassAntiLambda();
 379     }
 380
 381   } else if (fGenLcOption==kCountK0Sp) {
 382
 383     invMassV0 = v0toDaughters->MassK0Short();
 384
 385   }
 386
 387   vectorReco[0]  = pt;
 388   vectorReco[1]  = rapidity;
 389   vectorReco[2]  = phi;
 390   vectorReco[3]  = onTheFlyStatus;
 391   vectorReco[4]  = fzPrimVertex;
 392   vectorReco[5]  = fCentValue;
 393   vectorReco[6]  = fFake; // whether the reconstructed candidate was a fake (fFake = 0) or not (fFake = 2)
 394   vectorReco[7]  = fMultiplicity;
 395
 396   if (fConfiguration==AliCFTaskVertexingHF::kSnail) {
 397     //vectorReco[8]  = pTbachelor;
 398     vectorReco[8]  = pbachelor;
 399     vectorReco[9]  = v0toDaughters->Pt();
 400     if (fGenLcOption==kCountLambdapi) {
 401       vectorReco[10] = v0toDaughters->Y(3122);
 402     } else if (fGenLcOption==kCountK0Sp) {
 403       vectorReco[10] = v0toDaughters->Y(310);
 404     }
 405     vectorReco[11] = v0toDaughters->Phi();
 406     vectorReco[12] = invMassV0;
 407     vectorReco[13] = dcaV0;
 408     vectorReco[14] = cosPAwrtPrimVtxV0;
 409     vectorReco[15] = cosPointingAngleLc;
 410     //vectorReco[16] = cTV0*1.E4; // in micron
 411     //vectorReco[17] = cTLc*1.E4; // in micron
 412   }
 413
 414   bFillRecoValues = kTRUE;
 415
 416   return bFillRecoValues;
 417 }
 418
 419 //_____________________________________________________________
 420 Bool_t AliCFVertexingHFLctoV0bachelor::CheckMCChannelDecay() const
 421 {
 422   // check the required decay channel
 423
 424   Bool_t checkCD = kFALSE;
 425   Double_t sumPxDau=0.;
 426   Double_t sumPyDau=0.;
 427   Double_t sumPzDau=0.;
 428
 429
 430   if (fmcPartCandidate->GetNDaughters()!=2) {
 431     AliDebug(2, Form("The MC particle doesn't decay in 2 particles, skipping!!"));
 432     return checkCD;
 433   }
 434
 435   Int_t daughter0 = fmcPartCandidate->GetDaughter(0);
 436   Int_t daughter1 = fmcPartCandidate->GetDaughter(1);
 437   if (daughter0<=0 || daughter1<=0){
 438     AliDebug(2, Form("The MC particle doesn't have correct daughters, skipping!!"));
 439     return checkCD;
 440   }
 441   AliAODMCParticle* mcPartDaughter0 = dynamic_cast<AliAODMCParticle*>(fmcArray->At(daughter0));
 442   AliAODMCParticle* mcPartDaughter1 = dynamic_cast<AliAODMCParticle*>(fmcArray->At(daughter1));
 443   if (!mcPartDaughter0 || !mcPartDaughter1) {
 444     AliDebug(2,"Problems in the MC Daughters\n");
 445     return checkCD;
 446   }
 447
 448   // Lc -> Lambda + pion AND cc
 449   if (fGenLcOption==kCountLambdapi) {
 450
 451     if (!(TMath::Abs(mcPartDaughter0->GetPdgCode())==3122 &&
 452           TMath::Abs(mcPartDaughter1->GetPdgCode())==211) &&
 453         !(TMath::Abs(mcPartDaughter0->GetPdgCode())==211  &&
 454           TMath::Abs(mcPartDaughter1->GetPdgCode())==3122)) {
 455       AliDebug(2, "The Lc MC doesn't decay in Lambda+pion (or cc), skipping!!");
 456       return checkCD;
 457     }
 458
 459     if (TMath::Abs(mcPartDaughter0->GetPdgCode())==3122) {
 460       mcPartDaughter0 = dynamic_cast<AliAODMCParticle*>(fmcArray->At(daughter1)); // the bachelor
 461       mcPartDaughter1 = dynamic_cast<AliAODMCParticle*>(fmcArray->At(daughter0)); // the V0
 462     }
 463     if (!mcPartDaughter0 || !mcPartDaughter1) {
 464       AliDebug(2,"Problems in the MC Daughters\n");
 465       return checkCD;
 466     }
 467     // momentum of the bachelor
 468     sumPxDau+=mcPartDaughter0->Px();
 469     sumPyDau+=mcPartDaughter0->Py();
 470     sumPzDau+=mcPartDaughter0->Pz();
 471
 472     if (mcPartDaughter1->GetNDaughters()!=2) {
 473       AliDebug(2, "The Lambda MC particle doesn't decay in 2 particles, skipping!!");
 474       return checkCD;
 475     }
 476
 477     Int_t daughter1D0 = mcPartDaughter1->GetDaughter(0);
 478     Int_t daughter1D1 = mcPartDaughter1->GetDaughter(1);
 479     if (daughter1D0<=0 || daughter1D1<=0) {
 480       AliDebug(2, Form("The Lambda MC particle doesn't have correct daughters, skipping!!"));
 481       return checkCD;
 482     }
 483
 484     AliAODMCParticle* mcPartDaughter1D0 = dynamic_cast<AliAODMCParticle*>(fmcArray->At(daughter1D0));
 485     AliAODMCParticle* mcPartDaughter1D1 = dynamic_cast<AliAODMCParticle*>(fmcArray->At(daughter1D1));
 486     if(!mcPartDaughter1D0 || !mcPartDaughter1D1) {
 487       AliError("The Lambda daughter particle not found in MC array");
 488       return checkCD;
 489     }
 490
 491     if (!(TMath::Abs(mcPartDaughter1D0->GetPdgCode())==211   &&
 492           TMath::Abs(mcPartDaughter1D1->GetPdgCode())==2212) &&
 493         !(TMath::Abs(mcPartDaughter1D0->GetPdgCode())==2212  &&
 494           TMath::Abs(mcPartDaughter1D1->GetPdgCode())==211)) {
 495       AliDebug(2, "The Lambda MC doesn't decay in pi+proton (or cc), skipping!!");
 496       return checkCD;
 497     }
 498
 499     // momentum of the V0 daughters
 500     sumPxDau+=mcPartDaughter1D0->Px();
 501     sumPyDau+=mcPartDaughter1D0->Py();
 502     sumPzDau+=mcPartDaughter1D0->Pz();
 503     sumPxDau+=mcPartDaughter1D1->Px();
 504     sumPyDau+=mcPartDaughter1D1->Py();
 505     sumPzDau+=mcPartDaughter1D1->Pz();
 506
 507   } else if (fGenLcOption==kCountK0Sp) { // Lc -> K0bar + proton AND cc
 508
 509     if (!(TMath::Abs(mcPartDaughter0->GetPdgCode())==311   &&
 510           TMath::Abs(mcPartDaughter1->GetPdgCode())==2212) &&
 511         !(TMath::Abs(mcPartDaughter0->GetPdgCode())==2212  &&
 512           TMath::Abs(mcPartDaughter1->GetPdgCode())==311)) {
 513       AliDebug(2, "The Lc MC doesn't decay in K0+proton (or cc), skipping!!");
 514       return checkCD;
 515     }
 516
 517     if (TMath::Abs(mcPartDaughter0->GetPdgCode())==311) {
 518       mcPartDaughter0 = dynamic_cast<AliAODMCParticle*>(fmcArray->At(daughter1)); // the bachelor
 519       mcPartDaughter1 = dynamic_cast<AliAODMCParticle*>(fmcArray->At(daughter0)); // the V0
 520     }
 521     if (!mcPartDaughter0 || !mcPartDaughter1) {
 522       AliDebug(2,"Problems in the MC Daughters after swapping V0 and bachelor\n");
 523       return checkCD;
 524     }
 525
 526     // momentum of the bachelor
 527     sumPxDau+=mcPartDaughter0->Px();
 528     sumPyDau+=mcPartDaughter0->Py();
 529     sumPzDau+=mcPartDaughter0->Pz();
 530
 531     if (mcPartDaughter1->GetNDaughters()!=1) {
 532       AliDebug(2, "The K0/K0bar MC particle doesn't decay in 1 particles, skipping!!");
 533       return checkCD;
 534     }
 535
 536     Int_t daughter = mcPartDaughter1->GetDaughter(0);
 537     if (daughter<=0) {
 538       AliDebug(2, Form("The K0/K0bar MC particle doesn't have correct daughter, skipping!!"));
 539       return checkCD;
 540     }
 541
 542     AliAODMCParticle* mcPartDaughter = dynamic_cast<AliAODMCParticle*>(fmcArray->At(daughter));
 543     if(!mcPartDaughter){
 544       AliError("The K0/K0bar daughter particle not found in MC array");
 545       return checkCD;
 546     }
 547
 548     if (!(TMath::Abs(mcPartDaughter->GetPdgCode())==310)) {
 549       AliDebug(2, "The K0/K0bar MC doesn't go in K0S, skipping!!");
 550       return checkCD;
 551     }
 552
 553     if (mcPartDaughter->GetNDaughters()!=2) {
 554       AliDebug(2, "The K0S MC doesn't decay in 2 particles, skipping!!");
 555       return checkCD;
 556     }
 557
 558     Int_t daughterD0 = mcPartDaughter->GetDaughter(0);
 559     Int_t daughterD1 = mcPartDaughter->GetDaughter(1);
 560     if (daughterD0<=0 || daughterD1<=0) {
 561       AliDebug(2, Form("The K0S MC particle doesn't have correct daughters, skipping!!"));
 562       return checkCD;
 563     }
 564
 565     AliAODMCParticle* mcPartDaughterD0 = dynamic_cast<AliAODMCParticle*>(fmcArray->At(daughterD0));
 566     AliAODMCParticle* mcPartDaughterD1 = dynamic_cast<AliAODMCParticle*>(fmcArray->At(daughterD1));
 567     if (!mcPartDaughterD0 || !mcPartDaughterD1) {
 568       AliError("Daughter particle not found in MC array");
 569       return checkCD;
 570     }
 571
 572     if (! ( TMath::Abs(mcPartDaughterD0->GetPdgCode())==211 &&
 573             TMath::Abs(mcPartDaughterD1->GetPdgCode())==211 ) ) {
 574       AliDebug(2, "The K0S MC doesn't decay in pi+ pi-, skipping!!");
 575       return checkCD;
 576     }
 577
 578     // momentum of the V0 daughters
 579     sumPxDau+=mcPartDaughterD0->Px();
 580     sumPyDau+=mcPartDaughterD0->Py();
 581     sumPzDau+=mcPartDaughterD0->Pz();
 582     sumPxDau+=mcPartDaughterD1->Px();
 583     sumPyDau+=mcPartDaughterD1->Py();
 584     sumPzDau+=mcPartDaughterD1->Pz();
 585
 586   }
 587
 588   Double_t pxMother=fmcPartCandidate->Px();
 589   Double_t pyMother=fmcPartCandidate->Py();
 590   Double_t pzMother=fmcPartCandidate->Pz();
 591   if(TMath::Abs(pxMother-sumPxDau)/(TMath::Abs(pxMother)+1.e-13)>0.00001 ||
 592      TMath::Abs(pyMother-sumPyDau)/(TMath::Abs(pyMother)+1.e-13)>0.00001 ||
 593      TMath::Abs(pzMother-sumPzDau)/(TMath::Abs(pzMother)+1.e-13)>0.00001){
 594     AliDebug(2, "Momentum conservation violated, skipping!!");
 595     return checkCD;
 596   }
 597
 598   checkCD = kTRUE;
 599   return checkCD;
 600
 601 }
 602
 603 //_____________________________________________________________
 604 Double_t AliCFVertexingHFLctoV0bachelor::GetEtaProng(Int_t iProng) const
 605 {
 606   //
 607   // getting eta of the prong - overload the mother class method
 608   //
 609
 610   Double_t etaProng =-9999;
 611
 612   if (!fRecoCandidate) {
 613     AliDebug(2,"No reco candidate selected");
 614     return etaProng;
 615   }
 616
 617   AliAODRecoCascadeHF* lcV0bachelor = (AliAODRecoCascadeHF*)fRecoCandidate;
 618   AliAODTrack* bachelor = (AliAODTrack*)lcV0bachelor->GetBachelor();
 619   AliAODTrack* v0Pos = (AliAODTrack*)lcV0bachelor->Getv0PositiveTrack();
 620   AliAODTrack* v0Neg = (AliAODTrack*)lcV0bachelor->Getv0NegativeTrack();
 621   if (!(lcV0bachelor->Getv0()) || !bachelor || !v0Pos || !v0Neg) {
 622     AliDebug(2,"No V0 for this reco candidate selected");
 623     return etaProng;
 624   }
 625
 626   if (iProng==0) etaProng = bachelor->Eta();
 627   else if (iProng==1) etaProng = v0Pos->Eta();
 628   else if (iProng==2) etaProng = v0Neg->Eta();
 629
 630   AliDebug(4,Form("Eta value for prong number %1d = %f",iProng,etaProng));
 631
 632   return etaProng;
 633
 634 }
 635
 636 //_____________________________________________________________
 637
 638 Double_t AliCFVertexingHFLctoV0bachelor::GetPtProng(Int_t iProng) const
 639 {
 640   //
 641   // getting pt of the prong
 642   //
 643
 644   Double_t ptProng=-9999.;
 645
 646   if (!fRecoCandidate) {
 647     AliDebug(2,"No reco candidate selected");
 648     return ptProng;
 649   }
 650
 651   AliAODRecoCascadeHF* lcV0bachelor = (AliAODRecoCascadeHF*)fRecoCandidate;
 652   AliAODTrack* bachelor = (AliAODTrack*)lcV0bachelor->GetBachelor();
 653   AliAODTrack* v0Pos = (AliAODTrack*)lcV0bachelor->Getv0PositiveTrack();
 654   AliAODTrack* v0Neg = (AliAODTrack*)lcV0bachelor->Getv0NegativeTrack();
 655   if (!(lcV0bachelor->Getv0()) || !bachelor || !v0Pos || !v0Neg) {
 656     AliDebug(2,"No V0 for this reco candidate selected");
 657     return ptProng;
 658   }
 659
 660   if (iProng==0) ptProng = bachelor->Pt();
 661   else if (iProng==1) ptProng = v0Pos->Pt();
 662   else if (iProng==2) ptProng = v0Neg->Pt();
 663
 664   AliDebug(4,Form("Pt value for prong number %1d = %f",iProng,ptProng));
 665
 666   return ptProng;
 667
 668 }
 669
 670 //_____________________________________________________________
 671
 672 Double_t AliCFVertexingHFLctoV0bachelor::Ctau(AliAODMCParticle *mcPartCandidate)
 673 {
 674
 675   Double_t cTau = 999999.;
 676
 677   Int_t daughterD0 = mcPartCandidate->GetDaughter(0);
 678   Int_t daughterD1 = mcPartCandidate->GetDaughter(1);
 679   if (daughterD0<=0 || daughterD1<=0) {
 680     AliDebug(2, Form("The Lc MC particle doesn't have correct daughters, skipping!!"));
 681     return cTau;
 682   }
 683
 684   AliAODMCParticle *mcPartDaughter0 = dynamic_cast<AliAODMCParticle*>(fmcArray->At(daughterD0));
 685   AliAODMCParticle *mcPartDaughter1 = dynamic_cast<AliAODMCParticle*>(fmcArray->At(daughterD1));
 686   if (!mcPartDaughter0 || !mcPartDaughter1) {
 687     AliDebug(2,"The candidate daughter particles not found in MC array");
 688     return cTau;
 689   }
 690
 691   Double_t vtx1[3] = {0,0,0};   // primary vertex
 692   Bool_t hasProdVertex = mcPartCandidate->XvYvZv(vtx1);  // cm
 693
 694   Double_t vtx1daughter[3] = {0,0,0};   // secondary vertex
 695   Bool_t v0Vertex = mcPartDaughter0->XvYvZv(vtx1daughter);  //cm
 696   Double_t vtx2daughter[3] = {0,0,0};   // secondary vertex
 697   Bool_t bachVertex = hasProdVertex && mcPartDaughter1->XvYvZv(vtx2daughter);  //cm
 698
 699   if (!hasProdVertex || !v0Vertex || !bachVertex) {
 700     AliDebug(2,"At least one of Prim.vtx, V0vtx, BachelorVtx doesn't exist!");
 701     return cTau;
 702   }
 703
 704   if (TMath::Abs(vtx1daughter[0]-vtx2daughter[0])>1E-5 ||
 705       TMath::Abs(vtx1daughter[1]-vtx2daughter[1])>1E-5 ||
 706       TMath::Abs(vtx1daughter[2]-vtx2daughter[2])>1E-5) {
 707     AliDebug(2,"Bachelor and V0 haven't common vtx!");
 708     return cTau;
 709   }
 710
 711   Double_t decayLength = 0.;
 712   for (Int_t ii=0; ii<3; ii++) decayLength += (vtx1daughter[ii]-vtx1[ii])*(vtx1daughter[ii]-vtx1[ii]);
 713   decayLength = TMath::Sqrt(decayLength);
 714
 715   cTau = decayLength * mcPartCandidate->M()/mcPartCandidate->P();
 716
 717   AliDebug(2,Form(" cTau(4122)=%f",cTau));
 718
 719   return cTau;
 720
 721 }
 722
 723 //------------
 724 Bool_t AliCFVertexingHFLctoV0bachelor::SetLabelArray()
 725 {
 726   //
 727   // setting the label arrays
 728   //
 729
 730   Bool_t checkCD = kFALSE;
 731
 732   if (fmcPartCandidate->GetNDaughters()!=2) {
 733     AliDebug(2, Form("The MC particle have %d daughters (not 2), skipping!!",fmcPartCandidate->GetNDaughters()));
 734     //    fmcPartCandidate->Print();
 735     return checkCD;
 736   }
 737
 738   Int_t daughter0 = fmcPartCandidate->GetDaughter(0);
 739   Int_t daughter1 = fmcPartCandidate->GetDaughter(1);
 740   if (daughter0<=0 || daughter1<=0){
 741     AliDebug(2, Form("The MC particle doesn't have correct daughters, skipping!!"));
 742     return checkCD;
 743   }
 744
 745   AliAODMCParticle* mcPartDaughter0 = dynamic_cast<AliAODMCParticle*>(fmcArray->At(daughter0));
 746   AliAODMCParticle* mcPartDaughter1 = dynamic_cast<AliAODMCParticle*>(fmcArray->At(daughter1));
 747   if (!mcPartDaughter0 || !mcPartDaughter1) {
 748     AliDebug(2,"Problems in the MC Daughters\n");
 749     return checkCD;
 750   }
 751
 752
 753   fLabelArray = new Int_t[fProngs];
 754
 755   if (fGenLcOption==kCountLambdapi) { // Lc -> Lambda + pion OR cc
 756
 757     if (!(TMath::Abs(mcPartDaughter0->GetPdgCode())==3122 &&
 758           TMath::Abs(mcPartDaughter1->GetPdgCode())==211) &&
 759         !(TMath::Abs(mcPartDaughter0->GetPdgCode())==211  &&
 760           TMath::Abs(mcPartDaughter1->GetPdgCode())==3122)) {
 761       AliDebug(2, "The Lc MC doesn't decay in Lambda+pion (or cc), skipping!!");
 762       delete [] fLabelArray;
 763       fLabelArray = 0x0;
 764       return checkCD;
 765     }
 766
 767     // it is Lc -> Lambda + pion OR cc
 768     if (TMath::Abs(mcPartDaughter0->GetPdgCode())==3122) {
 769       mcPartDaughter0 = dynamic_cast<AliAODMCParticle*>(fmcArray->At(daughter1)); // the bachelor
 770       mcPartDaughter1 = dynamic_cast<AliAODMCParticle*>(fmcArray->At(daughter0)); // the V0
 771       Int_t daughterTemp = daughter0;
 772       daughter0 = daughter1; // the bachelor label
 773       daughter1 = daughterTemp; // the V0 label
 774     }
 775
 776     if (mcPartDaughter1->GetNDaughters()!=2) {
 777       AliDebug(2, "The Lambda MC particle doesn't decay in 2 particles, skipping!!");
 778       delete [] fLabelArray;
 779       fLabelArray = 0x0;
 780       return checkCD;
 781     }
 782
 783     Int_t daughter1D0 = mcPartDaughter1->GetDaughter(0);
 784     Int_t daughter1D1 = mcPartDaughter1->GetDaughter(1);
 785     if (daughter1D0<=0 || daughter1D1<=0) {
 786       AliDebug(2, Form("The Lambda MC particle doesn't have correct daughters, skipping!!"));
 787       delete [] fLabelArray;
 788       fLabelArray = 0x0;
 789       return checkCD;
 790     }
 791
 792     AliAODMCParticle* mcPartDaughter1D0 = dynamic_cast<AliAODMCParticle*>(fmcArray->At(daughter1D0));
 793     AliAODMCParticle* mcPartDaughter1D1 = dynamic_cast<AliAODMCParticle*>(fmcArray->At(daughter1D1));
 794     if (!mcPartDaughter1D0 || !mcPartDaughter1D1) {
 795       AliError("The Lambda daughter particles not found in MC array");
 796       delete [] fLabelArray;
 797       fLabelArray = 0x0;
 798       return checkCD;
 799     }
 800
 801     if (!(TMath::Abs(mcPartDaughter1D0->GetPdgCode())==211   &&
 802           TMath::Abs(mcPartDaughter1D1->GetPdgCode())==2212) &&
 803         !(TMath::Abs(mcPartDaughter1D0->GetPdgCode())==2212  &&
 804           TMath::Abs(mcPartDaughter1D1->GetPdgCode())==211)) {
 805       AliDebug(2, "The Lambda MC doesn't decay in pi+proton (or cc), skipping!!");
 806       delete [] fLabelArray;
 807       fLabelArray = 0x0;
 808       return checkCD;
 809     }
 810
 811     // Lambda -> p+pi OR cc
 812
 813     fLabelArray[0] = daughter0;//mcPartDaughter0->GetLabel(); // bachelor
 814
 815     if (fmcPartCandidate->Charge()>0) {
 816
 817       if (mcPartDaughter1D0->GetPdgCode()==2212) {
 818         fLabelArray[1] = daughter1D0;//mcPartDaughter1D0->GetLabel(); // proton
 819         fLabelArray[2] = daughter1D1;//mcPartDaughter1D1->GetLabel(); // pion
 820       } else if (mcPartDaughter1D1->GetPdgCode()==2212) {
 821         fLabelArray[1] = daughter1D1;//mcPartDaughter1D1->GetLabel(); // proton
 822         fLabelArray[2] = daughter1D0;//mcPartDaughter1D0->GetLabel(); // pion
 823       }
 824
 825     } else if (fmcPartCandidate->Charge()<0) {
 826
 827       if (mcPartDaughter1D0->GetPdgCode()==211) {
 828         fLabelArray[1] = daughter1D0;//mcPartDaughter1D0->GetLabel(); // pion
 829         fLabelArray[2] = daughter1D1;//mcPartDaughter1D1->GetLabel(); // proton
 830       } else if (mcPartDaughter1D1->GetPdgCode()==211) {
 831         fLabelArray[1] = daughter1D1;//mcPartDaughter1D1->GetLabel(); // pion
 832         fLabelArray[2] = daughter1D0;//mcPartDaughter1D0->GetLabel(); // proton
 833       }
 834
 835     }
 836
 837   } else if (fGenLcOption==kCountK0Sp) { // Lc -> K0bar + proton OR cc
 838
 839     if (!(TMath::Abs(mcPartDaughter0->GetPdgCode())==311   &&
 840           TMath::Abs(mcPartDaughter1->GetPdgCode())==2212) &&
 841         !(TMath::Abs(mcPartDaughter0->GetPdgCode())==2212  &&
 842           TMath::Abs(mcPartDaughter1->GetPdgCode())==311)) {
 843       AliDebug(2, "The Lc MC doesn't decay in K0bar+proton (or cc), skipping!!");
 844       delete [] fLabelArray;
 845       fLabelArray = 0x0;
 846       return checkCD;
 847     }
 848
 849     if (TMath::Abs(mcPartDaughter0->GetPdgCode())==311) {
 850       mcPartDaughter0 = dynamic_cast<AliAODMCParticle*>(fmcArray->At(daughter1)); // the bachelor
 851       mcPartDaughter1 = dynamic_cast<AliAODMCParticle*>(fmcArray->At(daughter0)); // the V0
 852       Int_t daughterTemp = daughter0;
 853       daughter0 = daughter1; // the bachelor label
 854       daughter1 = daughterTemp; // the V0 label
 855     }
 856
 857     if(!mcPartDaughter1){
 858       AliDebug(2, "Null pointer to K0 particle, skipping!!");
 859       delete [] fLabelArray;
 860       fLabelArray = 0x0;
 861       return checkCD;
 862     }
 863
 864     if (mcPartDaughter1->GetNDaughters()!=1) {
 865       AliDebug(2, "The K0/K0bar MC particle doesn't decay in 1 particles, skipping!!");
 866       delete [] fLabelArray;
 867       fLabelArray = 0x0;
 868       return checkCD;
 869     }
 870
 871     Int_t daughter = mcPartDaughter1->GetDaughter(0);
 872     if (daughter<=0) {
 873       AliDebug(2, Form("The K0/K0bar MC particle doesn't have correct daughter, skipping!!"));
 874       delete [] fLabelArray;
 875       fLabelArray = 0x0;
 876       return checkCD;
 877     }
 878
 879     AliAODMCParticle* mcPartDaughter = dynamic_cast<AliAODMCParticle*>(fmcArray->At(daughter));
 880     if (!mcPartDaughter) {
 881       AliError("The K0/K0bar daughter particle not found in MC array");
 882       delete [] fLabelArray;
 883       fLabelArray = 0x0;
 884       return checkCD;
 885     }
 886
 887     if (!(TMath::Abs(mcPartDaughter->GetPdgCode())==310)) {
 888       AliDebug(2, "The K0/K0bar MC doesn't go in K0S, skipping!!");
 889       delete [] fLabelArray;
 890       fLabelArray = 0x0;
 891       return checkCD;
 892     }
 893
 894     if (mcPartDaughter->GetNDaughters()!=2) {
 895       AliDebug(2, "The K0S MC doesn't decay in 2 particles, skipping!!");
 896       delete [] fLabelArray;
 897       fLabelArray = 0x0;
 898       return checkCD;
 899     }
 900
 901     Int_t daughterD0 = mcPartDaughter->GetDaughter(0);
 902     Int_t daughterD1 = mcPartDaughter->GetDaughter(1);
 903     if (daughterD0<=0 || daughterD1<=0) {
 904       AliDebug(2, Form("The K0S MC particle doesn't have correct daughters, skipping!!"));
 905       delete [] fLabelArray;
 906       fLabelArray = 0x0;
 907       return checkCD;
 908     }
 909
 910     AliAODMCParticle* mcPartDaughterD0 = dynamic_cast<AliAODMCParticle*>(fmcArray->At(daughterD0));
 911     AliAODMCParticle* mcPartDaughterD1 = dynamic_cast<AliAODMCParticle*>(fmcArray->At(daughterD1));
 912     if (!mcPartDaughterD0 || !mcPartDaughterD1) {
 913       AliError("The K0S daughter particles not found in MC array");
 914       delete [] fLabelArray;
 915       fLabelArray = 0x0;
 916       return checkCD;
 917     }
 918
 919     if (! ( TMath::Abs(mcPartDaughterD0->GetPdgCode())==211 &&
 920             TMath::Abs(mcPartDaughterD1->GetPdgCode())==211 ) ) {
 921       AliDebug(2, "The K0S MC doesn't decay in pi+ pi-, skipping!!");
 922       delete [] fLabelArray;
 923       fLabelArray = 0x0;
 924       return checkCD;
 925     }
 926
 927     // K0S -> pi+ pi-
 928
 929     fLabelArray[0] = daughter0;//mcPartDaughter0->GetLabel(); // bachelor
 930
 931     if (mcPartDaughterD0->GetPdgCode()==211) {
 932       fLabelArray[1] = daughterD0;//mcPartDaughterD0->GetLabel(); // pi+
 933       fLabelArray[2] = daughterD1;//mcPartDaughterD1->GetLabel(); // pi-
 934       AliDebug(2,Form(" daughter0=%d ------ daughter1=%d ------ dg0->GetLabel()=%d ------ dg1->GetLabel()=%d ",daughterD0,daughterD1,mcPartDaughterD0->GetLabel(),mcPartDaughterD1->GetLabel()));
 935     } else if (mcPartDaughterD1->GetPdgCode()==211) {
 936       fLabelArray[1] = daughterD1;//mcPartDaughterD1->GetLabel(); // pi+
 937       fLabelArray[2] = daughterD0;//mcPartDaughterD0->GetLabel(); // pi-
 938       AliDebug(2,Form(" daughter0=%d ------ daughter1=%d ------ dg0->GetLabel()=%d ------ dg1->GetLabel()=%d ",daughterD1,daughterD0,mcPartDaughterD1->GetLabel(),mcPartDaughterD0->GetLabel()));
 939     }
 940   }
 941
 942   AliDebug(2,Form(" label0=%d, label1=%d, label2=%d",fLabelArray[0],fLabelArray[1],fLabelArray[2]));
 943
 944   SetAccCut(); // setting the pt and eta acceptance cuts
 945
 946   checkCD = kTRUE;
 947   return checkCD;
 948
 949 }
 950 //____________________________________________
 951 Bool_t AliCFVertexingHFLctoV0bachelor::FillVectorFromMCarray(AliAODMCParticle *mcPartDaughterBachelor,
 952                                                              AliAODMCParticle *mcPartDaughterK0,
 953                                                              Double_t *vectorMC)
 954 {
 955   // fill the vector
 956
 957   Bool_t bGenValues = kFALSE;
 958
 959   AliAODMCParticle *mcPartV0DaughterPos = dynamic_cast<AliAODMCParticle*>(fmcArray->At(fLabelArray[1]));
 960   AliAODMCParticle *mcPartV0DaughterNeg = dynamic_cast<AliAODMCParticle*>(fmcArray->At(fLabelArray[2]));
 961   AliAODMCParticle *mcPartDaughterV0 = 0x0;
 962
 963   if(!mcPartV0DaughterPos || !mcPartV0DaughterNeg) return bGenValues;
 964
 965   if (TMath::Abs(mcPartDaughterK0->GetPdgCode())==311) {
 966     Int_t daughterK0 = mcPartDaughterK0->GetDaughter(0);
 967     if (daughterK0<=0) {
 968       AliDebug(2, Form("The K0/K0bar particle doesn't have correct daughter, skipping!!"));
 969       return bGenValues;
 970     }
 971     mcPartDaughterV0 = dynamic_cast<AliAODMCParticle*>(fmcArray->At(daughterK0));
 972     if (!mcPartDaughterV0) {
 973       AliDebug(2,"The K0/K0bar daughter particle not found in MC array");
 974       return bGenValues;
 975     }
 976     if (TMath::Abs(mcPartDaughterV0->GetPdgCode())!=310) {
 977       AliDebug(2,"The K0/K0bar daughter particle is not a K0S");
 978       return bGenValues;
 979     }
 980   } else if (TMath::Abs(mcPartDaughterK0->GetPdgCode())==3122) {
 981     mcPartDaughterV0 = dynamic_cast<AliAODMCParticle*>(mcPartDaughterK0);
 982     if (!mcPartDaughterV0) {
 983       AliDebug(2,"The Lambda particle not found in MC array");
 984       return bGenValues;
 985     }
 986   }
 987
 988   if (!mcPartDaughterV0) {
 989     AliDebug(2,"V0 particle not found in MC array");
 990     return bGenValues;
 991   }
 992
 993   //Double_t cTLc = Ctau(fmcPartCandidate); // by default wrt Primary Vtx
 994   //Double_t pTbach = mcPartDaughterBachelor->Pt(); // get the bachelor pT
 995   Double_t pbach = mcPartDaughterBachelor->P(); // get the bachelor p
 996
 997   Double_t vtx1[3] = {0,0,0};   // primary vertex
 998   Bool_t hasPrimVtx = fmcPartCandidate->XvYvZv(vtx1);  // cm
 999
1000   // getting vertex from daughters
1001   Double_t vtx1daughter0[3] = {0,0,0};   // secondary vertex from daughter 0
1002   Bool_t hasSecVtx1 = mcPartDaughterBachelor->XvYvZv(vtx1daughter0);  //cm
1003   Double_t vtx1daughter1[3] = {0,0,0};   // secondary vertex from daughter 1
1004   Bool_t hasSecVtx2 = mcPartDaughterV0->XvYvZv(vtx1daughter1);  //cm
1005   if (!hasPrimVtx || !hasSecVtx1 || !hasSecVtx2) {
1006     AliDebug(2,"At least one of Prim.vtx, V0vtx, BachelorVtx doesn't exist!");
1007     //return bGenValues;
1008   }
1009
1010   if (TMath::Abs(vtx1daughter0[0]-vtx1daughter1[0])>1E-5 ||
1011       TMath::Abs(vtx1daughter0[1]-vtx1daughter1[1])>1E-5 ||
1012       TMath::Abs(vtx1daughter0[2]-vtx1daughter1[2])>1E-5) {
1013     AliError("Daughters have different secondary vertex, skipping the track");
1014     //return bGenValues;
1015   }
1016
1017   // getting the momentum from the daughters
1018   Double_t px1[2] = {mcPartDaughterBachelor->Px(), mcPartDaughterV0->Px()};
1019   Double_t py1[2] = {mcPartDaughterBachelor->Py(), mcPartDaughterV0->Py()};
1020   Double_t pz1[2] = {mcPartDaughterBachelor->Pz(), mcPartDaughterV0->Pz()};
1021
1022   Int_t nprongs = 2;
1023   Short_t charge = mcPartDaughterBachelor->Charge();
1024   Double_t d0[2] = {0.,0.};
1025   AliAODRecoDecayHF* decayLc = new AliAODRecoDecayHF(vtx1,vtx1daughter0,nprongs,charge,px1,py1,pz1,d0);
1026   Double_t cosPAwrtPrimVtxLc = decayLc->CosPointingAngle();
1027   delete decayLc;
1028
1029   // getting vertex from daughters
1030   Double_t vtx2daughter0[3] = {0,0,0};   // secondary vertex from daughter 0
1031   Bool_t hasSecVtx3 = mcPartV0DaughterPos->XvYvZv(vtx2daughter0);  //cm
1032   Double_t vtx2daughter1[3] = {0,0,0};   // secondary vertex from daughter 1
1033   Bool_t hasSecVtx4 = mcPartV0DaughterNeg->XvYvZv(vtx2daughter1);  //cm
1034   if (!hasSecVtx3 || !hasSecVtx4) {
1035     AliDebug(2,"At least one of V0Posvtx, V0Negtx doesn't exist!");
1036     //return bGenValues;
1037   }
1038
1039   if (TMath::Abs(vtx2daughter0[0]-vtx2daughter1[0])>1E-5 ||
1040       TMath::Abs(vtx2daughter0[1]-vtx2daughter1[1])>1E-5 ||
1041       TMath::Abs(vtx2daughter0[2]-vtx2daughter1[2])>1E-5) {
1042     AliError("Daughters have different secondary vertex, skipping the track");
1043     //return bGenValues;
1044   }
1045
1046   // getting the momentum from the daughters
1047   Double_t px[2] = {mcPartV0DaughterPos->Px(), mcPartV0DaughterNeg->Px()};
1048   Double_t py[2] = {mcPartV0DaughterPos->Py(), mcPartV0DaughterNeg->Py()};
1049   Double_t pz[2] = {mcPartV0DaughterPos->Pz(), mcPartV0DaughterNeg->Pz()};
1050
1051   nprongs = 2;
1052   charge = 0;
1053   AliAODRecoDecayHF* decay = new AliAODRecoDecayHF(vtx1,vtx2daughter0,nprongs,charge,px,py,pz,d0);
1054   Double_t cosPAwrtPrimVtxV0 = decay->CosPointingAngle();
1055   /*
1056   Double_t cTV0 = 0.; //ct
1057   if (fGenLcOption==kCountK0Sp) {
1058     cTV0 = decay->Ct(310); // by default wrt Primary Vtx
1059   } else if (fGenLcOption==kCountLambdapi) {
1060     cTV0 = decay->Ct(3122); // by default wrt Primary Vtx
1061   }
1062   */
1063
1064   Double_t invMass = 0.; //invMass
1065   if (fGenLcOption==kCountK0Sp) {
1066     invMass = decay->InvMass2Prongs(0,1,211,211);
1067   } else if (fGenLcOption==kCountLambdapi) {
1068     if (fmcPartCandidate->GetPdgCode() == 4122)
1069       invMass = decay->InvMass2Prongs(0,1,2212,211);
1070     else if (fmcPartCandidate->GetPdgCode() ==-4122)
1071       invMass = decay->InvMass2Prongs(0,1,211,2212);
1072   }
1073   delete decay;
1074
1075   vectorMC[0]  = fmcPartCandidate->Pt();
1076   vectorMC[1]  = fmcPartCandidate->Y() ;
1077   vectorMC[2]  = fmcPartCandidate->Phi();
1078   vectorMC[3]  = 0; // dummy value x MC, onTheFlyStatus
1079   vectorMC[4]  = fzMCVertex;
1080   vectorMC[5]  = fCentValue; // reconstructed centrality
1081   vectorMC[6]  = 1; // dummy value x MC, fFake
1082   vectorMC[7]  = fMultiplicity; // reconstructed multiplicity
1083
1084   if (fConfiguration==AliCFTaskVertexingHF::kSnail) {
1085     //vectorMC[8]  = pTbach;
1086     vectorMC[8]  = pbach;
1087     vectorMC[9]  = mcPartDaughterV0->Pt();
1088     vectorMC[10] = mcPartDaughterV0->Y();
1089     vectorMC[11] = mcPartDaughterV0->Phi();
1090     vectorMC[12] = invMass;
1091     vectorMC[13] = 0; // dummy value x MC, V0 DCA
1092     vectorMC[14] = cosPAwrtPrimVtxV0;
1093     vectorMC[15] = cosPAwrtPrimVtxLc;
1094     //vectorMC[16] = cTV0*1.E4; // in micron
1095     //vectorMC[17] = cTLc*1.E4; // in micron
1096   }
1097
1098   bGenValues = kTRUE;
1099   return bGenValues;
1100
1101 }