[u/mrichter/AliRoot.git] / PWGHF / correlationHF / AliHFAssociatedTrackCuts.cxx

/**************************************************************************
 * Copyright(c) 1998-2010, 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.                  *
 **************************************************************************/

/* $Id$ */

////////////////////////////////////////////////////////////////////////
//
// Base class for cuts on Associated tracks for HF Correlation analysis
//
// Author: S.Bjelogrlic (Utrecht) sandro.bjelogrlic@cern.ch
////////////////////////////////////////////////////////////////////////
#include <Riostream.h>
#include "AliHFAssociatedTrackCuts.h"
#include "AliAODPidHF.h"
#include "AliESDtrackCuts.h"
#include "AliESDtrack.h"
#include "AliAODv0.h"
#include "AliAODVertex.h"
#include "AliAODMCParticle.h"
#include "TString.h"

using std::cout;
using std::endl;

ClassImp(AliHFAssociatedTrackCuts)

//--------------------------------------------------------------------------
AliHFAssociatedTrackCuts::AliHFAssociatedTrackCuts():
AliAnalysisCuts(),
fESDTrackCuts(0),
fPidObj(0),

fPoolMaxNEvents(0), 
fPoolMinNTracks(0), 
fMinEventsToMix(0),
fNzVtxBins(0), 
fNzVtxBinsDim(0), 
fZvtxBins(0), 
fNCentBins(0), 
fNCentBinsDim(0), 
fCentBins(0),

fNTrackCuts(0),
fAODTrackCuts(0),
fTrackCutsNames(0),
fNvZeroCuts(0),
fAODvZeroCuts(0),
fvZeroCutsNames(0),
fBit(0),
fCharge(0) 

{
	//
	//default constructor
	//
	//
	//default constructor
	//
	
}

//--------------------------------------------------------------------------
AliHFAssociatedTrackCuts::AliHFAssociatedTrackCuts(const char* name, const char* title):
AliAnalysisCuts(name,title),
fESDTrackCuts(0),
fPidObj(0),

fPoolMaxNEvents(0), 
fPoolMinNTracks(0), 
fMinEventsToMix(0),
fNzVtxBins(0), 
fNzVtxBinsDim(0), 
fZvtxBins(0), 
fNCentBins(0), 
fNCentBinsDim(0), 
fCentBins(0),

fNTrackCuts(0),
fAODTrackCuts(0),
fTrackCutsNames(0),
fNvZeroCuts(0),
fAODvZeroCuts(0),
fvZeroCutsNames(0),
fBit(0),
fCharge(0)

{
	//
	//default constructor
	//
	
}
//--------------------------------------------------------------------------
AliHFAssociatedTrackCuts::AliHFAssociatedTrackCuts(const AliHFAssociatedTrackCuts &source) :
AliAnalysisCuts(source),
fESDTrackCuts(source.fESDTrackCuts),
fPidObj(source.fPidObj),

fPoolMaxNEvents(source.fPoolMaxNEvents), 
fPoolMinNTracks(source.fPoolMinNTracks), 
fMinEventsToMix(source.fMinEventsToMix),
fNzVtxBins(source.fNzVtxBins), 
fNzVtxBinsDim(source.fNzVtxBinsDim), 
fZvtxBins(source.fZvtxBins), 
fNCentBins(source.fNCentBins), 
fNCentBinsDim(source.fNCentBinsDim), 
fCentBins(source.fCentBins),

fNTrackCuts(source.fNTrackCuts),
fAODTrackCuts(source.fAODTrackCuts),
fTrackCutsNames(source.fTrackCutsNames),
fNvZeroCuts(source.fNvZeroCuts),
fAODvZeroCuts(source.fAODvZeroCuts),
fvZeroCutsNames(source.fvZeroCutsNames),
fBit(source.fBit),
fCharge(source.fCharge)
{
	//
	// copy constructor
	//
	

        AliInfo("AliHFAssociatedTrackCuts::Copy constructor ");
	if(source.fESDTrackCuts) AddTrackCuts(source.fESDTrackCuts);
	if(source.fAODTrackCuts) SetAODTrackCuts(source.fAODTrackCuts);
	if(source.fAODvZeroCuts) SetAODvZeroCuts(source.fAODvZeroCuts);
	if(source.fPidObj) SetPidHF(source.fPidObj);
}
//--------------------------------------------------------------------------
AliHFAssociatedTrackCuts &AliHFAssociatedTrackCuts::operator=(const AliHFAssociatedTrackCuts &source)
{
	//
	// assignment operator
	//	
	if(&source == this) return *this;
	
	AliAnalysisCuts::operator=(source);
	fESDTrackCuts=source.fESDTrackCuts;
	fPidObj=source.fPidObj;
	fNTrackCuts=source.fNTrackCuts;
	fAODTrackCuts=source.fAODTrackCuts;
	fTrackCutsNames=source.fTrackCutsNames;
	fNvZeroCuts=source.fNvZeroCuts;
	fAODvZeroCuts=source.fAODvZeroCuts;
	fvZeroCutsNames=source.fvZeroCutsNames;
	fBit=source.fBit;
	fCharge=source.fCharge;
	
	return *this;

}


