[u/mrichter/AliRoot.git] / PWGHF / vertexingHF / AliCFVertexingHFCascade.cxx

/**************************************************************************
 * Copyright(c) 1998-2011, ALICE Experiment at CERN, All rights reserved. *
 *                                                                        *
 * Author: The ALICE Off-line Project.                                    *
 * Contributors are mentioned in the code where appropriate.              *
 *                                                                        *
 * Permission to use, copy, modify and distribute this software and its   *
 * documentation strictly for non-commercial purposes is hereby granted   *
 * without fee, provided that the above copyright notice appears in all   *
 * copies and that both the copyright notice and this permission notice   *
 * appear in the supporting documentation. The authors make no claims     *
 * about the suitability of this software for any purpose. It is          *
 * provided "as is" without express or implied warranty.                  *
 **************************************************************************/

//-----------------------------------------------------------------------
// Class for HF corrections as a function of many variables and steps
// For D* and other cascades
// 
// Author : A.Grelli a.grelli@uu.nl  UTECHT
//-----------------------------------------------------------------------

#include "AliAODRecoDecayHF2Prong.h"
#include "AliAODMCParticle.h"
#include "AliAODEvent.h"
#include "TClonesArray.h"
#include "AliCFVertexingHF.h"
#include "AliESDtrack.h"
#include "TDatabasePDG.h"
#include "AliAODRecoCascadeHF.h"
#include "AliCFVertexingHFCascade.h"
#include "AliCFContainer.h"
#include "AliCFTaskVertexingHF.h"

ClassImp(AliCFVertexingHFCascade)


//_________________________________________
  AliCFVertexingHFCascade::AliCFVertexingHFCascade(TClonesArray *mcArray, UShort_t originDselection):
  AliCFVertexingHF(mcArray, originDselection)
{	
  // standard constructor

  SetNProngs(3);
  fPtAccCut=new Float_t[fProngs];
  fEtaAccCut=new Float_t[fProngs];
  // element 0 in the cut arrays corresponds to the soft pion!!!!!!!! Careful when setting the values...
  fPtAccCut[0]=0.;
  fEtaAccCut[0]=0.;
  for(Int_t iP=1; iP<fProngs; iP++){
	  fPtAccCut[iP]=0.1;
	  fEtaAccCut[iP]=0.9;
  }

}


//_____________________________________
AliCFVertexingHFCascade& AliCFVertexingHFCascade::operator=(const AliCFVertexingHFCascade& c)
{
  // operator =
 
  if  (this != &c) {

    AliCFVertexingHF::operator=(c);
   
  }
  return *this;
}

//__________________________________________
Bool_t AliCFVertexingHFCascade::SetRecoCandidateParam(AliAODRecoDecayHF *recoCand)
{
  // set the AliAODRecoDecay candidate
  
  Bool_t bSignAssoc = kFALSE;
 
  fRecoCandidate = recoCand;
  AliAODRecoCascadeHF* dstarD0pi = (AliAODRecoCascadeHF*)recoCand;

  if (!fRecoCandidate) {
    AliError("fRecoCandidate not found, problem in assignement\n");
    return bSignAssoc;
  }
  
  if ( fRecoCandidate->GetPrimaryVtx()) AliDebug(3,"fReco Candidate has a pointer to PrimVtx\n");
  //if (recoCand->GetPrimaryVtx()) printf("Reco Cand has a pointer to PrimVtx\n");
  
  //Int_t pdgCand = 413;

  Int_t pdgDgDStartoD0pi[2]={421,211};
  Int_t pdgDgD0toKpi[2]={321,211};

  Int_t nentries = fmcArray->GetEntriesFast();

  AliDebug(3,Form("nentries = %d\n", nentries));
 
  Int_t mcLabel =  dstarD0pi->MatchToMC(413,421,pdgDgDStartoD0pi,pdgDgD0toKpi,fmcArray); 
  
  if (mcLabel == -1) return bSignAssoc;

  if (fRecoCandidate->NumberOfFakeDaughters()>0){
	  fFake = 0;    // fake candidate
	  if (fFakeSelection==1) return bSignAssoc;
  }
  if (fRecoCandidate->NumberOfFakeDaughters()==0){
	  fFake = 2;    // non-fake candidate
	  if (fFakeSelection==2) return bSignAssoc;
  }
  
  SetMCLabel(mcLabel);
  fmcPartCandidate = dynamic_cast<AliAODMCParticle*>(fmcArray->At(fmcLabel)); 

  if (!fmcPartCandidate){
    AliDebug(3,"No part candidate");
    return bSignAssoc;
  }
 
  bSignAssoc = kTRUE;
  return bSignAssoc;
}

