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

/**************************************************************************
 * Copyright(c) 1998-2009, 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.                  *
 **************************************************************************/
//
//
//             Base class for Heavy Flavour Correlations Analysis
//             Single Event and Mixed Event Analysis are implemented
//
//-----------------------------------------------------------------------
//          
//
//						   Author S.Bjelogrlic
//                         Utrecht University 
//                      sandro.bjelogrlic@cern.ch
//
//-----------------------------------------------------------------------

/* $Id$ */

#include <TParticle.h>
#include <TVector3.h>
#include <TChain.h>
#include "TROOT.h"
#include "AliHFCorrelator.h"
#include "AliRDHFCutsDStartoKpipi.h"
#include "AliHFAssociatedTrackCuts.h"
#include "AliEventPoolManager.h"
#include "AliReducedParticle.h"
#include "AliCentrality.h"

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

//_____________________________________________________
AliHFCorrelator::AliHFCorrelator() :
//
// default constructor
//
TNamed(),
fPoolMgr(0x0),         
fPool(0x0),
fhadcuts(0x0),
fAODEvent(0x0),
fAssociatedTracks(0x0),
fmcArray(0x0),
fReducedPart(0x0),
fD0cand(0x0), 
fhypD0(0), 
fDCharge(0),

fmixing(kFALSE),
fmontecarlo(kFALSE),
fsystem(kFALSE),
fselect(0),
fUseImpactParameter(0),
fPIDmode(0),

fNofTracks(0),
fPoolContent(0),

fPhiMin(0),
fPhiMax(0),

fPtTrigger(0),
fPhiTrigger(0),
fEtaTrigger(0),


fDeltaPhi(0),
fDeltaEta(0),
fk0InvMass(0)

{
	// default constructor	
}


//_____________________________________________________
AliHFCorrelator::AliHFCorrelator(const Char_t* name, AliHFAssociatedTrackCuts *cuts, Bool_t ppOrPbPb) :
TNamed(name,"title"),
fPoolMgr(0x0),         
fPool(0x0),
fhadcuts(0x0),
fAODEvent(0x0),
fAssociatedTracks(0x0),
fmcArray(0x0),
fReducedPart(0x0),
fD0cand(0x0), 
fhypD0(0),
fDCharge(0),

fmixing(kFALSE),
fmontecarlo(kFALSE),
fsystem(ppOrPbPb),
fselect(0),
fUseImpactParameter(0),
fPIDmode(0),

fNofTracks(0),
fPoolContent(0),

fPhiMin(0),
fPhiMax(0),

fPtTrigger(0),
fPhiTrigger(0),
fEtaTrigger(0),


fDeltaPhi(0),
fDeltaEta(0),
fk0InvMass(0)
{
	fhadcuts = cuts; 
}