//--------------------------------------------------------------------------
AliHFAssociatedTrackCuts::~AliHFAssociatedTrackCuts()
{
	if(fESDTrackCuts) {delete fESDTrackCuts; fESDTrackCuts = 0;}
	if(fPidObj) {delete fPidObj; fPidObj = 0;}
	if(fZvtxBins) {delete[] fZvtxBins; fZvtxBins=0;} 
	if(fCentBins) {delete[] fCentBins; fCentBins=0;}
	if(fAODTrackCuts) {delete[] fAODTrackCuts; fAODTrackCuts=0;}
	if(fTrackCutsNames) {delete[] fTrackCutsNames; fTrackCutsNames=0;}
	if(fAODvZeroCuts){delete[] fAODvZeroCuts; fAODvZeroCuts=0;}
	if(fvZeroCutsNames) {delete[] fvZeroCutsNames; fvZeroCutsNames=0;}
	
}
//--------------------------------------------------------------------------
Bool_t AliHFAssociatedTrackCuts::IsInAcceptance()
{
	printf("Careful: method AliHFAssociatedTrackCuts::IsInAcceptance is not implemented yet \n");
	return kFALSE;
}
//--------------------------------------------------------------------------
Bool_t AliHFAssociatedTrackCuts::IsHadronSelected(AliAODTrack * track)
{
	AliESDtrack esdtrack(track);
	if(!fESDTrackCuts->IsSelected(&esdtrack)) return kFALSE;
	
	if(fBit && !track->TestFilterBit(fBit)) return kFALSE; // check the filter bit
 
	return kTRUE;
	
}

//--------------------------------------------------------------------------
Bool_t AliHFAssociatedTrackCuts::CheckHadronKinematic(Double_t pt, Double_t d0) 
{
	
	
	if(pt < fAODTrackCuts[0]) return kFALSE;
	if(pt > fAODTrackCuts[1]) return kFALSE;
	if(d0 < fAODTrackCuts[2]) return kFALSE;
	if(d0 > fAODTrackCuts[3]) return kFALSE;
	
	return kTRUE;

	
}
//--------------------------------------------------------------------------

Bool_t AliHFAssociatedTrackCuts::Charge(Short_t charge, AliAODTrack* track) 
{// charge is the charge to compare to (for example, a daughter of a D meson)
	
	if(!fCharge) return kTRUE; // if fCharge is set to 0 (no selection on the charge), returns always true
	if(track->Charge()!= fCharge*charge) return kFALSE;
	return kTRUE;
}

//--------------------------------------------------------------------------
Bool_t AliHFAssociatedTrackCuts::CheckKaonCompatibility(AliAODTrack * track, Bool_t useMc, TClonesArray* mcArray, Int_t method)
{
	Bool_t isKaon = kFALSE;
	
	if(useMc) { // on MC
		Int_t hadLabel = track->GetLabel();
		if(hadLabel < 0) return kFALSE;
		AliAODMCParticle* hadron = dynamic_cast<AliAODMCParticle*>(mcArray->At(hadLabel));
		Int_t pdg = TMath::Abs(hadron->GetPdgCode()); 
		if (pdg == 321) isKaon = kTRUE;
	}
	
	if(!useMc) { // on DATA
          switch(method) {
	  case(1): {
		Bool_t isKTPC=kFALSE;
		Bool_t isPiTPC=kFALSE;
		Bool_t isPTPC=kFALSE;
		Bool_t isKTOF=kFALSE;
		Bool_t isPiTOF=kFALSE;
		Bool_t isPTOF=kFALSE;
		
		Bool_t KaonHyp = kFALSE;
		Bool_t PionHyp = kFALSE;
		Bool_t ProtonHyp = kFALSE;
		
		if(fPidObj->CheckStatus(track,"TOF")) {
			isKTOF=fPidObj->IsKaonRaw(track,"TOF");
			isPiTOF=fPidObj->IsPionRaw(track,"TOF");
			isPTOF=fPidObj->IsProtonRaw(track,"TOF");
		}
		if(fPidObj->CheckStatus(track,"TPC")){
			isKTPC=fPidObj->IsKaonRaw(track,"TPC");
			isPiTPC=fPidObj->IsPionRaw(track,"TPC");
			isPTPC=fPidObj->IsProtonRaw(track,"TPC");		
		}
		
		if (isKTOF && isKTPC) KaonHyp = kTRUE;
		if (isPiTOF && isPiTPC) PionHyp = kTRUE;
		if (isPTOF && isPTPC) ProtonHyp = kTRUE;
		
		if(KaonHyp && !PionHyp && !ProtonHyp) isKaon = kTRUE; 
	        break;
	      }
      case(2): {
		if(fPidObj->MakeRawPid(track,3)>=1) isKaon = kTRUE;
		break;
	      }
      }
	}
	
	return isKaon;
	
}
//--------------------------------------------------------------------------
Bool_t AliHFAssociatedTrackCuts::IsKZeroSelected(AliAODv0 *vzero, AliAODVertex *vtx1)
{
	
	if(vzero->DcaV0Daughters()>fAODvZeroCuts[0]) return kFALSE;
	if(vzero->Chi2V0()>fAODvZeroCuts[1]) return kFALSE;
	if(vzero->DecayLength(vtx1) < fAODvZeroCuts[2]) return kFALSE;
	if(vzero->DecayLength(vtx1) > fAODvZeroCuts[3]) return kFALSE;
	if(vzero->OpenAngleV0() > fAODvZeroCuts[4]) return kFALSE;
	if(vzero->Pt() < fAODvZeroCuts[5]) return kFALSE;
	if(TMath::Abs(vzero->Eta()) > fAODvZeroCuts[6]) return kFALSE;

	
	return kTRUE;
}
//--------------------------------------------------------------------------
Int_t AliHFAssociatedTrackCuts::IsMCpartFromHF(Int_t label, TClonesArray*mcArray){
  // Check origin in MC

  AliAODMCParticle* mcParticle;
  Int_t pdgCode = -1;
  if (label<0) return 0;
  Bool_t isCharmy = kFALSE;
  Bool_t isBeauty = kFALSE;
  Bool_t isD = kFALSE;
  Bool_t isB = kFALSE;

  while(pdgCode!=2212){ // loops back to the collision to check the particle source

    mcParticle = dynamic_cast<AliAODMCParticle*>(mcArray->At(label));
    if(!mcParticle) {AliError("NO MC PARTICLE"); break;}
    pdgCode =  TMath::Abs(mcParticle->GetPdgCode());

    label = mcParticle->GetMother();


    if((pdgCode>=400 && pdgCode <500) || (pdgCode>=4000 && pdgCode<5000 )) isD = kTRUE;
    if((pdgCode>=500 && pdgCode <600) || (pdgCode>=5000 && pdgCode<6000 )) {isD = kFALSE; isB = kTRUE;}


    if(pdgCode == 4) isCharmy = kTRUE;
    if(pdgCode == 5) {isBeauty = kTRUE; isCharmy = kFALSE;}
    if(label<0) break;

  }
  if (isCharmy && isD) return 44; // the first 4(5) indicates the charm/beauty quark, the second the charmy meson
  if (isBeauty && isD) return 54;
  if (isBeauty && isB) return 55;
  return 0;
}

