[u/mrichter/AliRoot.git] / PWGLF / SPECTRA / AntiprotonToProton / AliProtonAbsorptionCorrection.cxx

#include <Riostream.h>
#include <TFile.h>
#include <TSystem.h>
#include <TF1.h>
#include <TH2D.h>
#include <TH1D.h>
#include <TH1I.h>
#include <TParticle.h>
#include <string.h>

#include <AliExternalTrackParam.h>
#include <AliAODEvent.h>
#include <AliESDEvent.h>
//#include <AliLog.h>
#include <AliPID.h>
#include <AliStack.h>
#include <AliCFContainer.h>
#include <AliCFEffGrid.h>
#include <AliCFDataGrid.h>
#include <AliCFManager.h>
#include <AliCFTrackKineCuts.h>
#include <AliCFParticleGenCuts.h>
#include <AliCFAcceptanceCuts.h>
#include <AliMCEvent.h>
//#include <AliESDVertex.h>
class AliLog;
class AliESDVertex;

#include "AliProtonAbsorptionCorrection.h"
#include "AliProtonAnalysisBase.h"

ClassImp(AliProtonAbsorptionCorrection)

//____________________________________________________________________//
AliProtonAbsorptionCorrection::AliProtonAbsorptionCorrection() : 
  TObject(), fProtonAnalysisBase(0),
  fNBinsY(0), fMinY(0), fMaxY(0),
  fNBinsPt(0), fMinPt(0), fMaxPt(0),
  fCFManagerProtons(0), fCFManagerAntiProtons(0) {
  //fProtonContainer(0), fAntiProtonContainer(0) {
  //Default constructor
}

//____________________________________________________________________//
AliProtonAbsorptionCorrection::~AliProtonAbsorptionCorrection() {
  //Default destructor
  if(fCFManagerProtons) delete fCFManagerProtons;
  if(fCFManagerAntiProtons) delete fCFManagerAntiProtons;
  if(fProtonAnalysisBase) delete fProtonAnalysisBase;
  //if(fProtonContainer) delete fProtonContainer;
  //if(fAntiProtonContainer) delete fAntiProtonContainer;
}