//_____________________________________________________
AliHFCorrelator::~AliHFCorrelator() 
{
//
// destructor
//	
	
	if(fPoolMgr)  {delete fPoolMgr; fPoolMgr=0;}       
	if(fPool) {delete fPool; fPool=0;}
	if(fhadcuts) {delete fhadcuts; fhadcuts=0;}
	if(fAODEvent) {delete fAODEvent; fAODEvent=0;}
	if(fAssociatedTracks) {delete fAssociatedTracks; fAssociatedTracks=0;}
	if(fmcArray) {delete fmcArray; fmcArray=0;}
	if(fReducedPart) {delete fReducedPart; fReducedPart=0;}
	if(fD0cand) {delete fD0cand; fD0cand=0;}
	
	
	if(fNofTracks) fNofTracks = 0;
	
	if(fPhiMin) fPhiMin = 0;
	if(fPhiMax) fPhiMax = 0;
	
	if(fPtTrigger) fPtTrigger=0;
	if(fPhiTrigger) fPhiTrigger=0;
	if(fEtaTrigger) fEtaTrigger=0;
	
	if(fDeltaPhi) fDeltaPhi=0;
	if(fDeltaEta) fDeltaEta=0;
	
	if(fk0InvMass) fk0InvMass=0;
	
}
//_____________________________________________________
Bool_t AliHFCorrelator::DefineEventPool(){
	// definition of the Pool Manager for Event Mixing
	

	Int_t MaxNofEvents = fhadcuts->GetMaxNEventsInPool();
	Int_t MinNofTracks = fhadcuts->GetMinNTracksInPool();
	Int_t NofCentBins = fhadcuts->GetNCentPoolBins();
	Double_t * CentBins = fhadcuts->GetCentPoolBins();
	Int_t NofZVrtxBins = fhadcuts->GetNZvtxPoolBins();
	Double_t *ZVrtxBins = fhadcuts->GetZvtxPoolBins();
		
			
	fPoolMgr = new AliEventPoolManager(MaxNofEvents, MinNofTracks, NofCentBins, CentBins, NofZVrtxBins, ZVrtxBins);
	if(!fPoolMgr) return kFALSE;
	return kTRUE;
}
//_____________________________________________________
Bool_t AliHFCorrelator::Initialize(){
	
    //  std::cout << "AliHFCorrelator::Initialize"<< std::endl;
  AliInfo("AliHFCorrelator::Initialize") ;
	if(!fAODEvent)
    AliInfo("No AOD event") ;
    //std::cout << "No AOD event" << std::endl;
	
	AliCentrality *centralityObj = 0;
	Int_t multiplicity = -1;
	Double_t MultipOrCent = -1;
	
	// initialize the pool for event mixing
	if(!fsystem){ // pp
	multiplicity = fAODEvent->GetNTracks();
		MultipOrCent = multiplicity; // convert from Int_t to Double_t
	}
	if(fsystem){ // PbPb
		
		centralityObj = fAODEvent->GetHeader()->GetCentralityP();
		MultipOrCent = centralityObj->GetCentralityPercentileUnchecked("V0M");
		AliInfo(Form("Centrality is %f", MultipOrCent));
	}
	
	AliAODVertex *vtx = fAODEvent->GetPrimaryVertex();
	Double_t zvertex = vtx->GetZ(); // zvertex
	Double_t * CentBins = fhadcuts->GetCentPoolBins();
	Double_t poolmin=CentBins[0];
	Double_t poolmax=CentBins[fhadcuts->GetNCentPoolBins()];

	
		if(TMath::Abs(zvertex)>=10 || MultipOrCent>poolmax || MultipOrCent < poolmin) {
		if(!fsystem)AliInfo(Form("pp Event with Zvertex = %.2f cm and multiplicity = %.0f out of pool bounds, SKIPPING",zvertex,MultipOrCent));
		if(fsystem) AliInfo(Form("PbPb Event with Zvertex = %.2f cm and centrality = %.1f  out of pool bounds, SKIPPING",zvertex,MultipOrCent));

			return kFALSE;
		}
	
	fPool = fPoolMgr->GetEventPool(MultipOrCent, zvertex);
	
	if (!fPool){
		AliInfo(Form("No pool found for multiplicity = %f, zVtx = %f cm", MultipOrCent, zvertex));
	    return kFALSE;
	}
	//fPool->PrintInfo();
	return kTRUE;
}

//_____________________________________________________
Bool_t AliHFCorrelator::ProcessEventPool(){
	 // analysis on Mixed Events
	//cout << "AliHFCorrelator::ProcessEventPool"<< endl;
		if(!fmixing) return kFALSE;
		if(!fPool->IsReady()) return kFALSE;
		if(fPool->GetCurrentNEvents()<fhadcuts->GetMinEventsToMix()) return kFALSE;
	//	fPool->PrintInfo();
		fPoolContent = fPool->GetCurrentNEvents();
		
		return kTRUE;
	
}

//_____________________________________________________
Bool_t AliHFCorrelator::ProcessAssociatedTracks(Int_t EventLoopIndex){

	fAssociatedTracks = new TObjArray();
	
	if(!fmixing){ // analysis on Single Event
		
		
		if(fselect==1 || fselect ==2)	fAssociatedTracks = AcceptAndReduceTracks(fAODEvent);
		if(fselect==3) {fAssociatedTracks = AcceptAndReduceKZero(fAODEvent);}	
		
		
	}
	
	if(fmixing) { // analysis on Mixed Events
		
			
		fAssociatedTracks = fPool->GetEvent(EventLoopIndex);
				
				
	} // end if mixing
	
	if(!fAssociatedTracks) return kFALSE;
	
	fNofTracks = fAssociatedTracks->GetEntriesFast(); 
		
	return kTRUE;
	
}
//_____________________________________________________
Bool_t AliHFCorrelator::Correlate(Int_t loopindex){

	if(loopindex >= fNofTracks) return kFALSE;
	if(!fAssociatedTracks) return kFALSE;
	
	fReducedPart = (AliReducedParticle*)fAssociatedTracks->At(loopindex);
	

	fDeltaPhi = SetCorrectPhiRange(fPhiTrigger - fReducedPart->Phi());
	
	fDeltaEta = fEtaTrigger - fReducedPart->Eta();

	return kTRUE;
	
}
		