//--------------------------------------------------------------------------
Bool_t AliHFAssociatedTrackCuts::InvMassDstarRejection(AliAODRecoDecayHF2Prong* d, AliAODTrack *track, Int_t hypD0) const {
	//
	// Calculates invmass of track+D0 and rejects if compatible with D*
	// (to remove pions from D*)
	// 
	Double_t nsigma = 3.;
	
	Double_t mD0, mD0bar;
	d->InvMassD0(mD0,mD0bar);
	
	Double_t invmassDstar1 = 0, invmassDstar2 = 0; 
	Double_t e1Pi = d->EProng(0,211), e2K = d->EProng(1,321); //hyp 1 (pi,K) - D0
	Double_t e1K = d->EProng(0,321), e2Pi = d->EProng(1,211); //hyp 2 (K,pi) - D0bar
	Double_t psum2 = (d->Px()+track->Px())*(d->Px()+track->Px())
	+(d->Py()+track->Py())*(d->Py()+track->Py())
	+(d->Pz()+track->Pz())*(d->Pz()+track->Pz());
	
	switch(hypD0) {
		case 1:
			invmassDstar1 = TMath::Sqrt(pow(e1Pi+e2K+track->E(0.1396),2.)-psum2);
			if ((TMath::Abs(invmassDstar1-mD0)-0.14543) < nsigma*800.*pow(10.,-6.)) return kFALSE;
			break;
		case 2:
			invmassDstar2 = TMath::Sqrt(pow(e2Pi+e1K+track->E(0.1396),2.)-psum2);
			if ((TMath::Abs(invmassDstar2-mD0bar)-0.14543) < nsigma*800.*pow(10.,-6.)) return kFALSE;
			break;
		case 3:
			invmassDstar1 = TMath::Sqrt(pow(e1Pi+e2K+track->E(0.1396),2.)-psum2);
			invmassDstar2 = TMath::Sqrt(pow(e2Pi+e1K+track->E(0.1396),2.)-psum2);
			if ((TMath::Abs(invmassDstar1-mD0)-0.14543) < nsigma*800.*pow(10.,-6.)) return kFALSE;
			if ((TMath::Abs(invmassDstar2-mD0bar)-0.14543) < nsigma*800.*pow(10.,-6.)) return kFALSE;
			break;
	}
	
	return kTRUE;
}
//________________________________________________________
void AliHFAssociatedTrackCuts::SetAODTrackCuts(Float_t *cutsarray)
{
	cout << "Check 1" << endl;
	if(!fAODTrackCuts) fAODTrackCuts = new Float_t[fNTrackCuts];
	cout << "Check 2" << endl;
	for(Int_t i =0; i<fNTrackCuts; i++){
		cout << "Check 2." << i << endl;
		fAODTrackCuts[i] = cutsarray[i];
	}
	cout << "Check 3" << endl;
	SetTrackCutsNames();
	cout << "Check 4" << endl;
	return;
}
//________________________________________________________
void AliHFAssociatedTrackCuts::SetTrackCutsNames(/*TString *namearray*/){
	
	fTrackCutsNames = new TString[4];
	fTrackCutsNames[0]= "associated track:: pt min [GeV/c]................: ";
	fTrackCutsNames[1]= "associated track:: pt max [GeV/c]................: ";
	fTrackCutsNames[2]= "associated track:: d0 min [cm]...................: ";
	fTrackCutsNames[3]= "associated track:: d0 max [cm]...................: ";
	

	return;
}
//--------------------------------------------------------------------------
void AliHFAssociatedTrackCuts::SetAODvZeroCuts(Float_t *cutsarray)
{
	

	if(!fAODvZeroCuts) fAODvZeroCuts = new Float_t[fNvZeroCuts];
	for(Int_t i =0; i<fNvZeroCuts; i++){
		fAODvZeroCuts[i] = cutsarray[i];
	}
	SetvZeroCutsNames();
	return;
}
//--------------------------------------------------------------------------
void AliHFAssociatedTrackCuts::SetvZeroCutsNames(/*TString *namearray*/){
	
	fvZeroCutsNames = new TString[7];
	fvZeroCutsNames[0] = "vZero:: max DCA between two daughters [cm].......: ";
	fvZeroCutsNames[1] = "vZero:: max fit Chi Square between two daughters.: ";
	fvZeroCutsNames[2] = "vZero:: min decay length [cm]....................: ";
	fvZeroCutsNames[3] = "vZero:: max decay length [cm]....................: ";
	fvZeroCutsNames[4] = "vZero:: max opening angle between daughters [rad]: ";
	fvZeroCutsNames[5] = "vZero:: pt min [Gev/c]...........................: ";
	fvZeroCutsNames[6] = "vZero:: |Eta| range <............................: ";
	
	
	return;
}