//______________________________________________
Bool_t AliCFVertexingHFCascade::GetGeneratedValuesFromMCParticle(Double_t* vectorMC) 
{
  // 
  // collecting all the necessary info (pt, y, cosThetaStar, ptPi, ptKa, cT) from MC particle
  //
	
	Bool_t bGenValues = kFALSE;
	
	
	Int_t daughter0ds = fmcPartCandidate->GetDaughter(0);
	Int_t daughter1ds = fmcPartCandidate->GetDaughter(1);

        //the D0
	AliAODMCParticle* mcPartDaughterD0 = dynamic_cast<AliAODMCParticle*>(fmcArray->At(daughter0ds));
	AliAODMCParticle* mcPartDaughterPis = dynamic_cast<AliAODMCParticle*>(fmcArray->At(daughter1ds));

	if(!mcPartDaughterD0 || !mcPartDaughterPis) return kFALSE;

	Double_t vtx1[3] = {0,0,0};   // primary vertex		
	Double_t vtx2daughter0[3] = {0,0,0};   // secondary vertex from daughter 0
	Double_t vtx2daughter1[3] = {0,0,0};   // secondary vertex from daughter 1
	fmcPartCandidate->XvYvZv(vtx1);  // cm

	//D0 daughters
	Int_t daughter0 = mcPartDaughterD0->GetDaughter(0);
	Int_t daughter1 = mcPartDaughterD0->GetDaughter(1);

	AliAODMCParticle* mcPartDaughter0 = dynamic_cast<AliAODMCParticle*>(fmcArray->At(daughter0)); //D0
	AliAODMCParticle* mcPartDaughter1 = dynamic_cast<AliAODMCParticle*>(fmcArray->At(daughter1)); //pis

	if(!mcPartDaughter0 || !mcPartDaughter1) return kFALSE;

	// getting vertex from daughters
	mcPartDaughter0->XvYvZv(vtx2daughter0);  // cm
	mcPartDaughter1->XvYvZv(vtx2daughter1);  //cm
	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) {
	  AliError("Daughters have different secondary vertex, skipping the track");
	  return bGenValues;
	}
	
	Int_t nprongs = 2;
	Short_t charge = 0;
	// always instantiate the AliAODRecoDecay with the positive daughter first, the negative second
	AliAODMCParticle* positiveDaugh = mcPartDaughter0;
	AliAODMCParticle* negativeDaugh = mcPartDaughter1;
	if (mcPartDaughter0->GetPdgCode()<0 && mcPartDaughter1->GetPdgCode()>0){
		// inverting in case the positive daughter is the second one
		positiveDaugh = mcPartDaughter1;
		negativeDaugh = mcPartDaughter0;
	}
	// getting the momentum from the daughters
	Double_t px[2] = {positiveDaugh->Px(), negativeDaugh->Px()};		
	Double_t py[2] = {positiveDaugh->Py(), negativeDaugh->Py()};		
	Double_t pz[2] = {positiveDaugh->Pz(), negativeDaugh->Pz()};
	
	Double_t d0[2] = {0.,0.};		
	
	AliAODRecoDecayHF* decay = new AliAODRecoDecayHF(vtx1,vtx2daughter0,nprongs,charge,px,py,pz,d0);
	
	Double_t cosThetaStar = 0.;
	Double_t cosThetaStarD0 = 0.;
	Double_t cosThetaStarD0bar = 0.;
	cosThetaStarD0 = decay->CosThetaStar(1,421,211,321);
	cosThetaStarD0bar = decay->CosThetaStar(0,421,321,211);
	if (mcPartDaughterD0->GetPdgCode() == 421){  // D0
	  AliDebug(3, Form("D0, with pdgprong0 = %d, pdgprong1 = %d",mcPartDaughter0->GetPdgCode(),mcPartDaughter1->GetPdgCode()));
	  cosThetaStar = cosThetaStarD0;
	}
	else if (mcPartDaughterD0->GetPdgCode() == -421){  // D0bar{
	  AliDebug(3, Form("D0bar, with pdgprong0 = %d, pdgprong1 = %d",mcPartDaughter0->GetPdgCode(),mcPartDaughter1->GetPdgCode()));
	  cosThetaStar = cosThetaStarD0bar;
	}
	else{
	  AliWarning("There are problems!! particle was expected to be either a D0 or a D0bar, check...");
	  delete decay;
	  return vectorMC;
	}
	if (cosThetaStar < -1 || cosThetaStar > 1) {
	  AliWarning("Invalid value for cosine Theta star, returning");
	  delete decay;
	  return bGenValues;
	}
	
	Double_t vectorD0[2] ={0.,0.};

	// evaluate the correct cascade
	if (!EvaluateIfD0toKpi(mcPartDaughterD0,vectorD0)) {
	  AliDebug(2, "Error! the D0 MC doesn't have correct daughters!!");
	  delete decay;
	  return bGenValues;  
	}	

	//ct
	Double_t cT = decay->Ct(421);
	// get the pT of the daughters
	Double_t pTD0 = 0.;
	Double_t pTpi = 0.;
	
	if (TMath::Abs(fmcPartCandidate->GetPdgCode()) == 413) {
		pTD0 = mcPartDaughterD0->Pt();
		pTpi = mcPartDaughterPis->Pt();
	}

	
	switch (fConfiguration){
	case AliCFTaskVertexingHF::kSnail:
		vectorMC[0] = fmcPartCandidate->Pt();
		vectorMC[1] = fmcPartCandidate->Y() ;
		vectorMC[2] = cosThetaStar ;
		vectorMC[3] = vectorD0[0]; 
		vectorMC[4] = vectorD0[1];
		vectorMC[5] = cT*1.E4 ;  // in micron
		vectorMC[6] = 0.;   // dummy value, meaningless in MC
		vectorMC[7] = -100000.; // dummy value, meaningless in MC, in micron^2
		vectorMC[8] = 1.01;    // dummy value, meaningless in MC
		vectorMC[9] = fmcPartCandidate->Phi(); 
		vectorMC[10] = fzMCVertex;    // z of reconstructed of primary vertex
		vectorMC[11] = fCentValue; // reconstructed centrality
		vectorMC[12] = 1.;           // always filling with 1 at MC level 
		vectorMC[13] = 1.01; // dummy value for cosPointingXY  multiplicity
		vectorMC[14] = 0.; // dummy value for NormalizedDecayLengthXY multiplicity
		vectorMC[15] = fMultiplicity; // reconstructed multiplicity
		break;
	case AliCFTaskVertexingHF::kCheetah:
		vectorMC[0] = fmcPartCandidate->Pt();
		vectorMC[1] = fmcPartCandidate->Y() ;
		vectorMC[2] = cT*1.E4; // in micron
		vectorMC[3] = fmcPartCandidate->Phi();
		vectorMC[4] = fzMCVertex;
		vectorMC[5] = fCentValue;   // dummy value for dca, meaningless in MC
		vectorMC[6] = 1. ;  // fake: always filling with 1 at MC level 
		vectorMC[7] = fMultiplicity;   // dummy value for d0pi, meaningless in MC, in micron
		break;
	}

	delete decay;
	bGenValues = kTRUE;
	return bGenValues;
}
//____________________________________________
Bool_t AliCFVertexingHFCascade::GetRecoValuesFromCandidate(Double_t *vectorReco) const
{ 
  // read the variables for the container

  Bool_t bFillRecoValues=kFALSE;
  
  //Get the D* and the D0 from D*
  AliAODRecoCascadeHF* dstarD0pi = (AliAODRecoCascadeHF*)fRecoCandidate;
  AliAODRecoDecayHF2Prong* d0toKpi = (AliAODRecoDecayHF2Prong*)dstarD0pi->Get2Prong();
  

  //if (dstarD0pi->GetPrimaryVtx())printf("dstarD0pi has primary vtx\n");
  //if (fRecoCandidate->GetPrimaryVtx())printf("fRecoCandidateDstar has primary vtx\n");

  Double_t pt =  dstarD0pi->Pt();
  Double_t rapidity =  dstarD0pi->YDstar();
  Double_t invMass=0.;
  Double_t cosThetaStar = 9999.;
  Double_t pTpi = 0.;
  Double_t pTK = 0.;
  Double_t dca = d0toKpi->GetDCA();
  Double_t d0pi = 0.;
  Double_t d0K = 0.;
  Double_t d0xd0 = d0toKpi->Prodd0d0();
  Double_t cosPointingAngle = d0toKpi->CosPointingAngle();
  Double_t phi = dstarD0pi->Phi();
  Double_t cosPointingAngleXY = d0toKpi->CosPointingAngleXY();
  Double_t normDecayLengthXY = d0toKpi->NormalizedDecayLengthXY();

  Int_t pdgCode = fmcPartCandidate->GetPdgCode();
 
  if (pdgCode > 0){
    cosThetaStar = d0toKpi->CosThetaStarD0();
    pTpi = d0toKpi->PtProng(0);
    pTK = d0toKpi->PtProng(1);
    d0pi = d0toKpi->Getd0Prong(0);
    d0K = d0toKpi->Getd0Prong(1);
    invMass=d0toKpi->InvMassD0();
  }
  else {
    cosThetaStar = d0toKpi->CosThetaStarD0bar();
    pTpi = d0toKpi->PtProng(1);
    pTK = d0toKpi->PtProng(0);
    d0pi = d0toKpi->Getd0Prong(1);
    d0K = d0toKpi->Getd0Prong(0);
    invMass= d0toKpi->InvMassD0bar();
  }
  
  Double_t cT = d0toKpi->CtD0();
  
	switch (fConfiguration){
	case AliCFTaskVertexingHF::kSnail:
		vectorReco[0] = pt;
		vectorReco[1] = rapidity;
		vectorReco[2] = cosThetaStar;
		vectorReco[3] = pTpi;
		vectorReco[4] = pTK;
		vectorReco[5] = cT*1.E4;  // in micron
		vectorReco[6] = dca*1.E4;  // in micron
		vectorReco[7] = d0xd0*1.E8;  // in micron^2
		vectorReco[8] = cosPointingAngle;  // in micron
		vectorReco[9] = phi;  
		vectorReco[10] = fzPrimVertex;    // z of reconstructed of primary vertex
		vectorReco[11] = fCentValue;
		vectorReco[12] = fFake;      // whether the reconstructed candidate was a fake (fFake = 0) or not (fFake = 2) 
		vectorReco[13] = cosPointingAngleXY; 
		vectorReco[14] = normDecayLengthXY; // in cm
		vectorReco[15] = fMultiplicity; // reconstructed multiplicity
		break;
	case AliCFTaskVertexingHF::kCheetah:
		vectorReco[0] = pt;
		vectorReco[1] = rapidity ;
		vectorReco[2] = cT*1.E4; // in micron
		vectorReco[3] = phi; 
		vectorReco[4] = fzPrimVertex;
		vectorReco[5] = fCentValue;   
		vectorReco[6] = fFake ; 
		vectorReco[7] = fMultiplicity;  
		break;
	}

	bFillRecoValues = kTRUE;
		
  return bFillRecoValues;
}