//_____________________________________________________
Bool_t AliHFCorrelator::PoolUpdate(){

	if(!fmixing) return kFALSE;
	if(!fPool) return kFALSE;
	if(fmixing) { // update the pool for Event Mixing
		TObjArray* objArr = NULL;
		if(fselect==1 || fselect==2) objArr = (TObjArray*)AcceptAndReduceTracks(fAODEvent);
		if(fselect==3) objArr = (TObjArray*)AcceptAndReduceKZero(fAODEvent);
		if(objArr->GetEntriesFast()>0) fPool->UpdatePool(objArr); // updating the pool only if there are entries in the array
	}
		
	return kTRUE;
	
}
		
//_____________________________________________________
Double_t AliHFCorrelator::SetCorrectPhiRange(Double_t phi){
	Double_t pi = TMath::Pi();
	
	if(phi<fPhiMin) phi = phi + 2*pi;
	if(phi>fPhiMax) phi = phi - 2*pi;
	
	return phi;
}

//_____________________________________________________
TObjArray*  AliHFCorrelator::AcceptAndReduceTracks(AliAODEvent* inputEvent){
	
	Int_t nTracks = inputEvent->GetNTracks();
	AliAODVertex * vtx = inputEvent->GetPrimaryVertex();
	Double_t Bz = inputEvent->GetMagneticField();
	
	
	TObjArray* tracksClone = new TObjArray;
	tracksClone->SetOwner(kTRUE);
	for (Int_t iTrack=0; iTrack<nTracks; ++iTrack) {
		
		AliAODTrack* track = inputEvent->GetTrack(iTrack);
		if (!track) continue;
		if(!fhadcuts->IsHadronSelected(track)) continue; // apply ESD level selections
        if(!fhadcuts->Charge(fDCharge,track)) continue; // apply selection on charge, if required
		
		Double_t pT = track->Pt();
		
		//compute impact parameter
		Double_t d0z0[2],covd0z0[3];
		Double_t d0=-999999.;
		if(fUseImpactParameter) track->PropagateToDCA(vtx,Bz,100,d0z0,covd0z0);
		else d0z0[0] = 1. ; // random number - be careful with the cuts you applied
		
		if(fUseImpactParameter==1) d0 = TMath::Abs(d0z0[0]); // use impact parameter
		if(fUseImpactParameter==2) { // use impact parameter over resolution
		  if(TMath::Abs(covd0z0[0])>0.00000001) d0 = TMath::Abs(d0z0[0])/TMath::Sqrt(covd0z0[0]); 
		  else d0 = -1.; // if the resoultion is Zero, rejects the track - to be on the safe side
			
		}
		
        
		if(!fhadcuts->CheckHadronKinematic(pT,d0)) continue; // apply kinematic cuts
		Bool_t rejectsoftpi = kFALSE;
		if(fD0cand) rejectsoftpi = fhadcuts->InvMassDstarRejection(fD0cand,track,fhypD0);
		
		
		if(fselect ==2){	
			if(!fhadcuts->CheckKaonCompatibility(track,fmontecarlo,fmcArray,fPIDmode)) continue; // check if it is a Kaon - data and MC
		}
		
		tracksClone->Add(new AliReducedParticle(track->Eta(), track->Phi(), pT,track->GetLabel(),track->GetID(),d0,rejectsoftpi,track->Charge()));
	}
	
	
	return tracksClone;
}