//--------------------------------------------------------------------------
void AliHFAssociatedTrackCuts::PrintAll()
{
	printf("\nCuts for the associated track: \n \n");
	       
	printf("ITS Refit........................................: %s\n",fESDTrackCuts->GetRequireITSRefit() ? "Yes" : "No");
	printf("TPC Refit........................................: %s\n",fESDTrackCuts->GetRequireTPCRefit() ? "Yes" : "No");
	printf("ITS SA...........................................: %s\n",fESDTrackCuts->GetRequireITSStandAlone() ? "Yes" : "No");
	printf("TPC SA...........................................: %s\n",fESDTrackCuts->GetRequireTPCStandAlone() ? "Yes" : "No");
	printf("Min number of ITS clusters.......................: %d\n",fESDTrackCuts->GetMinNClustersITS());
	printf("Min number of TPC clusters.......................: %d\n",fESDTrackCuts->GetMinNClusterTPC());
	Int_t spd = fESDTrackCuts->GetClusterRequirementITS(AliESDtrackCuts::kSPD);
	if(spd==0) cout <<  "SPD..............................................: kOff"  << endl;
	if(spd==1) cout <<  "SPD..............................................: kNone"  << endl;
	if(spd==2) cout <<  "SPD..............................................: kAny"  << endl;
	if(spd==3) cout <<  "SPD..............................................: kFirst"  << endl;
	if(spd==4) cout <<  "SPD..............................................: kOnlyFirst"  << endl;
	if(spd==5) cout <<  "SPD..............................................: kSecond"  << endl;
	if(spd==6) cout <<  "SPD..............................................: kOnlySecond"  << endl;
	if(spd==7) cout <<  "SPD..............................................: kBoth"  << endl;
	
	cout <<  "Filter Bit.......................................: " << fBit  << endl;
	cout <<  "Charge...........................................: " << fCharge  << endl;
		
	for(Int_t j=0;j<fNTrackCuts;j++){
		cout << fTrackCutsNames[j] << fAODTrackCuts[j] << endl;
	}
	
	printf("\nCuts for the K0 candidates: \n \n");
	for(Int_t k=0;k<fNvZeroCuts;k++){
		cout << fvZeroCutsNames[k] <<  fAODvZeroCuts[k] << endl;
	}
	cout << " " << endl;
	PrintPoolParameters();

}

//--------------------------------------------------------------------------
void AliHFAssociatedTrackCuts::PrintPoolParameters()
{
	printf("\nEvent Pool settings: \n \n");
	
	printf("Number of zVtx Bins: %d\n", fNzVtxBins);
	printf("\nzVtx Bins:\n");
	//Double_t zVtxbinLims[fNzVtxBins+1] = fNzVtxBins;
	for(Int_t k=0; k<fNzVtxBins; k++){
		printf("Bin %d..............................................: %.1f - %.1f cm\n", k, fZvtxBins[k], fZvtxBins[k+1]);	
	}
	
	printf("\nNumber of Centrality(multiplicity) Bins: %d\n", fNCentBins);
	printf("\nCentrality(multiplicity) Bins:\n");
	for(Int_t k=0; k<fNCentBins; k++){
		printf("Bin %d..............................................: %.1f - %.1f\n", k, fCentBins[k], fCentBins[k+1]);
	}
	
	
}
Commit	Line	Data
815e6bab	1	/**************************************************************************
	2	* Copyright(c) 1998-2010, ALICE Experiment at CERN, All rights reserved. *
	3	* *
	4	* Author: The ALICE Off-line Project. *
	5	* Contributors are mentioned in the code where appropriate. *
	6	* *
	7	* Permission to use, copy, modify and distribute this software and its *
	8	* documentation strictly for non-commercial purposes is hereby granted *
	9	* without fee, provided that the above copyright notice appears in all *
	10	* copies and that both the copyright notice and this permission notice *
	11	* appear in the supporting documentation. The authors make no claims *
	12	* about the suitability of this software for any purpose. It is *
	13	* provided "as is" without express or implied warranty. *
	14	**************************************************************************/
	15
	16	/* $Id$ */
	17
	18	////////////////////////////////////////////////////////////////////////
	19	//
	20	// Base class for cuts on Associated tracks for HF Correlation analysis
	21	//
	22	// Author: S.Bjelogrlic (Utrecht) sandro.bjelogrlic@cern.ch
	23	////////////////////////////////////////////////////////////////////////
	24	#include <Riostream.h>
	25	#include "AliHFAssociatedTrackCuts.h"
	26	#include "AliAODPidHF.h"
	27	#include "AliESDtrackCuts.h"
	28	#include "AliESDtrack.h"
	29	#include "AliAODv0.h"
	30	#include "AliAODVertex.h"
	31	#include "AliAODMCParticle.h"
	32	#include "TString.h"
	33
14c8fedd	34	using std::cout;
14c8fedd	35	using std::endl;
815e6bab	36
	37	ClassImp(AliHFAssociatedTrackCuts)
	38
	39	//--------------------------------------------------------------------------
	40	AliHFAssociatedTrackCuts::AliHFAssociatedTrackCuts():
	41	AliAnalysisCuts(),
	42	fESDTrackCuts(0),
	43	fPidObj(0),