//_____________________________________________________________
Bool_t AliCFVertexingHFCascade::CheckMCChannelDecay() const
{ 
  // check the required decay channel

  Bool_t checkCD = kFALSE;

  
  Int_t daughter0 = fmcPartCandidate->GetDaughter(0);
  Int_t daughter1 = fmcPartCandidate->GetDaughter(1);
  AliAODMCParticle* mcPartDaughter0 = dynamic_cast<AliAODMCParticle*>(fmcArray->At(daughter0));
  AliAODMCParticle* mcPartDaughter1 = dynamic_cast<AliAODMCParticle*>(fmcArray->At(daughter1));

  if (!mcPartDaughter0 || !mcPartDaughter1) {
    AliDebug (2,"Problems in the MC Daughters\n");
    return checkCD;
  }

  if (!(TMath::Abs(mcPartDaughter0->GetPdgCode())==421 &&
	TMath::Abs(mcPartDaughter1->GetPdgCode())==211) && 
      !(TMath::Abs(mcPartDaughter0->GetPdgCode())==211 &&
	TMath::Abs(mcPartDaughter1->GetPdgCode())==421)) {
    AliDebug(2, "The D0 MC doesn't come from a Kpi decay, skipping!!");
    return checkCD;  
  }
  
  Double_t vectorD0[2] ={0.,0.};

  // D* is a cascade ...evaluate the correct cascade
  if (!EvaluateIfD0toKpi(mcPartDaughter0,vectorD0)) {
    AliDebug(2, "Error! the D0 MC doesn't have correct daughters!!");
    return checkCD;  
  }
   
  checkCD = kTRUE;
  return checkCD;
  
}

//__________________________________________
Bool_t AliCFVertexingHFCascade::EvaluateIfD0toKpi(AliAODMCParticle* neutralDaugh, Double_t* vectorD0)const
{  
  //
  // chack wether D0 is decaing into kpi
  //
  
  Bool_t isHadronic = kFALSE;
  
  Int_t daughterD00 = neutralDaugh->GetDaughter(0);
  Int_t daughterD01 = neutralDaugh->GetDaughter(1);
  
  AliDebug(2, Form("daughter0 = %d and daughter1 = %d",daughterD00,daughterD01));
  if (daughterD00 == 0 || daughterD01 == 0) {
    AliDebug(2, "Error! the D0 MC doesn't have correct daughters!!");
    return isHadronic;  
  }
  
  if (TMath::Abs(daughterD01 - daughterD00) != 1) { // should be everytime true - see PDGdatabooklet
    AliDebug(2, "The D0 MC doesn't come from a 2-prong decay, skipping!!");
    return isHadronic;  
  }
  
  AliAODMCParticle* mcPartDaughterD00 = dynamic_cast<AliAODMCParticle*>(fmcArray->At(daughterD00));
  AliAODMCParticle* mcPartDaughterD01 = dynamic_cast<AliAODMCParticle*>(fmcArray->At(daughterD01));
  if (!mcPartDaughterD00 || !mcPartDaughterD01) {
    AliWarning("D0 MC analysis: At least one Daughter Particle not found in tree, skipping"); 
    return isHadronic;  
  }
  
  if (!(TMath::Abs(mcPartDaughterD00->GetPdgCode())==321 &&
	TMath::Abs(mcPartDaughterD01->GetPdgCode())==211) && 
      !(TMath::Abs(mcPartDaughterD00->GetPdgCode())==211 &&
	TMath::Abs(mcPartDaughterD01->GetPdgCode())==321)) {
    AliDebug(2, "The D0 MC doesn't come from a Kpi decay, skipping!!");
    return isHadronic;  
  }
  
  Double_t sumPxDau=mcPartDaughterD00->Px()+mcPartDaughterD01->Px();
  Double_t sumPyDau=mcPartDaughterD00->Py()+mcPartDaughterD01->Py();
  Double_t sumPzDau=mcPartDaughterD00->Pz()+mcPartDaughterD01->Pz();
  Double_t pxMother=neutralDaugh->Px();
  Double_t pyMother=neutralDaugh->Py();
  Double_t pzMother=neutralDaugh->Pz();
  if(TMath::Abs(pxMother-sumPxDau)/(TMath::Abs(pxMother)+1.e-13)>0.00001 ||
     TMath::Abs(pyMother-sumPyDau)/(TMath::Abs(pyMother)+1.e-13)>0.00001 ||
     TMath::Abs(pzMother-sumPzDau)/(TMath::Abs(pzMother)+1.e-13)>0.00001){
    AliDebug(2, "Momentum conservation violated, skipping!!");
    return isHadronic;  
  }

  Double_t pTD0pi = 0;
  Double_t pTD0K = 0;
  
  
  if (TMath::Abs(mcPartDaughterD00->GetPdgCode()) == 211) {
    pTD0pi = mcPartDaughterD00->Pt();
    pTD0K = mcPartDaughterD01->Pt();
  }
  else {
    pTD0pi = mcPartDaughterD01->Pt();
    pTD0K  = mcPartDaughterD00->Pt();
  }
  
  isHadronic = kTRUE;
  
  vectorD0[0] = pTD0pi;
  vectorD0[1] = pTD0K;
 
  return isHadronic;

}