//_____________________________________________________
TObjArray*  AliHFCorrelator::AcceptAndReduceKZero(AliAODEvent* inputEvent){
	
	Int_t nOfVZeros = inputEvent->GetNumberOfV0s();
	TObjArray* KZeroClone = new TObjArray;
	AliAODVertex *vertex1 = (AliAODVertex*)inputEvent->GetPrimaryVertex();
	Int_t v0label = -1;
	Int_t pdgDgK0toPipi[2] = {211,211};
	Double_t mPDGK0=0.497614;//TDatabasePDG::Instance()->GetParticle(310)->Mass();
	const Int_t size = inputEvent->GetNumberOfV0s();
	Int_t prevnegID[size];
	Int_t prevposID[size];
	for(Int_t iv0 =0; iv0< nOfVZeros; iv0++){// loop on all v0 candidates
		AliAODv0 *v0 = (static_cast<AliAODEvent*> (inputEvent))->GetV0(iv0);
		if(!v0) {
      AliInfo(Form("is not a v0 at step, %d", iv0)) ;
        //cout << "is not a v0 at step " << iv0 << endl;
      continue;
    }
		
		if(!fhadcuts->IsKZeroSelected(v0,vertex1)) continue; // check if it is a k0
	    
		// checks to avoid double counting
		Int_t negID = -999;
		Int_t posID = -998;
		//Int_t a = 0;
		prevnegID[iv0] = -997;
		prevposID[iv0]= -996;
		negID = v0->GetNegID();
		posID = v0->GetPosID();
		Bool_t DoubleCounting = kFALSE;
		
		for(Int_t k=0; k<iv0; k++){
			if(((negID==prevnegID[k])&&(posID==prevposID[k]))||((negID==prevposID[k])&&(posID==prevnegID[k]))){
				DoubleCounting = kTRUE;
				//a=k;
				break;
			}//end if
		} // end for
		
		if(DoubleCounting) continue;
		else{ 
			prevposID[iv0] = posID; 
			prevnegID[iv0] = negID;
		}
		
		if(fmontecarlo)	v0label = v0->MatchToMC(310,fmcArray, 0, pdgDgK0toPipi); //match a K0 short
		Double_t k0pt = v0->Pt();
		Double_t k0eta = v0->Eta();
		Double_t k0Phi = v0->Phi();
	    fk0InvMass = v0->MassK0Short();	
		
		//if (loopindex == 0) {
		//	if(!plotassociation) ((TH2F*)fOutput->FindObject("KZeroSpectra"))->Fill(k0InvMass,k0pt); // spectra for all k0
		//	if(plotassociation) ((TH2F*)fOutput->FindObject("KZeroSpectraifHF"))->Fill(k0InvMass,k0pt); // spectra for k0 in association with a D*
		//}
		// if there are more D* candidates per event, loopindex == 0 makes sure you fill the mass spectra only once!
		
		if(TMath::Abs(fk0InvMass-mPDGK0)>3*0.004) continue; // select candidates within 3 sigma
		KZeroClone->Add(new AliReducedParticle(k0eta,k0Phi,k0pt,v0label));
		
	}
	
	return KZeroClone;
}
Commit	Line	Data
bce70c96	1	/**************************************************************************
	2	* Copyright(c) 1998-2009, ALICE Experiment at CERN, All rights reserved. *
	3	* *
	4	* Author: The ALICE Off-line Project. *
	5	* Contributors are mentioned in the code where appropriate. *
	6	* *
	7	* Permission to use, copy, modify and distribute this software and its *
	8	* documentation strictly for non-commercial purposes is hereby granted *
	9	* without fee, provided that the above copyright notice appears in all *
	10	* copies and that both the copyright notice and this permission notice *
	11	* appear in the supporting documentation. The authors make no claims *
	12	* about the suitability of this software for any purpose. It is *
	13	* provided "as is" without express or implied warranty. *
	14	**************************************************************************/
	15	//
	16	//
	17	// Base class for Heavy Flavour Correlations Analysis
	18	// Single Event and Mixed Event Analysis are implemented
	19	//
	20	//-----------------------------------------------------------------------
	21	//
	22	//
	23	// Author S.Bjelogrlic
	24	// Utrecht University
	25	// sandro.bjelogrlic@cern.ch
	26	//
	27	//-----------------------------------------------------------------------
	28
	29	/* $Id$ */
	30
	31	#include <TParticle.h>
	32	#include <TVector3.h>
	33	#include <TChain.h>
	34	#include "TROOT.h"
	35	#include "AliHFCorrelator.h"
	36	#include "AliRDHFCutsDStartoKpipi.h"
	37	#include "AliHFAssociatedTrackCuts.h"
	38	#include "AliEventPoolManager.h"
	39	#include "AliReducedParticle.h"
	40	#include "AliCentrality.h"
	41
e61afb80	42	using std::cout;
	43	using std::endl;
	44