bce70c96	44
	45	fPoolMaxNEvents(0),
	46	fPoolMinNTracks(0),
	47	fMinEventsToMix(0),
	48	fNzVtxBins(0),
	49	fNzVtxBinsDim(0),
	50	fZvtxBins(0),
	51	fNCentBins(0),
	52	fNCentBinsDim(0),
	53	fCentBins(0),
	54
815e6bab	55	fNTrackCuts(0),
	56	fAODTrackCuts(0),
	57	fTrackCutsNames(0),
	58	fNvZeroCuts(0),
	59	fAODvZeroCuts(0),
c84dbedf	60	fvZeroCutsNames(0),
	61	fBit(0),
	62	fCharge(0)
815e6bab	63
	64	{
	65	//
	66	//default constructor
	67	//
	68	//
	69	//default constructor
	70	//
	71
	72	}
	73
	74	//--------------------------------------------------------------------------
	75	AliHFAssociatedTrackCuts::AliHFAssociatedTrackCuts(const char* name, const char* title):
	76	AliAnalysisCuts(name,title),
	77	fESDTrackCuts(0),
	78	fPidObj(0),
bce70c96	79
	80	fPoolMaxNEvents(0),
	81	fPoolMinNTracks(0),
	82	fMinEventsToMix(0),
	83	fNzVtxBins(0),
	84	fNzVtxBinsDim(0),
	85	fZvtxBins(0),
	86	fNCentBins(0),
	87	fNCentBinsDim(0),
	88	fCentBins(0),
	89
815e6bab	90	fNTrackCuts(0),
	91	fAODTrackCuts(0),
	92	fTrackCutsNames(0),
	93	fNvZeroCuts(0),
	94	fAODvZeroCuts(0),
c84dbedf	95	fvZeroCutsNames(0),
	96	fBit(0),
	97	fCharge(0)
815e6bab	98
	99	{
	100	//
	101	//default constructor
	102	//
	103
	104	}
	105	//--------------------------------------------------------------------------
	106	AliHFAssociatedTrackCuts::AliHFAssociatedTrackCuts(const AliHFAssociatedTrackCuts &source) :
	107	AliAnalysisCuts(source),
	108	fESDTrackCuts(source.fESDTrackCuts),
	109	fPidObj(source.fPidObj),
bce70c96	110
	111	fPoolMaxNEvents(source.fPoolMaxNEvents),
	112	fPoolMinNTracks(source.fPoolMinNTracks),
	113	fMinEventsToMix(source.fMinEventsToMix),
	114	fNzVtxBins(source.fNzVtxBins),
	115	fNzVtxBinsDim(source.fNzVtxBinsDim),
	116	fZvtxBins(source.fZvtxBins),
	117	fNCentBins(source.fNCentBins),
	118	fNCentBinsDim(source.fNCentBinsDim),
	119	fCentBins(source.fCentBins),
	120
815e6bab	121	fNTrackCuts(source.fNTrackCuts),
	122	fAODTrackCuts(source.fAODTrackCuts),
	123	fTrackCutsNames(source.fTrackCutsNames),
	124	fNvZeroCuts(source.fNvZeroCuts),
	125	fAODvZeroCuts(source.fAODvZeroCuts),
c84dbedf	126	fvZeroCutsNames(source.fvZeroCutsNames),
	127	fBit(source.fBit),
	128	fCharge(source.fCharge)
815e6bab	129	{
	130	//
	131	// copy constructor
	132	//
	133
	134
14c8fedd	135	AliInfo("AliHFAssociatedTrackCuts::Copy constructor ");
815e6bab	136	if(source.fESDTrackCuts) AddTrackCuts(source.fESDTrackCuts);
	137	if(source.fAODTrackCuts) SetAODTrackCuts(source.fAODTrackCuts);
	138	if(source.fAODvZeroCuts) SetAODvZeroCuts(source.fAODvZeroCuts);
	139	if(source.fPidObj) SetPidHF(source.fPidObj);
	140	}
	141	//--------------------------------------------------------------------------
	142	AliHFAssociatedTrackCuts &AliHFAssociatedTrackCuts::operator=(const AliHFAssociatedTrackCuts &source)
	143	{
	144	//
	145	// assignment operator
	146	//
	147	if(&source == this) return *this;
	148
	149	AliAnalysisCuts::operator=(source);
	150	fESDTrackCuts=source.fESDTrackCuts;
	151	fPidObj=source.fPidObj;
	152	fNTrackCuts=source.fNTrackCuts;
	153	fAODTrackCuts=source.fAODTrackCuts;
	154	fTrackCutsNames=source.fTrackCutsNames;
	155	fNvZeroCuts=source.fNvZeroCuts;
	156	fAODvZeroCuts=source.fAODvZeroCuts;
	157	fvZeroCutsNames=source.fvZeroCutsNames;
c84dbedf	158	fBit=source.fBit;
c84dbedf	159	fCharge=source.fCharge;
815e6bab	160
	161	return *this;
	162
	163	}
	164
	165
	166	//--------------------------------------------------------------------------
	167	AliHFAssociatedTrackCuts::~AliHFAssociatedTrackCuts()
	168	{
	169	if(fESDTrackCuts) {delete fESDTrackCuts; fESDTrackCuts = 0;}
	170	if(fPidObj) {delete fPidObj; fPidObj = 0;}
bce70c96	171	if(fZvtxBins) {delete[] fZvtxBins; fZvtxBins=0;}
bce70c96	172	if(fCentBins) {delete[] fCentBins; fCentBins=0;}
815e6bab	173	if(fAODTrackCuts) {delete[] fAODTrackCuts; fAODTrackCuts=0;}
	174	if(fTrackCutsNames) {delete[] fTrackCutsNames; fTrackCutsNames=0;}
	175	if(fAODvZeroCuts){delete[] fAODvZeroCuts; fAODvZeroCuts=0;}
	176	if(fvZeroCutsNames) {delete[] fvZeroCutsNames; fvZeroCutsNames=0;}
	177
	178	}
	179	//--------------------------------------------------------------------------
	180	Bool_t AliHFAssociatedTrackCuts::IsInAcceptance()
	181	{
	182	printf("Careful: method AliHFAssociatedTrackCuts::IsInAcceptance is not implemented yet \n");
	183	return kFALSE;
	184	}
	185	//--------------------------------------------------------------------------