//____________________________________________________________________//
void AliProtonAbsorptionCorrection::InitAnalysisHistograms(Int_t nbinsY, 
							   Float_t fLowY, 
							   Float_t fHighY, 
							   Int_t nbinsPt, 
							   Float_t fLowPt, 
							   Float_t fHighPt) {
  //Initializes the histograms
  fNBinsY = nbinsY;
  fMinY = fLowY;
  fMaxY = fHighY;
  fNBinsPt = nbinsPt;
  fMinPt = fLowPt;
  fMaxPt = fHighPt;
  const Int_t    mintrackrefsTPC = 1;

  //=========================================================//
  //setting up the containers
  Int_t iBin[2];
  iBin[0] = nbinsY;
  iBin[1] = nbinsPt;
  Double_t *binLimY = new Double_t[nbinsY+1];
  Double_t *binLimPt = new Double_t[nbinsPt+1];
  //values for bin lower bounds
  for(Int_t i = 0; i <= nbinsY; i++) 
    binLimY[i]=(Double_t)fLowY  + (fHighY - fLowY)  /nbinsY*(Double_t)i;
  for(Int_t i = 0; i <= nbinsPt; i++) 
    binLimPt[i]=(Double_t)fLowPt  + (fHighPt - fLowPt)  /nbinsPt*(Double_t)i;

  //Proton container
  AliCFContainer *containerProtons = new AliCFContainer("containerProtons",
							"container for protons",
							kNSteps,2,iBin);
  containerProtons->SetBinLimits(0,binLimY); //rapidity or eta
  containerProtons->SetBinLimits(1,binLimPt); //pT

  //Anti-proton container
  AliCFContainer *containerAntiProtons = new AliCFContainer("containerAntiProtons",
							    "container for antiprotons",
							    kNSteps,2,iBin);
  containerAntiProtons->SetBinLimits(0,binLimY); //rapidity or eta
  containerAntiProtons->SetBinLimits(1,binLimPt); //pT
  
  //=========================================================//
  //Setting up the criteria for the generated particles
  Char_t titleCF[256] = "mcKineCutsProtons";
  Char_t nameCF[256] = "MC-level kinematic cuts";
  AliCFTrackKineCuts *mcKineCutsProtons = new AliCFTrackKineCuts(titleCF,nameCF);
  mcKineCutsProtons->SetPtRange(fMinPt,fMaxPt);
  if(fProtonAnalysisBase->GetEtaMode()) 
    mcKineCutsProtons->SetEtaRange(fMinY,fMaxY);
  else
    mcKineCutsProtons->SetRapidityRange(fMinY,fMaxY);
  mcKineCutsProtons->SetChargeMC(1.0);

  strncpy(titleCF,"mcKineCutsAntiProtons",22);
  AliCFTrackKineCuts *mcKineCutsAntiProtons = new AliCFTrackKineCuts(titleCF,nameCF);
  mcKineCutsAntiProtons->SetPtRange(fMinPt,fMaxPt);
  if(fProtonAnalysisBase->GetEtaMode()) 
    mcKineCutsAntiProtons->SetEtaRange(fMinY,fMaxY);
  else
    mcKineCutsAntiProtons->SetRapidityRange(fMinY,fMaxY);
  mcKineCutsAntiProtons->SetChargeMC(-1.0);

  strncpy(titleCF,"mcGenCutsProtons",17);
  strncpy(nameCF,"MC particle generation cuts",28);
  AliCFParticleGenCuts* mcGenCutsProtons = new AliCFParticleGenCuts(titleCF,nameCF);
  mcGenCutsProtons->SetRequireIsPrimary();
  mcGenCutsProtons->SetRequirePdgCode(2212);

  strncpy(titleCF,"mcGenCutsAntiProtons",21);
  AliCFParticleGenCuts* mcGenCutsAntiProtons = new AliCFParticleGenCuts(titleCF,nameCF);
  mcGenCutsAntiProtons->SetRequireIsPrimary();
  mcGenCutsAntiProtons->SetRequirePdgCode(-2212);

  TObjArray* mcListProtons = new TObjArray(0);
  mcListProtons->AddLast(mcKineCutsProtons);
  mcListProtons->AddLast(mcGenCutsProtons);
  TObjArray* mcListAntiProtons = new TObjArray(0);
  mcListAntiProtons->AddLast(mcKineCutsAntiProtons);
  mcListAntiProtons->AddLast(mcGenCutsAntiProtons);

  //Setting up the criteria for the reconstructible particles
  strncpy(titleCF,"mcAccCuts",10);
  strncpy(nameCF,"Acceptance cuts",16);
  AliCFAcceptanceCuts *mcAccCuts = new AliCFAcceptanceCuts(titleCF,nameCF);
  mcAccCuts->SetMinNHitTPC(mintrackrefsTPC);
  TObjArray* accList = new TObjArray(0);
  accList->AddLast(mcAccCuts);

  //____________________________________________//
  //Setting up the criteria for the reconstructed tracks
  //____________________________________________//
  strncpy(titleCF,"recKineCutsProtons",19);
  strncpy(nameCF,"rec-level kine cuts",20);
  AliCFTrackKineCuts *recKineCutsProtons = new AliCFTrackKineCuts(titleCF,nameCF);
  recKineCutsProtons->SetPtRange(fMinPt,fMaxPt);
  if(fProtonAnalysisBase->GetEtaMode()) 
    recKineCutsProtons->SetEtaRange(fMinY,fMaxY);
  else
    recKineCutsProtons->SetRapidityRange(fMinY,fMaxY);
  recKineCutsProtons->SetChargeRec(1.0);

  //____________________________________________//
  strncpy(titleCF,"recKineCutsAntiProtons",23);
  AliCFTrackKineCuts *recKineCutsAntiProtons = new AliCFTrackKineCuts(titleCF,nameCF);
  recKineCutsAntiProtons->SetPtRange(fMinPt,fMaxPt);
  if(fProtonAnalysisBase->GetEtaMode()) 
    recKineCutsAntiProtons->SetEtaRange(fMinY,fMinY);
  else
    recKineCutsAntiProtons->SetRapidityRange(fMinY,fMaxY);
  recKineCutsAntiProtons->SetChargeRec(-1.0);

  //____________________________________________//
  TObjArray* recListProtons = new TObjArray(0);
  recListProtons->AddLast(recKineCutsProtons);
  recListProtons->AddLast(mcGenCutsProtons);

  TObjArray* recListAntiProtons = new TObjArray(0);
  recListAntiProtons->AddLast(recKineCutsAntiProtons);
  recListAntiProtons->AddLast(mcGenCutsAntiProtons);

  //=========================================================//
  //CF manager - Protons
  fCFManagerProtons = new AliCFManager();
  fCFManagerProtons->SetParticleContainer(containerProtons);
  fCFManagerProtons->SetParticleCutsList(AliCFManager::kPartGenCuts,mcListProtons);
  fCFManagerProtons->SetParticleCutsList(AliCFManager::kPartAccCuts,accList);
  fCFManagerProtons->SetParticleCutsList(AliCFManager::kPartRecCuts,recListProtons);

  //CF manager - Protons
  fCFManagerAntiProtons = new AliCFManager();
  fCFManagerAntiProtons->SetParticleContainer(containerAntiProtons);
  fCFManagerAntiProtons->SetParticleCutsList(AliCFManager::kPartGenCuts,mcListAntiProtons);
  fCFManagerAntiProtons->SetParticleCutsList(AliCFManager::kPartAccCuts,accList);
  fCFManagerAntiProtons->SetParticleCutsList(AliCFManager::kPartRecCuts,recListAntiProtons);
}