bce70c96	45	//_____________________________________________________
	46	AliHFCorrelator::AliHFCorrelator() :
	47	//
	48	// default constructor
	49	//
	50	TNamed(),
	51	fPoolMgr(0x0),
	52	fPool(0x0),
	53	fhadcuts(0x0),
	54	fAODEvent(0x0),
	55	fAssociatedTracks(0x0),
	56	fmcArray(0x0),
	57	fReducedPart(0x0),
	58	fD0cand(0x0),
	59	fhypD0(0),
c84dbedf	60	fDCharge(0),
bce70c96	61
	62	fmixing(kFALSE),
	63	fmontecarlo(kFALSE),
	64	fsystem(kFALSE),
	65	fselect(0),
	66	fUseImpactParameter(0),
	67	fPIDmode(0),
	68
	69	fNofTracks(0),
	70	fPoolContent(0),
	71
	72	fPhiMin(0),
	73	fPhiMax(0),
	74
	75	fPtTrigger(0),
	76	fPhiTrigger(0),
	77	fEtaTrigger(0),
	78
	79
	80	fDeltaPhi(0),
	81	fDeltaEta(0),
	82	fk0InvMass(0)
	83
	84	{
	85	// default constructor
	86	}
	87
	88
	89
	90	//_____________________________________________________
	91	AliHFCorrelator::AliHFCorrelator(const Char_t* name, AliHFAssociatedTrackCuts *cuts, Bool_t ppOrPbPb) :
	92	TNamed(name,"title"),
	93	fPoolMgr(0x0),
	94	fPool(0x0),
	95	fhadcuts(0x0),
	96	fAODEvent(0x0),
	97	fAssociatedTracks(0x0),
	98	fmcArray(0x0),
	99	fReducedPart(0x0),
	100	fD0cand(0x0),
	101	fhypD0(0),