//___________________________________________________________

void AliCFVertexingHFCascade::SetPtAccCut(Float_t* ptAccCut)
{
	//
	// setting the pt cut to be used in the Acceptance steps (MC+Reco)
	//

	AliInfo("The 3rd element of the pt cut array will correspond to the cut applied to the soft pion - please check that it is correct");
	if (fProngs>0){
		for (Int_t iP=0; iP<fProngs; iP++){
			fPtAccCut[iP]=ptAccCut[iP];
		}
	}
	return;
}		


//___________________________________________________________

void AliCFVertexingHFCascade::SetEtaAccCut(Float_t* etaAccCut)
{
	//
	// setting the eta cut to be used in the Acceptance steps (MC+Reco)
	//

	AliInfo("The 3rd element of the eta cut array will correspond to the cut applied to the soft pion - please check that it is correct");
	if (fProngs>0){
		for (Int_t iP=0; iP<fProngs; iP++){
			fEtaAccCut[iP]=etaAccCut[iP];
		}
	}
	return;
}	
//___________________________________________________________

void AliCFVertexingHFCascade::SetAccCut(Float_t* ptAccCut, Float_t* etaAccCut)
{
	//
	// setting the pt and eta cut to be used in the Acceptance steps (MC+Reco)
	//

	AliInfo("The 3rd element of the pt and cut array will correspond to the cut applied to the soft pion - please check that they are correct");
	if (fProngs>0){
		for (Int_t iP=0; iP<fProngs; iP++){
			fPtAccCut[iP]=ptAccCut[iP];
			fEtaAccCut[iP]=etaAccCut[iP];
		}
	}
	return;
}		

//___________________________________________________________

void AliCFVertexingHFCascade::SetAccCut()
{
	//
	// setting the pt and eta cut to be used in the Acceptance steps (MC+Reco)
	//
  
  AliAODMCParticle* mcPartDaughter = dynamic_cast<AliAODMCParticle*>(fmcArray->At(fLabelArray[2]));  // should be the soft pion...  
  if(!mcPartDaughter) return;
  Int_t mother =  mcPartDaughter->GetMother();
  AliAODMCParticle* mcMother = dynamic_cast<AliAODMCParticle*>(fmcArray->At(mother)); 
  if(!mcMother) return;

  if (TMath::Abs(mcPartDaughter->GetPdgCode())!= 211 || TMath::Abs(mcMother->GetPdgCode())!=413){
    AliFatal("Apparently the soft pion is not in the third position, causing a crash!!");
  }	         
  if (fProngs>0){
    for (Int_t iP=0; iP<fProngs-1; iP++){
      fPtAccCut[iP]=0.1;
      fEtaAccCut[iP]=0.9;
    }
    fPtAccCut[2]=0.06;  // soft pion
    fEtaAccCut[2]=0.9;  // soft pion
  }
  return;
}		

//_____________________________________________________________
Double_t AliCFVertexingHFCascade::GetEtaProng(Int_t iProng) const 
{
	//
	// getting eta of the prong - overload the mother class method
	//

 if (fRecoCandidate){
   
   AliAODRecoCascadeHF* dstarD0pi = (AliAODRecoCascadeHF*)fRecoCandidate;

   Double_t etaProng =-9999;
    if(iProng==0) etaProng =dstarD0pi->Get2Prong()->EtaProng(0);
    if(iProng==1) etaProng =dstarD0pi->Get2Prong()->EtaProng(1);
    if(iProng==2) etaProng =dstarD0pi->EtaProng(1);
    
    return etaProng;
    
  }
  return 999999;    
}
//_____________________________________________________________
Double_t AliCFVertexingHFCascade::GetPtProng(Int_t iProng) const 
{
	//
	// getting pt of the prong
	//
  
  if (fRecoCandidate){

    AliAODRecoCascadeHF* dstarD0pi = (AliAODRecoCascadeHF*)fRecoCandidate;
    Double_t ptProng= -9999;
    if(iProng==0) ptProng =dstarD0pi->Get2Prong()->PtProng(0);
    if(iProng==1) ptProng =dstarD0pi->Get2Prong()->PtProng(1);
    if(iProng==2) ptProng =dstarD0pi->PtProng(1);
    
    //	Double_t ptProng = fRecoCandidate->PtProng(iProng);  
    return ptProng;
    
  }
  return 999999;  
  
}
Commit	Line	Data
a9f238cf	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	//-----------------------------------------------------------------------
	17	// Class for HF corrections as a function of many variables and steps
	18	// For D* and other cascades
	19	//
	20	// Author : A.Grelli a.grelli@uu.nl UTECHT
	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	#include "AliAODRecoCascadeHF.h"
	31	#include "AliCFVertexingHFCascade.h"
	32	#include "AliCFContainer.h"