//_________________________________________________________________________//
void AliProtonAbsorptionCorrection::FillAbsorptionMaps(AliESDEvent *esd, 
						       //const AliESDVertex *vertex,
						       AliMCEvent *mcEvent) {	
  //Function to fill the correction containers
  fCFManagerProtons->SetMCEventInfo(mcEvent);
  fCFManagerAntiProtons->SetMCEventInfo(mcEvent);
 
  //Dummy container
  Double_t containerInput[2];
  //__________________________________________________________//
  //loop on the MC event
  for (Int_t ipart = 0; ipart < mcEvent->GetNumberOfTracks(); ipart++) { 
    AliMCParticle *mcPart  = dynamic_cast<AliMCParticle *>(mcEvent->GetTrack(ipart));
    if(!mcPart) continue;

    Double_t vz = mcPart->Zv();
    if (TMath::Abs(vz) > 50.) continue;//exclude particles generated out of the acceptance
    if(TMath::Abs(mcPart->Eta()) > 1.0) continue;

    //Protons
    //check the MC-level cuts
    if (fCFManagerProtons->CheckParticleCuts(AliCFManager::kPartGenCuts,
					     mcPart)) {
      if(fProtonAnalysisBase->GetEtaMode()) 
	containerInput[0] = (Float_t)mcPart->Eta();
      else 
	containerInput[0] = (Float_t)fProtonAnalysisBase->Rapidity(mcPart->Px(),mcPart->Pz(),mcPart->Pz());
      containerInput[1] = (Float_t)mcPart->Pt();
      //fill the container for Gen-level selection
      fCFManagerProtons->GetParticleContainer()->Fill(containerInput,
						      kStepGenerated);
      //check the Acceptance-level cuts
      if (fCFManagerProtons->CheckParticleCuts(AliCFManager::kPartAccCuts,
					       mcPart)) {
	//fill the container for Acceptance-level selection
	fCFManagerProtons->GetParticleContainer()->Fill(containerInput,
							kStepReconstructible);
      }//acceptance cuts - protons
    }//MC level cuts - protons

    //Antiprotons
    //check the MC-level cuts
    if (fCFManagerAntiProtons->CheckParticleCuts(AliCFManager::kPartGenCuts,
						 mcPart)) {
      if(fProtonAnalysisBase->GetEtaMode()) 
	containerInput[0] = (Float_t)mcPart->Eta();
      else
	containerInput[0] = (Float_t)fProtonAnalysisBase->Rapidity(mcPart->Px(),mcPart->Pz(),mcPart->Pz());
      containerInput[1] = (Float_t)mcPart->Pt();
      //fill the container for Gen-level selection
      fCFManagerAntiProtons->GetParticleContainer()->Fill(containerInput,
							  kStepGenerated);
      //check the Acceptance-level cuts
      if (fCFManagerAntiProtons->CheckParticleCuts(AliCFManager::kPartAccCuts,
						   mcPart)) {
	//fill the container for Acceptance-level selection
	fCFManagerAntiProtons->GetParticleContainer()->Fill(containerInput,
							    kStepReconstructible);
      }//acceptance cuts - antiprotons
    }//MC level cuts - antiprotons
  }//loop over MC particles
  
  
  //__________________________________________________________//
  //ESD track loop
  for (Int_t iTrack = 0; iTrack < esd->GetNumberOfTracks(); iTrack++) {
    AliESDtrack *track = dynamic_cast<AliESDtrack *>(esd->GetTrack(iTrack));
    if(!track) continue;

    // is track associated to particle ?
    Int_t label = track->GetLabel();
    if (label < 0) continue;
    if(label > mcEvent->GetNumberOfTracks()) continue;

    AliMCParticle *mcPart  = dynamic_cast<AliMCParticle *>(mcEvent->GetTrack(label));
    if(!mcPart) continue;

    Double_t vz = mcPart->Zv();
    if (TMath::Abs(vz) > 50.) continue;//exclude particles generated out of the acceptance
    if(TMath::Abs(mcPart->Eta()) > 1.0) continue;

    if((fProtonAnalysisBase->GetAnalysisMode()==AliProtonAnalysisBase::kTPC)||(fProtonAnalysisBase->GetAnalysisMode()==AliProtonAnalysisBase::kHybrid)) {
      AliExternalTrackParam *tpcTrack = (AliExternalTrackParam *)track->GetTPCInnerParam();
      if(!tpcTrack) continue;
    }//Hybrid or standalone TPC
    
    //Protons
    // check if this track was part of the signal
    if (fCFManagerProtons->CheckParticleCuts(AliCFManager::kPartGenCuts,
					     mcPart)) {
      //fill the container
      if(fProtonAnalysisBase->GetEtaMode()) 
	containerInput[0] = mcPart->Eta();
      else 
	containerInput[0] = (Float_t)fProtonAnalysisBase->Rapidity(mcPart->Px(),mcPart->Pz(),mcPart->Pz());
      containerInput[1] = mcPart->Pt();
      fCFManagerProtons->GetParticleContainer()->Fill(containerInput,
						      kStepReconstructed) ;   
    }//MC primaries - protons
    
    //Antiprotons    
    // check if this track was part of the signal
    if (fCFManagerAntiProtons->CheckParticleCuts(AliCFManager::kPartGenCuts,
						 mcPart)) {
      //fill the container
      if(fProtonAnalysisBase->GetEtaMode()) 
	containerInput[0] = mcPart->Eta();
      else 
	containerInput[0] = (Float_t)fProtonAnalysisBase->Rapidity(mcPart->Px(),mcPart->Pz(),mcPart->Pz());
      containerInput[1] = mcPart->Pt();
      fCFManagerAntiProtons->GetParticleContainer()->Fill(containerInput,
							  kStepReconstructed);
    }//MC primaries - antiprotons
  }//track loop
  
  //if(fProtonAnalysisBase->GetDebugMode())
  //Printf("Initial number of tracks: %d | Identified (anti)protons: %d - %d | Survived (anti)protons: %d - %d",nTracks,nIdentifiedProtons,nIdentifiedAntiProtons,nSurvivedProtons,nSurvivedAntiProtons);
}