c84dbedf	102	fDCharge(0),
bce70c96	103
	104	fmixing(kFALSE),
	105	fmontecarlo(kFALSE),
	106	fsystem(ppOrPbPb),
	107	fselect(0),
	108	fUseImpactParameter(0),
	109	fPIDmode(0),
	110
	111	fNofTracks(0),
	112	fPoolContent(0),
	113
	114	fPhiMin(0),
	115	fPhiMax(0),
	116
	117	fPtTrigger(0),
	118	fPhiTrigger(0),
	119	fEtaTrigger(0),
	120
	121
	122	fDeltaPhi(0),
	123	fDeltaEta(0),
	124	fk0InvMass(0)
	125	{
	126	fhadcuts = cuts;
	127	}
	128
	129	//_____________________________________________________
	130	AliHFCorrelator::~AliHFCorrelator()
	131	{
	132	//
	133	// destructor
	134	//
	135
	136	if(fPoolMgr) {delete fPoolMgr; fPoolMgr=0;}
	137	if(fPool) {delete fPool; fPool=0;}
	138	if(fhadcuts) {delete fhadcuts; fhadcuts=0;}
	139	if(fAODEvent) {delete fAODEvent; fAODEvent=0;}
	140	if(fAssociatedTracks) {delete fAssociatedTracks; fAssociatedTracks=0;}
	141	if(fmcArray) {delete fmcArray; fmcArray=0;}
	142	if(fReducedPart) {delete fReducedPart; fReducedPart=0;}
	143	if(fD0cand) {delete fD0cand; fD0cand=0;}
	144
	145
	146	if(fNofTracks) fNofTracks = 0;
	147
	148	if(fPhiMin) fPhiMin = 0;
	149	if(fPhiMax) fPhiMax = 0;
	150
	151	if(fPtTrigger) fPtTrigger=0;
	152	if(fPhiTrigger) fPhiTrigger=0;
	153	if(fEtaTrigger) fEtaTrigger=0;
	154
	155	if(fDeltaPhi) fDeltaPhi=0;
	156	if(fDeltaEta) fDeltaEta=0;
	157
	158	if(fk0InvMass) fk0InvMass=0;
	159
	160	}
	161	//_____________________________________________________
	162	Bool_t AliHFCorrelator::DefineEventPool(){
	163	// definition of the Pool Manager for Event Mixing
	164
	165
	166	Int_t MaxNofEvents = fhadcuts->GetMaxNEventsInPool();
167	Int_t MinNofTracks = fhadcuts->GetMinNTracksInPool();
168	Int_t NofCentBins = fhadcuts->GetNCentPoolBins();
169	Double_t * CentBins = fhadcuts->GetCentPoolBins();
170	Int_t NofZVrtxBins = fhadcuts->GetNZvtxPoolBins();
171	Double_t *ZVrtxBins = fhadcuts->GetZvtxPoolBins();
172
173
174	fPoolMgr = new AliEventPoolManager(MaxNofEvents, MinNofTracks, NofCentBins, CentBins, NofZVrtxBins, ZVrtxBins);
175	if(!fPoolMgr) return kFALSE;
176	return kTRUE;
177	}
178	//_____________________________________________________
179	Bool_t AliHFCorrelator::Initialize(){
180
2ad98fc5	181	// std::cout << "AliHFCorrelator::Initialize"<< std::endl;
	182	AliInfo("AliHFCorrelator::Initialize") ;
	183	if(!fAODEvent)
	184	AliInfo("No AOD event") ;
	185	//std::cout << "No AOD event" << std::endl;
bce70c96	186
	187	AliCentrality *centralityObj = 0;
	188	Int_t multiplicity = -1;
	189	Double_t MultipOrCent = -1;
	190
	191	// initialize the pool for event mixing
	192	if(!fsystem){ // pp
	193	multiplicity = fAODEvent->GetNTracks();
	194	MultipOrCent = multiplicity; // convert from Int_t to Double_t
	195	}
	196	if(fsystem){ // PbPb
	197
	198	centralityObj = fAODEvent->GetHeader()->GetCentralityP();
	199	MultipOrCent = centralityObj->GetCentralityPercentileUnchecked("V0M");
	200	AliInfo(Form("Centrality is %f", MultipOrCent));
	201	}
	202
	203	AliAODVertex *vtx = fAODEvent->GetPrimaryVertex();
	204	Double_t zvertex = vtx->GetZ(); // zvertex
	205	Double_t * CentBins = fhadcuts->GetCentPoolBins();
	206	Double_t poolmin=CentBins[0];
	207	Double_t poolmax=CentBins[fhadcuts->GetNCentPoolBins()];
c84dbedf	208
bce70c96	209
	210	if(TMath::Abs(zvertex)>=10 \|\| MultipOrCent>poolmax \|\| MultipOrCent < poolmin) {
	211	if(!fsystem)AliInfo(Form("pp Event with Zvertex = %.2f cm and multiplicity = %.0f out of pool bounds, SKIPPING",zvertex,MultipOrCent));
c84dbedf	212	if(fsystem) AliInfo(Form("PbPb Event with Zvertex = %.2f cm and centrality = %.1f out of pool bounds, SKIPPING",zvertex,MultipOrCent));
bce70c96	213
	214	return kFALSE;
	215	}
	216
	217	fPool = fPoolMgr->GetEventPool(MultipOrCent, zvertex);
	218
	219	if (!fPool){
	220	AliInfo(Form("No pool found for multiplicity = %f, zVtx = %f cm", MultipOrCent, zvertex));
	221	return kFALSE;
	222	}
c84dbedf	223	//fPool->PrintInfo();
bce70c96	224	return kTRUE;
	225	}
	226
	227	//_____________________________________________________
	228	Bool_t AliHFCorrelator::ProcessEventPool(){
	229	// analysis on Mixed Events
c84dbedf	230	//cout << "AliHFCorrelator::ProcessEventPool"<< endl;
bce70c96	231	if(!fmixing) return kFALSE;
	232	if(!fPool->IsReady()) return kFALSE;
	233	if(fPool->GetCurrentNEvents()<fhadcuts->GetMinEventsToMix()) return kFALSE;
c84dbedf	234	// fPool->PrintInfo();
bce70c96	235	fPoolContent = fPool->GetCurrentNEvents();
	236
	237	return kTRUE;
	238
	239	}
	240
	241	//_____________________________________________________
	242	Bool_t AliHFCorrelator::ProcessAssociatedTracks(Int_t EventLoopIndex){
	243
	244	fAssociatedTracks = new TObjArray();
	245
	246	if(!fmixing){ // analysis on Single Event
	247
	248
	249	if(fselect==1 \|\| fselect ==2) fAssociatedTracks = AcceptAndReduceTracks(fAODEvent);
	250	if(fselect==3) {fAssociatedTracks = AcceptAndReduceKZero(fAODEvent);}
	251
	252
	253	}
	254
	255	if(fmixing) { // analysis on Mixed Events
	256
	257
	258	fAssociatedTracks = fPool->GetEvent(EventLoopIndex);
	259
	260
	261
	262
	263	} // end if mixing
	264
	265	if(!fAssociatedTracks) return kFALSE;
	266
	267	fNofTracks = fAssociatedTracks->GetEntriesFast();
	268
	269	return kTRUE;
	270
	271	}
	272	//_____________________________________________________
	273	Bool_t AliHFCorrelator::Correlate(Int_t loopindex){
	274
	275	if(loopindex >= fNofTracks) return kFALSE;
	276	if(!fAssociatedTracks) return kFALSE;
	277
	278	fReducedPart = (AliReducedParticle*)fAssociatedTracks->At(loopindex);
	279
	280
	281	fDeltaPhi = SetCorrectPhiRange(fPhiTrigger - fReducedPart->Phi());
	282
	283	fDeltaEta = fEtaTrigger - fReducedPart->Eta();
	284
	285	return kTRUE;
	286
	287	}
	288
	289	//_____________________________________________________
	290	Bool_t AliHFCorrelator::PoolUpdate(){
	291
	292	if(!fmixing) return kFALSE;
	293	if(!fPool) return kFALSE;
	294	if(fmixing) { // update the pool for Event Mixing
c84dbedf	295	TObjArray* objArr = NULL;
	296	if(fselect==1 \|\| fselect==2) objArr = (TObjArray*)AcceptAndReduceTracks(fAODEvent);
	297	if(fselect==3) objArr = (TObjArray*)AcceptAndReduceKZero(fAODEvent);
	298	if(objArr->GetEntriesFast()>0) fPool->UpdatePool(objArr); // updating the pool only if there are entries in the array
bce70c96	299	}
c84dbedf	300
bce70c96	301	return kTRUE;
c84dbedf	302
bce70c96	303	}
	304
	305	//_____________________________________________________
	306	Double_t AliHFCorrelator::SetCorrectPhiRange(Double_t phi){
	307	Double_t pi = TMath::Pi();
	308
	309	if(phi<fPhiMin) phi = phi + 2*pi;
	310	if(phi>fPhiMax) phi = phi - 2*pi;
	311
	312	return phi;
	313	}
	314
	315	//_____________________________________________________
	316	TObjArray* AliHFCorrelator::AcceptAndReduceTracks(AliAODEvent* inputEvent){
	317
	318	Int_t nTracks = inputEvent->GetNTracks();
	319	AliAODVertex * vtx = inputEvent->GetPrimaryVertex();
	320	Double_t Bz = inputEvent->GetMagneticField();
	321
	322
	323	TObjArray* tracksClone = new TObjArray;
	324	tracksClone->SetOwner(kTRUE);
	325	for (Int_t iTrack=0; iTrack<nTracks; ++iTrack) {
	326
	327	AliAODTrack* track = inputEvent->GetTrack(iTrack);
	328	if (!track) continue;
	329	if(!fhadcuts->IsHadronSelected(track)) continue; // apply ESD level selections
c84dbedf	330	if(!fhadcuts->Charge(fDCharge,track)) continue; // apply selection on charge, if required
c84dbedf	331
bce70c96	332	Double_t pT = track->Pt();
	333
	334	//compute impact parameter
	335	Double_t d0z0[2],covd0z0[3];
	336	Double_t d0=-999999.;
	337	if(fUseImpactParameter) track->PropagateToDCA(vtx,Bz,100,d0z0,covd0z0);
	338	else d0z0[0] = 1. ; // random number - be careful with the cuts you applied
	339
	340	if(fUseImpactParameter==1) d0 = TMath::Abs(d0z0[0]); // use impact parameter
	341	if(fUseImpactParameter==2) { // use impact parameter over resolution
bec72d8c	342	if(TMath::Abs(covd0z0[0])>0.00000001) d0 = TMath::Abs(d0z0[0])/TMath::Sqrt(covd0z0[0]);
bce70c96	343	else d0 = -1.; // if the resoultion is Zero, rejects the track - to be on the safe side
	344
	345	}
	346
c84dbedf	347
bce70c96	348
bce70c96	349
c84dbedf	350	if(!fhadcuts->CheckHadronKinematic(pT,d0)) continue; // apply kinematic cuts
bce70c96	351	Bool_t rejectsoftpi = kFALSE;
	352	if(fD0cand) rejectsoftpi = fhadcuts->InvMassDstarRejection(fD0cand,track,fhypD0);
	353
	354
	355	if(fselect ==2){
	356	if(!fhadcuts->CheckKaonCompatibility(track,fmontecarlo,fmcArray,fPIDmode)) continue; // check if it is a Kaon - data and MC
	357	}
	358
c84dbedf	359	tracksClone->Add(new AliReducedParticle(track->Eta(), track->Phi(), pT,track->GetLabel(),track->GetID(),d0,rejectsoftpi,track->Charge()));
bce70c96	360	}
	361
	362
	363	return tracksClone;
	364	}
	365
	366	//_____________________________________________________
	367	TObjArray* AliHFCorrelator::AcceptAndReduceKZero(AliAODEvent* inputEvent){
	368
	369	Int_t nOfVZeros = inputEvent->GetNumberOfV0s();
	370	TObjArray* KZeroClone = new TObjArray;
	371	AliAODVertex vertex1 = (AliAODVertex)inputEvent->GetPrimaryVertex();
	372	Int_t v0label = -1;
	373	Int_t pdgDgK0toPipi[2] = {211,211};
	374	Double_t mPDGK0=0.497614;//TDatabasePDG::Instance()->GetParticle(310)->Mass();
	375	const Int_t size = inputEvent->GetNumberOfV0s();
	376	Int_t prevnegID[size];
	377	Int_t prevposID[size];
	378	for(Int_t iv0 =0; iv0< nOfVZeros; iv0++){// loop on all v0 candidates
	379	AliAODv0 v0 = (static_cast<AliAODEvent> (inputEvent))->GetV0(iv0);
2ad98fc5	380	if(!v0) {
	381	AliInfo(Form("is not a v0 at step, %d", iv0)) ;
	382	//cout << "is not a v0 at step " << iv0 << endl;
	383	continue;
	384	}
bce70c96	385
	386	if(!fhadcuts->IsKZeroSelected(v0,vertex1)) continue; // check if it is a k0
	387
	388	// checks to avoid double counting
	389	Int_t negID = -999;
	390	Int_t posID = -998;
	391	//Int_t a = 0;
	392	prevnegID[iv0] = -997;
	393	prevposID[iv0]= -996;
	394	negID = v0->GetNegID();
	395	posID = v0->GetPosID();
	396	Bool_t DoubleCounting = kFALSE;
	397
	398	for(Int_t k=0; k<iv0; k++){
	399	if(((negID==prevnegID[k])&&(posID==prevposID[k]))\|\|((negID==prevposID[k])&&(posID==prevnegID[k]))){
	400	DoubleCounting = kTRUE;
	401	//a=k;
	402	break;
	403	}//end if
	404	} // end for
	405
c84dbedf	406	if(DoubleCounting) continue;
	407	else{
	408	prevposID[iv0] = posID;
	409	prevnegID[iv0] = negID;
	410	}
bce70c96	411
	412	if(fmontecarlo) v0label = v0->MatchToMC(310,fmcArray, 0, pdgDgK0toPipi); //match a K0 short
	413	Double_t k0pt = v0->Pt();
	414	Double_t k0eta = v0->Eta();
	415	Double_t k0Phi = v0->Phi();
	416	fk0InvMass = v0->MassK0Short();
	417
	418	//if (loopindex == 0) {
	419	// if(!plotassociation) ((TH2F*)fOutput->FindObject("KZeroSpectra"))->Fill(k0InvMass,k0pt); // spectra for all k0
	420	// if(plotassociation) ((TH2F)fOutput->FindObject("KZeroSpectraifHF"))->Fill(k0InvMass,k0pt); // spectra for k0 in association with a D
	421	//}
	422	// if there are more D* candidates per event, loopindex == 0 makes sure you fill the mass spectra only once!
	423
	424	if(TMath::Abs(fk0InvMass-mPDGK0)>3*0.004) continue; // select candidates within 3 sigma
	425	KZeroClone->Add(new AliReducedParticle(k0eta,k0Phi,k0pt,v0label));
	426
	427	}
	428
	429	return KZeroClone;
	430	}