bce70c96	186	Bool_t AliHFAssociatedTrackCuts::IsHadronSelected(AliAODTrack * track)
815e6bab	187	{
815e6bab	188	AliESDtrack esdtrack(track);
bce70c96	189	if(!fESDTrackCuts->IsSelected(&esdtrack)) return kFALSE;
c84dbedf	190
	191	if(fBit && !track->TestFilterBit(fBit)) return kFALSE; // check the filter bit
	192
bce70c96	193	return kTRUE;
815e6bab	194
bce70c96	195	}
	196
	197	//--------------------------------------------------------------------------
	198	Bool_t AliHFAssociatedTrackCuts::CheckHadronKinematic(Double_t pt, Double_t d0)
	199	{
815e6bab	200
815e6bab	201
815e6bab	202
bce70c96	203	if(pt < fAODTrackCuts[0]) return kFALSE;
	204	if(pt > fAODTrackCuts[1]) return kFALSE;
	205	if(d0 < fAODTrackCuts[2]) return kFALSE;
	206	if(d0 > fAODTrackCuts[3]) return kFALSE;
815e6bab	207
	208	return kTRUE;
	209
	210
	211	}
c84dbedf	212	//--------------------------------------------------------------------------
	213
	214	Bool_t AliHFAssociatedTrackCuts::Charge(Short_t charge, AliAODTrack* track)
	215	{// charge is the charge to compare to (for example, a daughter of a D meson)
	216
	217	if(!fCharge) return kTRUE; // if fCharge is set to 0 (no selection on the charge), returns always true
	218	if(track->Charge()!= fCharge*charge) return kFALSE;
	219	return kTRUE;
	220	}
	221
815e6bab	222	//--------------------------------------------------------------------------
14c8fedd	223	Bool_t AliHFAssociatedTrackCuts::CheckKaonCompatibility(AliAODTrack * track, Bool_t useMc, TClonesArray* mcArray, Int_t method)
815e6bab	224	{
	225	Bool_t isKaon = kFALSE;
	226
	227	if(useMc) { // on MC
	228	Int_t hadLabel = track->GetLabel();
	229	if(hadLabel < 0) return kFALSE;
	230	AliAODMCParticle* hadron = dynamic_cast<AliAODMCParticle*>(mcArray->At(hadLabel));
	231	Int_t pdg = TMath::Abs(hadron->GetPdgCode());
	232	if (pdg == 321) isKaon = kTRUE;
	233	}
	234
	235	if(!useMc) { // on DATA
14c8fedd	236	switch(method) {
14c8fedd	237	case(1): {
815e6bab	238	Bool_t isKTPC=kFALSE;
	239	Bool_t isPiTPC=kFALSE;
	240	Bool_t isPTPC=kFALSE;
	241	Bool_t isKTOF=kFALSE;
	242	Bool_t isPiTOF=kFALSE;
	243	Bool_t isPTOF=kFALSE;
	244
	245	Bool_t KaonHyp = kFALSE;
	246	Bool_t PionHyp = kFALSE;
	247	Bool_t ProtonHyp = kFALSE;
	248
	249	if(fPidObj->CheckStatus(track,"TOF")) {
	250	isKTOF=fPidObj->IsKaonRaw(track,"TOF");
	251	isPiTOF=fPidObj->IsPionRaw(track,"TOF");
	252	isPTOF=fPidObj->IsProtonRaw(track,"TOF");
	253	}
	254	if(fPidObj->CheckStatus(track,"TPC")){
	255	isKTPC=fPidObj->IsKaonRaw(track,"TPC");
	256	isPiTPC=fPidObj->IsPionRaw(track,"TPC");
	257	isPTPC=fPidObj->IsProtonRaw(track,"TPC");
	258	}
	259
	260	if (isKTOF && isKTPC) KaonHyp = kTRUE;
	261	if (isPiTOF && isPiTPC) PionHyp = kTRUE;
	262	if (isPTOF && isPTPC) ProtonHyp = kTRUE;
	263
	264	if(KaonHyp && !PionHyp && !ProtonHyp) isKaon = kTRUE;
14c8fedd	265	break;
14c8fedd	266	}
b9f81146	267	case(2): {
b9f81146	268	if(fPidObj->MakeRawPid(track,3)>=1) isKaon = kTRUE;
14c8fedd	269	break;
14c8fedd	270	}
b9f81146	271	}
815e6bab	272	}
	273
	274	return isKaon;
	275
	276	}
	277	//--------------------------------------------------------------------------
	278	Bool_t AliHFAssociatedTrackCuts::IsKZeroSelected(AliAODv0 vzero, AliAODVertex vtx1)
	279	{
	280
	281	if(vzero->DcaV0Daughters()>fAODvZeroCuts[0]) return kFALSE;
	282	if(vzero->Chi2V0()>fAODvZeroCuts[1]) return kFALSE;
	283	if(vzero->DecayLength(vtx1) < fAODvZeroCuts[2]) return kFALSE;
	284	if(vzero->DecayLength(vtx1) > fAODvZeroCuts[3]) return kFALSE;
	285	if(vzero->OpenAngleV0() > fAODvZeroCuts[4]) return kFALSE;
	286	if(vzero->Pt() < fAODvZeroCuts[5]) return kFALSE;
	287	if(TMath::Abs(vzero->Eta()) > fAODvZeroCuts[6]) return kFALSE;
	288
	289
	290	return kTRUE;
	291	}
	292	//--------------------------------------------------------------------------