//_________________________________________________________________________//
void AliProtonAbsorptionCorrection::FillAbsorptionMaps(AliAODEvent *fAOD) {
  // Analysis from AOD: to be implemented...
  // in the near future.
  fAOD->Print();
}
Commit	Line	Data
8c0ff853	1	#include <Riostream.h>
	2	#include <TFile.h>
	3	#include <TSystem.h>
	4	#include <TF1.h>
	5	#include <TH2D.h>
	6	#include <TH1D.h>
	7	#include <TH1I.h>
	8	#include <TParticle.h>
d7a38708	9	#include <string.h>
8c0ff853	10
	11	#include <AliExternalTrackParam.h>
	12	#include <AliAODEvent.h>
	13	#include <AliESDEvent.h>
	14	//#include <AliLog.h>
	15	#include <AliPID.h>
	16	#include <AliStack.h>
	17	#include <AliCFContainer.h>
	18	#include <AliCFEffGrid.h>
	19	#include <AliCFDataGrid.h>
	20	#include <AliCFManager.h>
	21	#include <AliCFTrackKineCuts.h>
	22	#include <AliCFParticleGenCuts.h>
	23	#include <AliCFAcceptanceCuts.h>
	24	#include <AliMCEvent.h>
	25	//#include <AliESDVertex.h>
	26	class AliLog;
	27	class AliESDVertex;
	28
	29	#include "AliProtonAbsorptionCorrection.h"
	30	#include "AliProtonAnalysisBase.h"
	31
	32	ClassImp(AliProtonAbsorptionCorrection)
	33
	34	//____________________________________________________________________//
	35	AliProtonAbsorptionCorrection::AliProtonAbsorptionCorrection() :
	36	TObject(), fProtonAnalysisBase(0),
	37	fNBinsY(0), fMinY(0), fMaxY(0),
	38	fNBinsPt(0), fMinPt(0), fMaxPt(0),
	39	fCFManagerProtons(0), fCFManagerAntiProtons(0) {
	40	//fProtonContainer(0), fAntiProtonContainer(0) {
	41	//Default constructor
	42	}
	43
	44	//____________________________________________________________________//
	45	AliProtonAbsorptionCorrection::~AliProtonAbsorptionCorrection() {
	46	//Default destructor
	47	if(fCFManagerProtons) delete fCFManagerProtons;
	48	if(fCFManagerAntiProtons) delete fCFManagerAntiProtons;
	49	if(fProtonAnalysisBase) delete fProtonAnalysisBase;
	50	//if(fProtonContainer) delete fProtonContainer;
	51	//if(fAntiProtonContainer) delete fAntiProtonContainer;
	52	}
	53
	54	//____________________________________________________________________//
	55	void AliProtonAbsorptionCorrection::InitAnalysisHistograms(Int_t nbinsY,
	56	Float_t fLowY,
	57	Float_t fHighY,
	58	Int_t nbinsPt,
	59	Float_t fLowPt,
	60	Float_t fHighPt) {
	61	//Initializes the histograms
	62	fNBinsY = nbinsY;
	63	fMinY = fLowY;
	64	fMaxY = fHighY;
	65	fNBinsPt = nbinsPt;
	66	fMinPt = fLowPt;
	67	fMaxPt = fHighPt;
	68	const Int_t mintrackrefsTPC = 1;
	69
	70	//=========================================================//
	71	//setting up the containers
	72	Int_t iBin[2];
	73	iBin[0] = nbinsY;
74	iBin[1] = nbinsPt;
75	Double_t *binLimY = new Double_t[nbinsY+1];
76	Double_t *binLimPt = new Double_t[nbinsPt+1];
77	//values for bin lower bounds
78	for(Int_t i = 0; i <= nbinsY; i++)
79	binLimY[i]=(Double_t)fLowY + (fHighY - fLowY) /nbinsY*(Double_t)i;
80	for(Int_t i = 0; i <= nbinsPt; i++)
81	binLimPt[i]=(Double_t)fLowPt + (fHighPt - fLowPt) /nbinsPt*(Double_t)i;
82
83	//Proton container
84	AliCFContainer *containerProtons = new AliCFContainer("containerProtons",
85	"container for protons",
86	kNSteps,2,iBin);
87	containerProtons->SetBinLimits(0,binLimY); //rapidity or eta
88	containerProtons->SetBinLimits(1,binLimPt); //pT
89
90	//Anti-proton container
91	AliCFContainer *containerAntiProtons = new AliCFContainer("containerAntiProtons",
92	"container for antiprotons",
93	kNSteps,2,iBin);
94	containerAntiProtons->SetBinLimits(0,binLimY); //rapidity or eta
95	containerAntiProtons->SetBinLimits(1,binLimPt); //pT
96
97	//=========================================================//
98	//Setting up the criteria for the generated particles
d7a38708	99	Char_t titleCF[256] = "mcKineCutsProtons";
d7a38708	100	Char_t nameCF[256] = "MC-level kinematic cuts";
86014430	101	AliCFTrackKineCuts *mcKineCutsProtons = new AliCFTrackKineCuts(titleCF,nameCF);
8c0ff853	102	mcKineCutsProtons->SetPtRange(fMinPt,fMaxPt);
	103	if(fProtonAnalysisBase->GetEtaMode())
	104	mcKineCutsProtons->SetEtaRange(fMinY,fMaxY);
	105	else
	106	mcKineCutsProtons->SetRapidityRange(fMinY,fMaxY);
	107	mcKineCutsProtons->SetChargeMC(1.0);
	108
93a33311	109	strncpy(titleCF,"mcKineCutsAntiProtons",22);
86014430	110	AliCFTrackKineCuts *mcKineCutsAntiProtons = new AliCFTrackKineCuts(titleCF,nameCF);
8c0ff853	111	mcKineCutsAntiProtons->SetPtRange(fMinPt,fMaxPt);
	112	if(fProtonAnalysisBase->GetEtaMode())
	113	mcKineCutsAntiProtons->SetEtaRange(fMinY,fMaxY);
	114	else
	115	mcKineCutsAntiProtons->SetRapidityRange(fMinY,fMaxY);
	116	mcKineCutsAntiProtons->SetChargeMC(-1.0);
	117
93a33311	118	strncpy(titleCF,"mcGenCutsProtons",17);
93a33311	119	strncpy(nameCF,"MC particle generation cuts",28);
86014430	120	AliCFParticleGenCuts* mcGenCutsProtons = new AliCFParticleGenCuts(titleCF,nameCF);
8c0ff853	121	mcGenCutsProtons->SetRequireIsPrimary();
8c0ff853	122	mcGenCutsProtons->SetRequirePdgCode(2212);
86014430	123
93a33311	124	strncpy(titleCF,"mcGenCutsAntiProtons",21);
86014430	125	AliCFParticleGenCuts* mcGenCutsAntiProtons = new AliCFParticleGenCuts(titleCF,nameCF);
8c0ff853	126	mcGenCutsAntiProtons->SetRequireIsPrimary();
	127	mcGenCutsAntiProtons->SetRequirePdgCode(-2212);
	128
	129	TObjArray* mcListProtons = new TObjArray(0);
	130	mcListProtons->AddLast(mcKineCutsProtons);
	131	mcListProtons->AddLast(mcGenCutsProtons);
	132	TObjArray* mcListAntiProtons = new TObjArray(0);
	133	mcListAntiProtons->AddLast(mcKineCutsAntiProtons);
	134	mcListAntiProtons->AddLast(mcGenCutsAntiProtons);
	135
	136	//Setting up the criteria for the reconstructible particles
93a33311	137	strncpy(titleCF,"mcAccCuts",10);
93a33311	138	strncpy(nameCF,"Acceptance cuts",16);
86014430	139	AliCFAcceptanceCuts *mcAccCuts = new AliCFAcceptanceCuts(titleCF,nameCF);
8c0ff853	140	mcAccCuts->SetMinNHitTPC(mintrackrefsTPC);
	141	TObjArray* accList = new TObjArray(0);
	142	accList->AddLast(mcAccCuts);
	143
	144	//____________________________________________//
	145	//Setting up the criteria for the reconstructed tracks
	146	//____________________________________________//
93a33311	147	strncpy(titleCF,"recKineCutsProtons",19);
93a33311	148	strncpy(nameCF,"rec-level kine cuts",20);
86014430	149	AliCFTrackKineCuts *recKineCutsProtons = new AliCFTrackKineCuts(titleCF,nameCF);
8c0ff853	150	recKineCutsProtons->SetPtRange(fMinPt,fMaxPt);
	151	if(fProtonAnalysisBase->GetEtaMode())
	152	recKineCutsProtons->SetEtaRange(fMinY,fMaxY);
	153	else
	154	recKineCutsProtons->SetRapidityRange(fMinY,fMaxY);
	155	recKineCutsProtons->SetChargeRec(1.0);
	156
	157	//____________________________________________//
93a33311	158	strncpy(titleCF,"recKineCutsAntiProtons",23);
86014430	159	AliCFTrackKineCuts *recKineCutsAntiProtons = new AliCFTrackKineCuts(titleCF,nameCF);
8c0ff853	160	recKineCutsAntiProtons->SetPtRange(fMinPt,fMaxPt);
	161	if(fProtonAnalysisBase->GetEtaMode())
	162	recKineCutsAntiProtons->SetEtaRange(fMinY,fMinY);
	163	else
	164	recKineCutsAntiProtons->SetRapidityRange(fMinY,fMaxY);
	165	recKineCutsAntiProtons->SetChargeRec(-1.0);
	166
	167	//____________________________________________//
	168	TObjArray* recListProtons = new TObjArray(0);
	169	recListProtons->AddLast(recKineCutsProtons);
	170	recListProtons->AddLast(mcGenCutsProtons);
	171
	172	TObjArray* recListAntiProtons = new TObjArray(0);
	173	recListAntiProtons->AddLast(recKineCutsAntiProtons);
	174	recListAntiProtons->AddLast(mcGenCutsAntiProtons);
	175
	176	//=========================================================//
	177	//CF manager - Protons
	178	fCFManagerProtons = new AliCFManager();
	179	fCFManagerProtons->SetParticleContainer(containerProtons);
	180	fCFManagerProtons->SetParticleCutsList(AliCFManager::kPartGenCuts,mcListProtons);
	181	fCFManagerProtons->SetParticleCutsList(AliCFManager::kPartAccCuts,accList);
	182	fCFManagerProtons->SetParticleCutsList(AliCFManager::kPartRecCuts,recListProtons);
	183
	184	//CF manager - Protons
	185	fCFManagerAntiProtons = new AliCFManager();
	186	fCFManagerAntiProtons->SetParticleContainer(containerAntiProtons);
	187	fCFManagerAntiProtons->SetParticleCutsList(AliCFManager::kPartGenCuts,mcListAntiProtons);
	188	fCFManagerAntiProtons->SetParticleCutsList(AliCFManager::kPartAccCuts,accList);
	189	fCFManagerAntiProtons->SetParticleCutsList(AliCFManager::kPartRecCuts,recListAntiProtons);
	190	}
	191
	192	//_________________________________________________________________________//
	193	void AliProtonAbsorptionCorrection::FillAbsorptionMaps(AliESDEvent *esd,
224f441a	194	//const AliESDVertex *vertex,
8c0ff853	195	AliMCEvent *mcEvent) {
8c0ff853	196	//Function to fill the correction containers
0b98b55a	197	fCFManagerProtons->SetMCEventInfo(mcEvent);
0b98b55a	198	fCFManagerAntiProtons->SetMCEventInfo(mcEvent);
8c0ff853	199
	200	//Dummy container
	201	Double_t containerInput[2];
	202	//__________________________________________________________//
	203	//loop on the MC event
	204	for (Int_t ipart = 0; ipart < mcEvent->GetNumberOfTracks(); ipart++) {
	205	AliMCParticle mcPart = dynamic_cast<AliMCParticle >(mcEvent->GetTrack(ipart));
3ef6ad84	206	if(!mcPart) continue;
	207
	208	Double_t vz = mcPart->Zv();
	209	if (TMath::Abs(vz) > 50.) continue;//exclude particles generated out of the acceptance
	210	if(TMath::Abs(mcPart->Eta()) > 1.0) continue;
8c0ff853	211
	212	//Protons
	213	//check the MC-level cuts
	214	if (fCFManagerProtons->CheckParticleCuts(AliCFManager::kPartGenCuts,
	215	mcPart)) {
6c04258d	216	if(fProtonAnalysisBase->GetEtaMode())
	217	containerInput[0] = (Float_t)mcPart->Eta();
	218	else
	219	containerInput[0] = (Float_t)fProtonAnalysisBase->Rapidity(mcPart->Px(),mcPart->Pz(),mcPart->Pz());
8c0ff853	220	containerInput[1] = (Float_t)mcPart->Pt();
	221	//fill the container for Gen-level selection
	222	fCFManagerProtons->GetParticleContainer()->Fill(containerInput,
	223	kStepGenerated);
	224	//check the Acceptance-level cuts
	225	if (fCFManagerProtons->CheckParticleCuts(AliCFManager::kPartAccCuts,
	226	mcPart)) {
	227	//fill the container for Acceptance-level selection
	228	fCFManagerProtons->GetParticleContainer()->Fill(containerInput,
	229	kStepReconstructible);
	230	}//acceptance cuts - protons
	231	}//MC level cuts - protons
	232
	233	//Antiprotons
	234	//check the MC-level cuts
	235	if (fCFManagerAntiProtons->CheckParticleCuts(AliCFManager::kPartGenCuts,
	236	mcPart)) {
6c04258d	237	if(fProtonAnalysisBase->GetEtaMode())
	238	containerInput[0] = (Float_t)mcPart->Eta();
	239	else
	240	containerInput[0] = (Float_t)fProtonAnalysisBase->Rapidity(mcPart->Px(),mcPart->Pz(),mcPart->Pz());
8c0ff853	241	containerInput[1] = (Float_t)mcPart->Pt();
	242	//fill the container for Gen-level selection
	243	fCFManagerAntiProtons->GetParticleContainer()->Fill(containerInput,
	244	kStepGenerated);
	245	//check the Acceptance-level cuts
	246	if (fCFManagerAntiProtons->CheckParticleCuts(AliCFManager::kPartAccCuts,
	247	mcPart)) {
	248	//fill the container for Acceptance-level selection
	249	fCFManagerAntiProtons->GetParticleContainer()->Fill(containerInput,
	250	kStepReconstructible);
	251	}//acceptance cuts - antiprotons
	252	}//MC level cuts - antiprotons
	253	}//loop over MC particles
	254
	255
	256	//__________________________________________________________//
	257	//ESD track loop
	258	for (Int_t iTrack = 0; iTrack < esd->GetNumberOfTracks(); iTrack++) {
	259	AliESDtrack track = dynamic_cast<AliESDtrack >(esd->GetTrack(iTrack));
	260	if(!track) continue;
	261
	262	// is track associated to particle ?
	263	Int_t label = track->GetLabel();
6c04258d	264	if (label < 0) continue;
	265	if(label > mcEvent->GetNumberOfTracks()) continue;
	266
8c0ff853	267	AliMCParticle mcPart = dynamic_cast<AliMCParticle >(mcEvent->GetTrack(label));
3ef6ad84	268	if(!mcPart) continue;
	269
	270	Double_t vz = mcPart->Zv();
	271	if (TMath::Abs(vz) > 50.) continue;//exclude particles generated out of the acceptance
	272	if(TMath::Abs(mcPart->Eta()) > 1.0) continue;
8c0ff853	273
	274	if((fProtonAnalysisBase->GetAnalysisMode()==AliProtonAnalysisBase::kTPC)\|\|(fProtonAnalysisBase->GetAnalysisMode()==AliProtonAnalysisBase::kHybrid)) {
	275	AliExternalTrackParam tpcTrack = (AliExternalTrackParam )track->GetTPCInnerParam();
	276	if(!tpcTrack) continue;
	277	}//Hybrid or standalone TPC
	278
	279	//Protons
	280	// check if this track was part of the signal
	281	if (fCFManagerProtons->CheckParticleCuts(AliCFManager::kPartGenCuts,
	282	mcPart)) {
	283	//fill the container
6c04258d	284	if(fProtonAnalysisBase->GetEtaMode())
	285	containerInput[0] = mcPart->Eta();
	286	else
	287	containerInput[0] = (Float_t)fProtonAnalysisBase->Rapidity(mcPart->Px(),mcPart->Pz(),mcPart->Pz());
8c0ff853	288	containerInput[1] = mcPart->Pt();
	289	fCFManagerProtons->GetParticleContainer()->Fill(containerInput,
	290	kStepReconstructed) ;
	291	}//MC primaries - protons
	292
	293	//Antiprotons
	294	// check if this track was part of the signal
	295	if (fCFManagerAntiProtons->CheckParticleCuts(AliCFManager::kPartGenCuts,
	296	mcPart)) {
	297	//fill the container
6c04258d	298	if(fProtonAnalysisBase->GetEtaMode())
	299	containerInput[0] = mcPart->Eta();
	300	else
	301	containerInput[0] = (Float_t)fProtonAnalysisBase->Rapidity(mcPart->Px(),mcPart->Pz(),mcPart->Pz());
8c0ff853	302	containerInput[1] = mcPart->Pt();
	303	fCFManagerAntiProtons->GetParticleContainer()->Fill(containerInput,
	304	kStepReconstructed);
	305	}//MC primaries - antiprotons
	306	}//track loop
	307
	308	//if(fProtonAnalysisBase->GetDebugMode())
	309	//Printf("Initial number of tracks: %d \| Identified (anti)protons: %d - %d \| Survived (anti)protons: %d - %d",nTracks,nIdentifiedProtons,nIdentifiedAntiProtons,nSurvivedProtons,nSurvivedAntiProtons);
	310	}
	311
	312	//_________________________________________________________________________//
	313	void AliProtonAbsorptionCorrection::FillAbsorptionMaps(AliAODEvent *fAOD) {
	314	// Analysis from AOD: to be implemented...
	315	// in the near future.
	316	fAOD->Print();
	317	}