ec5288c3	33	#include "AliCFTaskVertexingHF.h"
a9f238cf	34
	35	ClassImp(AliCFVertexingHFCascade)
	36
	37
	38	//_________________________________________
	39	AliCFVertexingHFCascade::AliCFVertexingHFCascade(TClonesArray *mcArray, UShort_t originDselection):
	40	AliCFVertexingHF(mcArray, originDselection)
	41	{
	42	// standard constructor
	43
	44	SetNProngs(3);
2bf2e62b	45	fPtAccCut=new Float_t[fProngs];
	46	fEtaAccCut=new Float_t[fProngs];
	47	// element 0 in the cut arrays corresponds to the soft pion!!!!!!!! Careful when setting the values...
	48	fPtAccCut[0]=0.;
	49	fEtaAccCut[0]=0.;
	50	for(Int_t iP=1; iP<fProngs; iP++){
	51	fPtAccCut[iP]=0.1;
	52	fEtaAccCut[iP]=0.9;
	53	}
	54
a9f238cf	55	}
	56
	57
	58	//_____________________________________
	59	AliCFVertexingHFCascade& AliCFVertexingHFCascade::operator=(const AliCFVertexingHFCascade& c)
	60	{
	61	// operator =
	62
	63	if (this != &c) {
	64
	65	AliCFVertexingHF::operator=(c);
	66
	67	}
	68	return *this;
	69	}
	70
	71	//__________________________________________
	72	Bool_t AliCFVertexingHFCascade::SetRecoCandidateParam(AliAODRecoDecayHF *recoCand)
	73	{
	74	// set the AliAODRecoDecay candidate
	75
	76	Bool_t bSignAssoc = kFALSE;
	77
	78	fRecoCandidate = recoCand;
	79	AliAODRecoCascadeHF* dstarD0pi = (AliAODRecoCascadeHF*)recoCand;
	80
	81	if (!fRecoCandidate) {
	82	AliError("fRecoCandidate not found, problem in assignement\n");
	83	return bSignAssoc;
	84	}
	85
	86	if ( fRecoCandidate->GetPrimaryVtx()) AliDebug(3,"fReco Candidate has a pointer to PrimVtx\n");
	87	//if (recoCand->GetPrimaryVtx()) printf("Reco Cand has a pointer to PrimVtx\n");
	88
	89	//Int_t pdgCand = 413;
	90
	91	Int_t pdgDgDStartoD0pi[2]={421,211};
	92	Int_t pdgDgD0toKpi[2]={321,211};
	93
	94	Int_t nentries = fmcArray->GetEntriesFast();
	95
	96	AliDebug(3,Form("nentries = %d\n", nentries));
	97
	98	Int_t mcLabel = dstarD0pi->MatchToMC(413,421,pdgDgDStartoD0pi,pdgDgD0toKpi,fmcArray);
	99
fbec9fa9	100	if (mcLabel == -1) return bSignAssoc;
fbec9fa9	101
1f780958	102	if (fRecoCandidate->NumberOfFakeDaughters()>0){
	103	fFake = 0; // fake candidate
	104	if (fFakeSelection==1) return bSignAssoc;
	105	}
	106	if (fRecoCandidate->NumberOfFakeDaughters()==0){
	107	fFake = 2; // non-fake candidate
	108	if (fFakeSelection==2) return bSignAssoc;
	109	}
fbec9fa9	110
a9f238cf	111	SetMCLabel(mcLabel);
	112	fmcPartCandidate = dynamic_cast<AliAODMCParticle*>(fmcArray->At(fmcLabel));
	113
	114	if (!fmcPartCandidate){
	115	AliDebug(3,"No part candidate");
	116	return bSignAssoc;
	117	}
	118
	119	bSignAssoc = kTRUE;
	120	return bSignAssoc;
	121	}
	122
	123	//______________________________________________
	124	Bool_t AliCFVertexingHFCascade::GetGeneratedValuesFromMCParticle(Double_t* vectorMC)
	125	{
	126	//
	127	// collecting all the necessary info (pt, y, cosThetaStar, ptPi, ptKa, cT) from MC particle
	128	//
	129
	130	Bool_t bGenValues = kFALSE;
	131
	132
	133	Int_t daughter0ds = fmcPartCandidate->GetDaughter(0);
	134	Int_t daughter1ds = fmcPartCandidate->GetDaughter(1);
	135
	136	//the D0
	137	AliAODMCParticle* mcPartDaughterD0 = dynamic_cast<AliAODMCParticle*>(fmcArray->At(daughter0ds));
	138	AliAODMCParticle* mcPartDaughterPis = dynamic_cast<AliAODMCParticle*>(fmcArray->At(daughter1ds));
	139
a828d135	140	if(!mcPartDaughterD0 \|\| !mcPartDaughterPis) return kFALSE;
a828d135	141
a9f238cf	142	Double_t vtx1[3] = {0,0,0}; // primary vertex
	143	Double_t vtx2daughter0[3] = {0,0,0}; // secondary vertex from daughter 0
	144	Double_t vtx2daughter1[3] = {0,0,0}; // secondary vertex from daughter 1
	145	fmcPartCandidate->XvYvZv(vtx1); // cm
	146
	147	//D0 daughters
	148	Int_t daughter0 = mcPartDaughterD0->GetDaughter(0);
	149	Int_t daughter1 = mcPartDaughterD0->GetDaughter(1);
	150
	151	AliAODMCParticle* mcPartDaughter0 = dynamic_cast<AliAODMCParticle*>(fmcArray->At(daughter0)); //D0
	152	AliAODMCParticle* mcPartDaughter1 = dynamic_cast<AliAODMCParticle*>(fmcArray->At(daughter1)); //pis
	153
a828d135	154	if(!mcPartDaughter0 \|\| !mcPartDaughter1) return kFALSE;
a828d135	155
a9f238cf	156	// getting vertex from daughters
	157	mcPartDaughter0->XvYvZv(vtx2daughter0); // cm
	158	mcPartDaughter1->XvYvZv(vtx2daughter1); //cm
	159	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) {
	160	AliError("Daughters have different secondary vertex, skipping the track");
	161	return bGenValues;
	162	}
	163
	164	Int_t nprongs = 2;
	165	Short_t charge = 0;
	166	// always instantiate the AliAODRecoDecay with the positive daughter first, the negative second
	167	AliAODMCParticle* positiveDaugh = mcPartDaughter0;
	168	AliAODMCParticle* negativeDaugh = mcPartDaughter1;
	169	if (mcPartDaughter0->GetPdgCode()<0 && mcPartDaughter1->GetPdgCode()>0){
	170	// inverting in case the positive daughter is the second one
	171	positiveDaugh = mcPartDaughter1;
	172	negativeDaugh = mcPartDaughter0;
	173	}
	174	// getting the momentum from the daughters
	175	Double_t px[2] = {positiveDaugh->Px(), negativeDaugh->Px()};
	176	Double_t py[2] = {positiveDaugh->Py(), negativeDaugh->Py()};
	177	Double_t pz[2] = {positiveDaugh->Pz(), negativeDaugh->Pz()};
	178
	179	Double_t d0[2] = {0.,0.};
	180
	181	AliAODRecoDecayHF* decay = new AliAODRecoDecayHF(vtx1,vtx2daughter0,nprongs,charge,px,py,pz,d0);
	182
	183	Double_t cosThetaStar = 0.;
	184	Double_t cosThetaStarD0 = 0.;
	185	Double_t cosThetaStarD0bar = 0.;
	186	cosThetaStarD0 = decay->CosThetaStar(1,421,211,321);
	187	cosThetaStarD0bar = decay->CosThetaStar(0,421,321,211);
	188	if (mcPartDaughterD0->GetPdgCode() == 421){ // D0
	189	AliDebug(3, Form("D0, with pdgprong0 = %d, pdgprong1 = %d",mcPartDaughter0->GetPdgCode(),mcPartDaughter1->GetPdgCode()));
	190	cosThetaStar = cosThetaStarD0;
	191	}
	192	else if (mcPartDaughterD0->GetPdgCode() == -421){ // D0bar{
	193	AliDebug(3, Form("D0bar, with pdgprong0 = %d, pdgprong1 = %d",mcPartDaughter0->GetPdgCode(),mcPartDaughter1->GetPdgCode()));
	194	cosThetaStar = cosThetaStarD0bar;
	195	}
	196	else{
	197	AliWarning("There are problems!! particle was expected to be either a D0 or a D0bar, check...");
b477193a	198	delete decay;
a9f238cf	199	return vectorMC;
	200	}
	201	if (cosThetaStar < -1 \|\| cosThetaStar > 1) {
	202	AliWarning("Invalid value for cosine Theta star, returning");
b477193a	203	delete decay;
a9f238cf	204	return bGenValues;
	205	}
	206
	207	Double_t vectorD0[2] ={0.,0.};
	208
	209	// evaluate the correct cascade
	210	if (!EvaluateIfD0toKpi(mcPartDaughterD0,vectorD0)) {
	211	AliDebug(2, "Error! the D0 MC doesn't have correct daughters!!");
b477193a	212	delete decay;
a9f238cf	213	return bGenValues;
	214	}
	215
	216	//ct
	217	Double_t cT = decay->Ct(421);
	218	// get the pT of the daughters
	219	Double_t pTD0 = 0.;
	220	Double_t pTpi = 0.;
	221
	222	if (TMath::Abs(fmcPartCandidate->GetPdgCode()) == 413) {
	223	pTD0 = mcPartDaughterD0->Pt();
	224	pTpi = mcPartDaughterPis->Pt();
	225	}
	226
	227