02f91926	293	Int_t AliHFAssociatedTrackCuts::IsMCpartFromHF(Int_t label, TClonesArray*mcArray){
02f91926	294	// Check origin in MC
815e6bab	295
02f91926	296	AliAODMCParticle* mcParticle;
	297	Int_t pdgCode = -1;
	298	if (label<0) return 0;
	299	Bool_t isCharmy = kFALSE;
	300	Bool_t isBeauty = kFALSE;
	301	Bool_t isD = kFALSE;
	302	Bool_t isB = kFALSE;
815e6bab	303
02f91926	304	while(pdgCode!=2212){ // loops back to the collision to check the particle source
	305
	306	mcParticle = dynamic_cast<AliAODMCParticle*>(mcArray->At(label));
bf658d5c	307	if(!mcParticle) {AliError("NO MC PARTICLE"); break;}
02f91926	308	pdgCode = TMath::Abs(mcParticle->GetPdgCode());
	309
	310	label = mcParticle->GetMother();
	311
	312
	313	if((pdgCode>=400 && pdgCode <500) \|\| (pdgCode>=4000 && pdgCode<5000 )) isD = kTRUE;
	314	if((pdgCode>=500 && pdgCode <600) \|\| (pdgCode>=5000 && pdgCode<6000 )) {isD = kFALSE; isB = kTRUE;}
	315
	316
	317	if(pdgCode == 4) isCharmy = kTRUE;
	318	if(pdgCode == 5) {isBeauty = kTRUE; isCharmy = kFALSE;}
	319	if(label<0) break;
	320
	321	}
	322	if (isCharmy && isD) return 44; // the first 4(5) indicates the charm/beauty quark, the second the charmy meson
	323	if (isBeauty && isD) return 54;
	324	if (isBeauty && isB) return 55;
	325	return 0;
815e6bab	326	}
bce70c96	327
	328	//--------------------------------------------------------------------------
	329	Bool_t AliHFAssociatedTrackCuts::InvMassDstarRejection(AliAODRecoDecayHF2Prong* d, AliAODTrack *track, Int_t hypD0) const {
	330	//
	331	// Calculates invmass of track+D0 and rejects if compatible with D*
	332	// (to remove pions from D*)
	333	//
	334	Double_t nsigma = 3.;
	335
	336	Double_t mD0, mD0bar;
	337	d->InvMassD0(mD0,mD0bar);
	338
	339	Double_t invmassDstar1 = 0, invmassDstar2 = 0;
	340	Double_t e1Pi = d->EProng(0,211), e2K = d->EProng(1,321); //hyp 1 (pi,K) - D0
	341	Double_t e1K = d->EProng(0,321), e2Pi = d->EProng(1,211); //hyp 2 (K,pi) - D0bar
	342	Double_t psum2 = (d->Px()+track->Px())*(d->Px()+track->Px())
	343	+(d->Py()+track->Py())*(d->Py()+track->Py())
	344	+(d->Pz()+track->Pz())*(d->Pz()+track->Pz());
	345
	346	switch(hypD0) {
	347	case 1:
	348	invmassDstar1 = TMath::Sqrt(pow(e1Pi+e2K+track->E(0.1396),2.)-psum2);
	349	if ((TMath::Abs(invmassDstar1-mD0)-0.14543) < nsigma800.pow(10.,-6.)) return kFALSE;
	350	break;
	351	case 2:
	352	invmassDstar2 = TMath::Sqrt(pow(e2Pi+e1K+track->E(0.1396),2.)-psum2);
	353	if ((TMath::Abs(invmassDstar2-mD0bar)-0.14543) < nsigma800.pow(10.,-6.)) return kFALSE;
	354	break;
	355	case 3:
	356	invmassDstar1 = TMath::Sqrt(pow(e1Pi+e2K+track->E(0.1396),2.)-psum2);
	357	invmassDstar2 = TMath::Sqrt(pow(e2Pi+e1K+track->E(0.1396),2.)-psum2);
	358	if ((TMath::Abs(invmassDstar1-mD0)-0.14543) < nsigma800.pow(10.,-6.)) return kFALSE;
	359	if ((TMath::Abs(invmassDstar2-mD0bar)-0.14543) < nsigma800.pow(10.,-6.)) return kFALSE;
	360	break;
	361	}
	362
	363	return kTRUE;
	364	}
815e6bab	365	//________________________________________________________
	366	void AliHFAssociatedTrackCuts::SetAODTrackCuts(Float_t *cutsarray)
	367	{
bce70c96	368	cout << "Check 1" << endl;
815e6bab	369	if(!fAODTrackCuts) fAODTrackCuts = new Float_t[fNTrackCuts];
bce70c96	370	cout << "Check 2" << endl;
815e6bab	371	for(Int_t i =0; i<fNTrackCuts; i++){
bce70c96	372	cout << "Check 2." << i << endl;
815e6bab	373	fAODTrackCuts[i] = cutsarray[i];
815e6bab	374	}
bce70c96	375	cout << "Check 3" << endl;
815e6bab	376	SetTrackCutsNames();
bce70c96	377	cout << "Check 4" << endl;
815e6bab	378	return;
	379	}
	380	//________________________________________________________
	381	void AliHFAssociatedTrackCuts::SetTrackCutsNames(/TString namearray*/){
	382
	383	fTrackCutsNames = new TString[4];
	384	fTrackCutsNames[0]= "associated track:: pt min [GeV/c]................: ";
	385	fTrackCutsNames[1]= "associated track:: pt max [GeV/c]................: ";
	386	fTrackCutsNames[2]= "associated track:: d0 min [cm]...................: ";
	387	fTrackCutsNames[3]= "associated track:: d0 max [cm]...................: ";
	388
	389
	390
	391	return;
	392	}
	393	//--------------------------------------------------------------------------
	394	void AliHFAssociatedTrackCuts::SetAODvZeroCuts(Float_t *cutsarray)
	395	{
	396
	397
	398	if(!fAODvZeroCuts) fAODvZeroCuts = new Float_t[fNvZeroCuts];
	399	for(Int_t i =0; i<fNvZeroCuts; i++){
	400	fAODvZeroCuts[i] = cutsarray[i];
	401	}
	402	SetvZeroCutsNames();
	403	return;
	404	}
	405	//--------------------------------------------------------------------------
	406	void AliHFAssociatedTrackCuts::SetvZeroCutsNames(/TString namearray*/){
	407
	408	fvZeroCutsNames = new TString[7];
	409	fvZeroCutsNames[0] = "vZero:: max DCA between two daughters [cm].......: ";
	410	fvZeroCutsNames[1] = "vZero:: max fit Chi Square between two daughters.: ";
	411	fvZeroCutsNames[2] = "vZero:: min decay length [cm]....................: ";
	412	fvZeroCutsNames[3] = "vZero:: max decay length [cm]....................: ";
	413	fvZeroCutsNames[4] = "vZero:: max opening angle between daughters [rad]: ";
	414	fvZeroCutsNames[5] = "vZero:: pt min [Gev/c]...........................: ";
	415	fvZeroCutsNames[6] = "vZero:: \|Eta\| range <............................: ";
	416
	417
	418	return;
	419	}
	420
	421	//--------------------------------------------------------------------------
	422	void AliHFAssociatedTrackCuts::PrintAll()
	423	{
	424	printf("\nCuts for the associated track: \n \n");
	425
	426	printf("ITS Refit........................................: %s\n",fESDTrackCuts->GetRequireITSRefit() ? "Yes" : "No");
	427	printf("TPC Refit........................................: %s\n",fESDTrackCuts->GetRequireTPCRefit() ? "Yes" : "No");
	428	printf("ITS SA...........................................: %s\n",fESDTrackCuts->GetRequireITSStandAlone() ? "Yes" : "No");
	429	printf("TPC SA...........................................: %s\n",fESDTrackCuts->GetRequireTPCStandAlone() ? "Yes" : "No");
	430	printf("Min number of ITS clusters.......................: %d\n",fESDTrackCuts->GetMinNClustersITS());
	431	printf("Min number of TPC clusters.......................: %d\n",fESDTrackCuts->GetMinNClusterTPC());
	432	Int_t spd = fESDTrackCuts->GetClusterRequirementITS(AliESDtrackCuts::kSPD);
	433	if(spd==0) cout << "SPD..............................................: kOff" << endl;
	434	if(spd==1) cout << "SPD..............................................: kNone" << endl;
	435	if(spd==2) cout << "SPD..............................................: kAny" << endl;
	436	if(spd==3) cout << "SPD..............................................: kFirst" << endl;
	437	if(spd==4) cout << "SPD..............................................: kOnlyFirst" << endl;
	438	if(spd==5) cout << "SPD..............................................: kSecond" << endl;
	439	if(spd==6) cout << "SPD..............................................: kOnlySecond" << endl;
	440	if(spd==7) cout << "SPD..............................................: kBoth" << endl;
c84dbedf	441
	442	cout << "Filter Bit.......................................: " << fBit << endl;
	443	cout << "Charge...........................................: " << fCharge << endl;
815e6bab	444
	445	for(Int_t j=0;j<fNTrackCuts;j++){
	446	cout << fTrackCutsNames[j] << fAODTrackCuts[j] << endl;
	447	}
	448
	449	printf("\nCuts for the K0 candidates: \n \n");
	450	for(Int_t k=0;k<fNvZeroCuts;k++){
	451	cout << fvZeroCutsNames[k] << fAODvZeroCuts[k] << endl;
	452	}
	453	cout << " " << endl;
bce70c96	454	PrintPoolParameters();
bce70c96	455
815e6bab	456	}
14c8fedd	457
14c8fedd	458	//--------------------------------------------------------------------------
bce70c96	459	void AliHFAssociatedTrackCuts::PrintPoolParameters()
	460	{
	461	printf("\nEvent Pool settings: \n \n");
	462
	463	printf("Number of zVtx Bins: %d\n", fNzVtxBins);
	464	printf("\nzVtx Bins:\n");
	465	//Double_t zVtxbinLims[fNzVtxBins+1] = fNzVtxBins;
	466	for(Int_t k=0; k<fNzVtxBins; k++){
c84dbedf	467	printf("Bin %d..............................................: %.1f - %.1f cm\n", k, fZvtxBins[k], fZvtxBins[k+1]);
bce70c96	468	}
	469
	470	printf("\nNumber of Centrality(multiplicity) Bins: %d\n", fNCentBins);
	471	printf("\nCentrality(multiplicity) Bins:\n");
	472	for(Int_t k=0; k<fNCentBins; k++){
c84dbedf	473	printf("Bin %d..............................................: %.1f - %.1f\n", k, fCentBins[k], fCentBins[k+1]);
bce70c96	474	}
	475
	476
	477
14c8fedd	478	}
bce70c96	479
bce70c96	480