ec5288c3	228	switch (fConfiguration){
	229	case AliCFTaskVertexingHF::kSnail:
	230	vectorMC[0] = fmcPartCandidate->Pt();
	231	vectorMC[1] = fmcPartCandidate->Y() ;
	232	vectorMC[2] = cosThetaStar ;
	233	vectorMC[3] = vectorD0[0];
	234	vectorMC[4] = vectorD0[1];
	235	vectorMC[5] = cT*1.E4 ; // in micron
	236	vectorMC[6] = 0.; // dummy value, meaningless in MC
	237	vectorMC[7] = -100000.; // dummy value, meaningless in MC, in micron^2
	238	vectorMC[8] = 1.01; // dummy value, meaningless in MC
	239	vectorMC[9] = fmcPartCandidate->Phi();
	240	vectorMC[10] = fzMCVertex; // z of reconstructed of primary vertex
	241	vectorMC[11] = fCentValue; // reconstructed centrality
	242	vectorMC[12] = 1.; // always filling with 1 at MC level
	243	vectorMC[13] = 1.01; // dummy value for cosPointingXY multiplicity
	244	vectorMC[14] = 0.; // dummy value for NormalizedDecayLengthXY multiplicity
	245	vectorMC[15] = fMultiplicity; // reconstructed multiplicity
	246	break;
	247	case AliCFTaskVertexingHF::kCheetah:
	248	vectorMC[0] = fmcPartCandidate->Pt();
	249	vectorMC[1] = fmcPartCandidate->Y() ;
	250	vectorMC[2] = cT*1.E4; // in micron
	251	vectorMC[3] = fmcPartCandidate->Phi();
	252	vectorMC[4] = fzMCVertex;
	253	vectorMC[5] = fCentValue; // dummy value for dca, meaningless in MC
	254	vectorMC[6] = 1. ; // fake: always filling with 1 at MC level
	255	vectorMC[7] = fMultiplicity; // dummy value for d0pi, meaningless in MC, in micron
	256	break;
	257	}
a9f238cf	258
	259	delete decay;
	260	bGenValues = kTRUE;
	261	return bGenValues;
	262	}
	263	//____________________________________________
	264	Bool_t AliCFVertexingHFCascade::GetRecoValuesFromCandidate(Double_t *vectorReco) const
	265	{
	266	// read the variables for the container
	267
	268	Bool_t bFillRecoValues=kFALSE;
	269
	270	//Get the D* and the D0 from D*
	271	AliAODRecoCascadeHF* dstarD0pi = (AliAODRecoCascadeHF*)fRecoCandidate;
	272	AliAODRecoDecayHF2Prong* d0toKpi = (AliAODRecoDecayHF2Prong*)dstarD0pi->Get2Prong();
	273
	274
c9b41591	275	//if (dstarD0pi->GetPrimaryVtx())printf("dstarD0pi has primary vtx\n");
c9b41591	276	//if (fRecoCandidate->GetPrimaryVtx())printf("fRecoCandidateDstar has primary vtx\n");
a9f238cf	277
	278	Double_t pt = dstarD0pi->Pt();
	279	Double_t rapidity = dstarD0pi->YDstar();
	280	Double_t invMass=0.;
	281	Double_t cosThetaStar = 9999.;
	282	Double_t pTpi = 0.;
	283	Double_t pTK = 0.;
	284	Double_t dca = d0toKpi->GetDCA();
	285	Double_t d0pi = 0.;
	286	Double_t d0K = 0.;
	287	Double_t d0xd0 = d0toKpi->Prodd0d0();
	288	Double_t cosPointingAngle = d0toKpi->CosPointingAngle();
	289	Double_t phi = dstarD0pi->Phi();
ec5288c3	290	Double_t cosPointingAngleXY = d0toKpi->CosPointingAngleXY();
ec5288c3	291	Double_t normDecayLengthXY = d0toKpi->NormalizedDecayLengthXY();
a9f238cf	292
	293	Int_t pdgCode = fmcPartCandidate->GetPdgCode();
	294
	295	if (pdgCode > 0){
	296	cosThetaStar = d0toKpi->CosThetaStarD0();
	297	pTpi = d0toKpi->PtProng(0);
	298	pTK = d0toKpi->PtProng(1);
	299	d0pi = d0toKpi->Getd0Prong(0);
	300	d0K = d0toKpi->Getd0Prong(1);
	301	invMass=d0toKpi->InvMassD0();
	302	}
	303	else {
	304	cosThetaStar = d0toKpi->CosThetaStarD0bar();
	305	pTpi = d0toKpi->PtProng(1);
	306	pTK = d0toKpi->PtProng(0);
	307	d0pi = d0toKpi->Getd0Prong(1);
	308	d0K = d0toKpi->Getd0Prong(0);
	309	invMass= d0toKpi->InvMassD0bar();
	310	}
	311
	312	Double_t cT = d0toKpi->CtD0();
	313
ec5288c3	314	switch (fConfiguration){
	315	case AliCFTaskVertexingHF::kSnail:
	316	vectorReco[0] = pt;
	317	vectorReco[1] = rapidity;
	318	vectorReco[2] = cosThetaStar;
	319	vectorReco[3] = pTpi;
	320	vectorReco[4] = pTK;
	321	vectorReco[5] = cT*1.E4; // in micron
	322	vectorReco[6] = dca*1.E4; // in micron
	323	vectorReco[7] = d0xd0*1.E8; // in micron^2
	324	vectorReco[8] = cosPointingAngle; // in micron
	325	vectorReco[9] = phi;
	326	vectorReco[10] = fzPrimVertex; // z of reconstructed of primary vertex
	327	vectorReco[11] = fCentValue;
	328	vectorReco[12] = fFake; // whether the reconstructed candidate was a fake (fFake = 0) or not (fFake = 2)
	329	vectorReco[13] = cosPointingAngleXY;
	330	vectorReco[14] = normDecayLengthXY; // in cm
	331	vectorReco[15] = fMultiplicity; // reconstructed multiplicity
	332	break;
	333	case AliCFTaskVertexingHF::kCheetah:
	334	vectorReco[0] = pt;
	335	vectorReco[1] = rapidity ;
	336	vectorReco[2] = cT*1.E4; // in micron
	337	vectorReco[3] = phi;
	338	vectorReco[4] = fzPrimVertex;
	339	vectorReco[5] = fCentValue;
	340	vectorReco[6] = fFake ;
	341	vectorReco[7] = fMultiplicity;
	342	break;
	343	}
a9f238cf	344
ec5288c3	345	bFillRecoValues = kTRUE;
ec5288c3	346
a9f238cf	347	return bFillRecoValues;
	348	}
	349
	350
	351	//_____________________________________________________________
	352	Bool_t AliCFVertexingHFCascade::CheckMCChannelDecay() const
	353	{
	354	// check the required decay channel
	355
	356	Bool_t checkCD = kFALSE;
	357
	358
	359	Int_t daughter0 = fmcPartCandidate->GetDaughter(0);
	360	Int_t daughter1 = fmcPartCandidate->GetDaughter(1);
	361	AliAODMCParticle* mcPartDaughter0 = dynamic_cast<AliAODMCParticle*>(fmcArray->At(daughter0));
	362	AliAODMCParticle* mcPartDaughter1 = dynamic_cast<AliAODMCParticle*>(fmcArray->At(daughter1));
	363
	364	if (!mcPartDaughter0 \|\| !mcPartDaughter1) {
	365	AliDebug (2,"Problems in the MC Daughters\n");
	366	return checkCD;
	367	}
	368
	369	if (!(TMath::Abs(mcPartDaughter0->GetPdgCode())==421 &&
	370	TMath::Abs(mcPartDaughter1->GetPdgCode())==211) &&
	371	!(TMath::Abs(mcPartDaughter0->GetPdgCode())==211 &&
	372	TMath::Abs(mcPartDaughter1->GetPdgCode())==421)) {
	373	AliDebug(2, "The D0 MC doesn't come from a Kpi decay, skipping!!");
	374	return checkCD;
	375	}
	376
	377	Double_t vectorD0[2] ={0.,0.};
	378
	379	// D* is a cascade ...evaluate the correct cascade
	380	if (!EvaluateIfD0toKpi(mcPartDaughter0,vectorD0)) {
	381	AliDebug(2, "Error! the D0 MC doesn't have correct daughters!!");
	382	return checkCD;
	383	}
	384
	385	checkCD = kTRUE;
	386	return checkCD;
	387
	388	}
	389
	390	//__________________________________________
	391	Bool_t AliCFVertexingHFCascade::EvaluateIfD0toKpi(AliAODMCParticle* neutralDaugh, Double_t* vectorD0)const
	392	{
	393	//
	394	// chack wether D0 is decaing into kpi
	395	//
	396
	397	Bool_t isHadronic = kFALSE;
	398
	399	Int_t daughterD00 = neutralDaugh->GetDaughter(0);
	400	Int_t daughterD01 = neutralDaugh->GetDaughter(1);
	401
	402	AliDebug(2, Form("daughter0 = %d and daughter1 = %d",daughterD00,daughterD01));
	403	if (daughterD00 == 0 \|\| daughterD01 == 0) {
	404	AliDebug(2, "Error! the D0 MC doesn't have correct daughters!!");
	405	return isHadronic;
	406	}
	407
	408	if (TMath::Abs(daughterD01 - daughterD00) != 1) { // should be everytime true - see PDGdatabooklet
	409	AliDebug(2, "The D0 MC doesn't come from a 2-prong decay, skipping!!");
	410	return isHadronic;
411	}
412
413	AliAODMCParticle* mcPartDaughterD00 = dynamic_cast<AliAODMCParticle*>(fmcArray->At(daughterD00));
414	AliAODMCParticle* mcPartDaughterD01 = dynamic_cast<AliAODMCParticle*>(fmcArray->At(daughterD01));
415	if (!mcPartDaughterD00 \|\| !mcPartDaughterD01) {
416	AliWarning("D0 MC analysis: At least one Daughter Particle not found in tree, skipping");
417	return isHadronic;
418	}
419
420	if (!(TMath::Abs(mcPartDaughterD00->GetPdgCode())==321 &&
421	TMath::Abs(mcPartDaughterD01->GetPdgCode())==211) &&
422	!(TMath::Abs(mcPartDaughterD00->GetPdgCode())==211 &&
423	TMath::Abs(mcPartDaughterD01->GetPdgCode())==321)) {
424	AliDebug(2, "The D0 MC doesn't come from a Kpi decay, skipping!!");
425	return isHadronic;
426	}
427
88bc191b	428	Double_t sumPxDau=mcPartDaughterD00->Px()+mcPartDaughterD01->Px();
	429	Double_t sumPyDau=mcPartDaughterD00->Py()+mcPartDaughterD01->Py();
	430	Double_t sumPzDau=mcPartDaughterD00->Pz()+mcPartDaughterD01->Pz();
	431	Double_t pxMother=neutralDaugh->Px();
	432	Double_t pyMother=neutralDaugh->Py();
	433	Double_t pzMother=neutralDaugh->Pz();
	434	if(TMath::Abs(pxMother-sumPxDau)/(TMath::Abs(pxMother)+1.e-13)>0.00001 \|\|
	435	TMath::Abs(pyMother-sumPyDau)/(TMath::Abs(pyMother)+1.e-13)>0.00001 \|\|
	436	TMath::Abs(pzMother-sumPzDau)/(TMath::Abs(pzMother)+1.e-13)>0.00001){
	437	AliDebug(2, "Momentum conservation violated, skipping!!");
	438	return isHadronic;
	439	}
	440
a9f238cf	441	Double_t pTD0pi = 0;
	442	Double_t pTD0K = 0;
	443
	444
	445	if (TMath::Abs(mcPartDaughterD00->GetPdgCode()) == 211) {
	446	pTD0pi = mcPartDaughterD00->Pt();
	447	pTD0K = mcPartDaughterD01->Pt();
	448	}
	449	else {
	450	pTD0pi = mcPartDaughterD01->Pt();
	451	pTD0K = mcPartDaughterD00->Pt();
	452	}
	453
	454	isHadronic = kTRUE;
	455
	456	vectorD0[0] = pTD0pi;
	457	vectorD0[1] = pTD0K;
	458
	459	return isHadronic;
	460
	461	}
	462
2bf2e62b	463	//___________________________________________________________
	464
	465	void AliCFVertexingHFCascade::SetPtAccCut(Float_t* ptAccCut)
	466	{
	467	//
	468	// setting the pt cut to be used in the Acceptance steps (MC+Reco)
	469	//
	470
	471	AliInfo("The 3rd element of the pt cut array will correspond to the cut applied to the soft pion - please check that it is correct");
	472	if (fProngs>0){
	473	for (Int_t iP=0; iP<fProngs; iP++){
	474	fPtAccCut[iP]=ptAccCut[iP];
	475	}
	476	}
	477	return;
	478	}
	479
	480
	481
	482	//___________________________________________________________
	483
	484	void AliCFVertexingHFCascade::SetEtaAccCut(Float_t* etaAccCut)
	485	{
	486	//
	487	// setting the eta cut to be used in the Acceptance steps (MC+Reco)
	488	//
	489
	490	AliInfo("The 3rd element of the eta cut array will correspond to the cut applied to the soft pion - please check that it is correct");
	491	if (fProngs>0){
	492	for (Int_t iP=0; iP<fProngs; iP++){
	493	fEtaAccCut[iP]=etaAccCut[iP];
	494	}
	495	}
	496	return;
	497	}
	498	//___________________________________________________________
	499
	500	void AliCFVertexingHFCascade::SetAccCut(Float_t* ptAccCut, Float_t* etaAccCut)
	501	{
	502	//
	503	// setting the pt and eta cut to be used in the Acceptance steps (MC+Reco)
	504	//
	505
	506	AliInfo("The 3rd element of the pt and cut array will correspond to the cut applied to the soft pion - please check that they are correct");
	507	if (fProngs>0){
	508	for (Int_t iP=0; iP<fProngs; iP++){
	509	fPtAccCut[iP]=ptAccCut[iP];
	510	fEtaAccCut[iP]=etaAccCut[iP];
	511	}
	512	}
	513	return;
	514	}
	515
	516	//___________________________________________________________
	517
	518	void AliCFVertexingHFCascade::SetAccCut()
	519	{
	520	//
	521	// setting the pt and eta cut to be used in the Acceptance steps (MC+Reco)
	522	//
367e9aa3	523
	524	AliAODMCParticle* mcPartDaughter = dynamic_cast<AliAODMCParticle*>(fmcArray->At(fLabelArray[2])); // should be the soft pion...
	525	if(!mcPartDaughter) return;
	526	Int_t mother = mcPartDaughter->GetMother();
	527	AliAODMCParticle* mcMother = dynamic_cast<AliAODMCParticle*>(fmcArray->At(mother));
	528	if(!mcMother) return;
	529
	530	if (TMath::Abs(mcPartDaughter->GetPdgCode())!= 211 \|\| TMath::Abs(mcMother->GetPdgCode())!=413){
	531	AliFatal("Apparently the soft pion is not in the third position, causing a crash!!");
	532	}
	533	if (fProngs>0){
	534	for (Int_t iP=0; iP<fProngs-1; iP++){
	535	fPtAccCut[iP]=0.1;
	536	fEtaAccCut[iP]=0.9;
	537	}
	538	fPtAccCut[2]=0.06; // soft pion
	539	fEtaAccCut[2]=0.9; // soft pion
	540	}
	541	return;
2bf2e62b	542	}
a9f238cf	543
367e9aa3	544	//_____________________________________________________________
	545	Double_t AliCFVertexingHFCascade::GetEtaProng(Int_t iProng) const
	546	{
	547	//
	548	// getting eta of the prong - overload the mother class method
	549	//
a9f238cf	550
367e9aa3	551	if (fRecoCandidate){
	552
	553	AliAODRecoCascadeHF* dstarD0pi = (AliAODRecoCascadeHF*)fRecoCandidate;
	554
	555	Double_t etaProng =-9999;
	556	if(iProng==0) etaProng =dstarD0pi->Get2Prong()->EtaProng(0);
	557	if(iProng==1) etaProng =dstarD0pi->Get2Prong()->EtaProng(1);
	558	if(iProng==2) etaProng =dstarD0pi->EtaProng(1);
	559
	560	return etaProng;
	561
	562	}
	563	return 999999;
	564	}
	565	//_____________________________________________________________
	566	Double_t AliCFVertexingHFCascade::GetPtProng(Int_t iProng) const
	567	{
	568	//
	569	// getting pt of the prong
	570	//
	571
	572	if (fRecoCandidate){
	573
	574	AliAODRecoCascadeHF* dstarD0pi = (AliAODRecoCascadeHF*)fRecoCandidate;
	575	Double_t ptProng= -9999;
	576	if(iProng==0) ptProng =dstarD0pi->Get2Prong()->PtProng(0);
	577	if(iProng==1) ptProng =dstarD0pi->Get2Prong()->PtProng(1);
	578	if(iProng==2) ptProng =dstarD0pi->PtProng(1);
	579
	580	// Double_t ptProng = fRecoCandidate->PtProng(iProng);
	581	return ptProng;
	582
	583	}
	584	return 999999;
	585
	586	}