[u/mrichter/AliRoot.git] / PWGLF / SPECTRA / PiKaPr / TestAOD / AliSpectraBothPID.cxx

#include "AliSpectraBothPID.h"
#include "AliAODEvent.h"      
#include "TH1F.h"             
#include "TH2F.h"             
#include "TList.h"            
#include "AliAODTrack.h"      
#include "AliAODMCParticle.h" 
#include "AliPIDResponse.h"   
#include "AliAnalysisManager.h"
#include "AliInputEventHandler.h"
#include "AliSpectraBothTrackCuts.h"

ClassImp(AliSpectraBothPID)

AliSpectraBothPID::AliSpectraBothPID() : TNamed("PID", "PID object"), fPIDType(kNSigmacircleTPCTOF), fNSigmaPID(3), fPIDResponse(0),fshiftTPC(0),fshiftTOF(0) 
{

}

AliSpectraBothPID::AliSpectraBothPID(BothPIDType_t pidType) : TNamed("PID", "PID object"), fPIDType(pidType), fNSigmaPID(3), fPIDResponse(0), fshiftTPC(0),fshiftTOF(0)
{


}


void AliSpectraBothPID::FillQAHistos(AliSpectraBothHistoManager * hman, AliVTrack * track, AliSpectraBothTrackCuts * trackCuts) 
{

  // fill a bunch of QA histos


  // Get PID response object, if needed
  	if(!fPIDResponse) 
  	{
    		AliAnalysisManager *man = AliAnalysisManager::GetAnalysisManager();
   		 AliInputEventHandler* inputHandler = (AliInputEventHandler*)(man->GetInputEventHandler());
    		fPIDResponse = inputHandler->GetPIDResponse();
  	}
	Bool_t ycut[3]={false,false,false};
   	ycut[0]=trackCuts->CheckYCut ((BothParticleSpecies_t) kSpPion);
  	ycut[1]=trackCuts->CheckYCut ((BothParticleSpecies_t) kSpKaon);
	ycut[2]=trackCuts->CheckYCut ((BothParticleSpecies_t) kSpProton);
 
 
  //Response
  	AliVParticle *inEvHMain = dynamic_cast<AliVParticle *>(track);
  
  	hman->GetPIDHistogram(kHistPIDTPC)->Fill(track->GetTPCmomentum(), track->GetTPCsignal()*track->Charge()); // PID histo
  	hman->GetPIDHistogram(kHistPIDTPCPion)->Fill(track->GetTPCmomentum(),fPIDResponse->GetTPCResponse().GetExpectedSignal(track->GetTPCmomentum(),AliPID::kPion)*track->Charge()); // Expected PIDPion histo
  	hman->GetPIDHistogram(kHistPIDTPCKaon)->Fill(track->GetTPCmomentum(),fPIDResponse->GetTPCResponse().GetExpectedSignal(track->GetTPCmomentum(),AliPID::kKaon)*track->Charge()); // Expected PIDKaon histo
  	hman->GetPIDHistogram(kHistPIDTPCProton)->Fill(track->GetTPCmomentum(),fPIDResponse->GetTPCResponse().GetExpectedSignal(track->GetTPCmomentum(),AliPID::kProton)*track->Charge()); // Expected PIDProton histo
  
  	Double_t nsigmaTPCkProton = fPIDResponse->NumberOfSigmasTPC(inEvHMain, AliPID::kProton)+fshiftTPC;
  	Double_t nsigmaTPCkKaon = fPIDResponse->NumberOfSigmasTPC(inEvHMain, AliPID::kKaon)+fshiftTPC; 
  	Double_t nsigmaTPCkPion = fPIDResponse->NumberOfSigmasTPC(inEvHMain, AliPID::kPion)+fshiftTPC; 
  	Double_t nsigmaTOFkProton=100.0,nsigmaTOFkKaon=100.0,nsigmaTOFkPion=100.0;

  	Double_t nsigmaTPCTOFkProton = TMath::Abs(nsigmaTPCkProton);
  	Double_t nsigmaTPCTOFkKaon   = TMath::Abs(nsigmaTPCkKaon);
  	Double_t nsigmaTPCTOFkPion   = TMath::Abs(nsigmaTPCkPion);
	if(fPIDType==kNSigmaTOF)
	{
		nsigmaTPCTOFkProton = 100.0;
  		nsigmaTPCTOFkKaon   = 100.0;
  		nsigmaTPCTOFkPion   =100.0;

	}
	Float_t tmppt=TMath::Min(trackCuts->GetPtTOFMatchingPion(),TMath::Min(trackCuts->GetPtTOFMatchingKaon(),trackCuts->GetPtTOFMatchingProton()));
 	if(track->Pt()>=tmppt&&trackCuts->CheckTOFMatchingParticleType(kSpProton)&&trackCuts->CheckTOFMatchingParticleType(kSpKaon)&&trackCuts->CheckTOFMatchingParticleType(kSpPion))
  	{
    
   		 hman->GetPIDHistogram(kHistPIDTOF)->Fill(track->P(),(track->GetTOFsignal()/100)*track->Charge()); // PID histo
    
   		 nsigmaTOFkProton = fPIDResponse->NumberOfSigmasTOF(inEvHMain, AliPID::kProton)+fshiftTOF;
  		 nsigmaTOFkKaon = fPIDResponse->NumberOfSigmasTOF(inEvHMain, AliPID::kKaon)+fshiftTOF; 
   		 nsigmaTOFkPion = fPIDResponse->NumberOfSigmasTOF(inEvHMain, AliPID::kPion)+fshiftTOF; 
    
    //TOF	 	 
		if(ycut[0])
    		{
    			hman->GetPtHistogram(kHistNSigPionTOF)->Fill(track->P(),nsigmaTOFkPion );
    			hman->GetPtHistogram(kHistNSigPionPtTOF)->Fill(track->Pt(),nsigmaTOFkPion );
    		}
    		if(ycut[1])
  	 	{	
    			hman->GetPtHistogram(kHistNSigKaonTOF)->Fill(track->P(),nsigmaTOFkKaon);
   	 		hman->GetPtHistogram(kHistNSigKaonPtTOF)->Fill(track->Pt(),nsigmaTOFkKaon );	 
   		}
		if(ycut[2])
		{					
   			hman->GetPtHistogram(kHistNSigProtonTOF)->Fill(track->P(),nsigmaTOFkProton ); 
    			hman->GetPtHistogram(kHistNSigProtonPtTOF)->Fill(track->Pt(),nsigmaTOFkProton );
		}
	
    
   		 nsigmaTPCTOFkProton = TMath::Sqrt((nsigmaTPCkProton*nsigmaTPCkProton+nsigmaTOFkProton*nsigmaTOFkProton)/2.0);
    		nsigmaTPCTOFkKaon = TMath::Sqrt((nsigmaTPCkKaon*nsigmaTPCkKaon+nsigmaTOFkKaon*nsigmaTOFkKaon)/2.0);
    		nsigmaTPCTOFkPion = TMath::Sqrt((nsigmaTPCkPion*nsigmaTPCkPion+nsigmaTOFkPion*nsigmaTOFkPion)/2.0);
  
  	}

   	if(ycut[0])
    	{
    		hman->GetPtHistogram(kHistNSigPionTPC)->Fill(track->P(),nsigmaTPCkPion );
    		hman->GetPtHistogram(kHistNSigPionPtTPC)->Fill(track->Pt(),nsigmaTPCkPion );
  		hman->GetPtHistogram(kHistNSigPionTPCTOF)->Fill(track->P(),nsigmaTPCTOFkPion);
  		hman->GetPtHistogram(kHistNSigPionPtTPCTOF)->Fill(track->Pt(),nsigmaTPCTOFkPion);
	}
    	if(ycut[1])
  	 {	
    		hman->GetPtHistogram(kHistNSigKaonTPC)->Fill(track->P(),nsigmaTPCkKaon );
   	 	hman->GetPtHistogram(kHistNSigKaonPtTPC)->Fill(track->Pt(),nsigmaTPCkKaon );	 
  		hman->GetPtHistogram(kHistNSigKaonTPCTOF)->Fill(track->P(),nsigmaTPCTOFkKaon);
  		hman->GetPtHistogram(kHistNSigKaonPtTPCTOF)->Fill(track->Pt(),nsigmaTPCTOFkKaon);
		
	}
	if(ycut[2])
	{					
   		hman->GetPtHistogram(kHistNSigProtonTPC)->Fill(track->P(),nsigmaTPCkProton ); 
    		hman->GetPtHistogram(kHistNSigProtonPtTPC)->Fill(track->Pt(),nsigmaTPCkProton );
  		hman->GetPtHistogram(kHistNSigProtonTPCTOF)->Fill(track->P(),nsigmaTPCTOFkProton);
  		hman->GetPtHistogram(kHistNSigProtonPtTPCTOF)->Fill(track->Pt(),nsigmaTPCTOFkProton);
	}

}

Int_t AliSpectraBothPID::GetParticleSpecie(AliAODMCParticle * part) 
{
  // return PID according to MC truth
  	switch(TMath::Abs(part->PdgCode()))
	{
  		case 2212:
   		return kSpProton;
   		 break;
  		case 321:
    		return kSpKaon;
    		break;
  		case 211:
    		return kSpPion;
    		break;
  		default:
    		return kSpUndefined;
  	} 
}

Int_t AliSpectraBothPID::GetParticleSpecie(TParticle* part) 
{
  // return PID according to MC truth
  switch(TMath::Abs(part->GetPdgCode()))
  {
  	case 2212:
    	return kSpProton;
   	 break;
  	case 321:
    	return kSpKaon;
    	break;
  	case 211:
    	return kSpPion;
    	break;
  	default:
    	return kSpUndefined;
  } 
}

Int_t AliSpectraBothPID::GetParticleSpecie(AliSpectraBothHistoManager * hman,AliVTrack      * trk, AliSpectraBothTrackCuts * trackCuts, Bool_t* rec) 
{
  // return PID according to detectors
  // Get PID response object, if needed

	  // guess the particle based on the smaller nsigma
  	rec[kSpPion]=false;
 	 rec[kSpKaon]=false;
  	rec[kSpProton]=false;

 	if(!fPIDResponse) 
 	 {
    		AliAnalysisManager *man = AliAnalysisManager::GetAnalysisManager();
    		AliInputEventHandler* inputHandler = (AliInputEventHandler*)(man->GetInputEventHandler());
   		 fPIDResponse = inputHandler->GetPIDResponse();
  	}

  	if(!fPIDResponse) 
	{
    		AliFatal("Cannot get pid response");
    	return 0;
  	}


  // Compute nsigma for each hypthesis
  	AliVParticle *inEvHMain = dynamic_cast<AliVParticle *>(trk);

  // --- TPC
  	Double_t nsigmaTPCkProton = TMath::Abs(fPIDResponse->NumberOfSigmasTPC(inEvHMain, AliPID::kProton)+fshiftTPC);
  	Double_t nsigmaTPCkKaon   = TMath::Abs(fPIDResponse->NumberOfSigmasTPC(inEvHMain, AliPID::kKaon)+fshiftTPC); 
 	 Double_t nsigmaTPCkPion   = TMath::Abs(fPIDResponse->NumberOfSigmasTPC(inEvHMain, AliPID::kPion)+fshiftTPC);
  
 
	Double_t nsigmaTOFkProton=0.0;
	Double_t nsigmaTOFkKaon=0.0;
	Double_t nsigmaTOFkPion=0.0;


  	Double_t nsigmaTPCTOFkProton = TMath::Abs(nsigmaTPCkProton);
  	Double_t nsigmaTPCTOFkKaon   = TMath::Abs(nsigmaTPCkKaon);
  	Double_t nsigmaTPCTOFkPion   = TMath::Abs(nsigmaTPCkPion);
	if(fPIDType==kNSigmaTOF)
	{
		nsigmaTPCTOFkProton = 10000.0;
  		nsigmaTPCTOFkKaon   = 10000.0;
  		nsigmaTPCTOFkPion   =10000.0;
		nsigmaTOFkProton = 10000.0;
  		nsigmaTOFkKaon   = 10000.0;
  		nsigmaTOFkPion   =10000.0;

	}
	

	if(trackCuts->GetUseTypeDependedTOFCut())
  	{
		if(trackCuts->CheckTOFMatchingParticleType(kSpPion)&&trk->Pt()>=trackCuts->GetPtTOFMatchingPion())
		{	
			nsigmaTOFkPion   = TMath::Abs(fPIDResponse->NumberOfSigmasTOF(inEvHMain, AliPID::kPion)+fshiftTOF); 
			nsigmaTPCTOFkPion   = TMath::Sqrt((nsigmaTPCkPion*nsigmaTPCkPion+nsigmaTOFkPion*nsigmaTOFkPion)/2.);
		}
		else if (trk->Pt()>=trackCuts->GetPtTOFMatchingPion()) // tracks without proper flags
		{
			nsigmaTOFkPion   = 10000.0; 
			nsigmaTPCTOFkPion   =10000.0;

		}	
			
		if(trackCuts->CheckTOFMatchingParticleType(kSpKaon)&&trk->Pt()>=trackCuts->GetPtTOFMatchingKaon())
		{	
			nsigmaTOFkKaon   = TMath::Abs(fPIDResponse->NumberOfSigmasTOF(inEvHMain, AliPID::kKaon)+fshiftTOF); 
			nsigmaTPCTOFkKaon   = TMath::Sqrt((nsigmaTPCkKaon*nsigmaTPCkKaon+nsigmaTOFkKaon*nsigmaTOFkKaon)/2.);
		}
		else if (trk->Pt()>=trackCuts->GetPtTOFMatchingKaon()) // tracks without proper flags
		{
			nsigmaTOFkKaon   = 10000.0; 
			nsigmaTPCTOFkKaon   =10000.0;

		}	

		if(trackCuts->CheckTOFMatchingParticleType(kSpProton)&&trk->Pt()>=trackCuts->GetPtTOFMatchingProton())
		{	
			nsigmaTOFkProton   = TMath::Abs(fPIDResponse->NumberOfSigmasTOF(inEvHMain, AliPID::kProton)+fshiftTOF); 
			nsigmaTPCTOFkProton   = TMath::Sqrt((nsigmaTPCkProton*nsigmaTPCkProton+nsigmaTOFkProton*nsigmaTOFkProton)/2.);
		}
		else if (trk->Pt()>=trackCuts->GetPtTOFMatchingProton()) // tracks without proper flags
		{
			nsigmaTOFkProton   = 10000.0; 
			nsigmaTPCTOFkProton   =10000.0;

		}	


	}
	else if(trk->Pt()>=trackCuts->GetPtTOFMatching())
  	{
    		nsigmaTOFkProton = TMath::Abs(fPIDResponse->NumberOfSigmasTOF(inEvHMain, AliPID::kProton)+fshiftTOF);
   	 	nsigmaTOFkKaon   = TMath::Abs(fPIDResponse->NumberOfSigmasTOF(inEvHMain, AliPID::kKaon)+fshiftTOF); 
    		nsigmaTOFkPion   = TMath::Abs(fPIDResponse->NumberOfSigmasTOF(inEvHMain, AliPID::kPion)+fshiftTOF); 
    		// the TOF info is used in combined
    		nsigmaTPCTOFkProton = TMath::Sqrt((nsigmaTPCkProton*nsigmaTPCkProton+nsigmaTOFkProton*nsigmaTOFkProton)/2.);
		nsigmaTPCTOFkKaon   = TMath::Sqrt((nsigmaTPCkKaon*nsigmaTPCkKaon+nsigmaTOFkKaon*nsigmaTOFkKaon)/2.);
		nsigmaTPCTOFkPion   = TMath::Sqrt((nsigmaTPCkPion*nsigmaTPCkPion+nsigmaTOFkPion*nsigmaTOFkPion)/2.);
  	}
	 else
	{ 
		if (fPIDType==kNSigmaTOF)
			return kSpUndefined;

	}
	
  	// select the nsigma to be used for the actual PID
  	Double_t nsigmaPion=10000.0; 
	Double_t nsigmaKaon=10000.0; 
	Double_t nsigmaProton=10000.0;

  	switch (fPIDType) 
	{
  		case kNSigmaTPC:
   		nsigmaProton  =  nsigmaTPCkProton;
    		nsigmaKaon	  =  nsigmaTPCkKaon  ;
    		nsigmaPion    =  nsigmaTPCkPion  ;
    		break;
  		case kNSigmaTOF:
    		nsigmaProton  =  nsigmaTOFkProton;
    		nsigmaKaon	  =  nsigmaTOFkKaon  ;
    		nsigmaPion    =  nsigmaTOFkPion  ;
    		break;
  		case kNSigmacircleTPCTOF:
    		nsigmaProton  =  nsigmaTPCTOFkProton;
    		nsigmaKaon	  =  nsigmaTPCTOFkKaon  ;
    		nsigmaPion    =  nsigmaTPCTOFkPion  ;
    		break;
		case kNSigmasquareTPCTOF:
		nsigmaProton  =  TMath::Max(nsigmaTPCkProton,nsigmaTOFkProton);
    		nsigmaKaon	  =  TMath::Max(nsigmaTPCkKaon,nsigmaTOFkKaon);
    		nsigmaPion    =  TMath::Max(nsigmaTPCkPion,nsigmaTOFkPion)  ;
    		break;

  	}	
  
  	if(nsigmaPion   < fNSigmaPID)
		rec[kSpPion]=true;
  	if(nsigmaKaon   < fNSigmaPID)
		rec[kSpKaon]=true;
 	 if(nsigmaProton   < fNSigmaPID)
		rec[kSpProton]=true;
	
  	if( ( nsigmaKaon==nsigmaPion ) && ( nsigmaKaon==nsigmaProton )) 
		return kSpUndefined;//if is the default value for the three
  	if( ( nsigmaKaon   < nsigmaPion ) && ( nsigmaKaon < nsigmaProton ) && (nsigmaKaon   < fNSigmaPID))
	{
		if(hman->GetPIDHistogram(kHistPIDTPCKaonRec))
    			hman->GetPIDHistogram(kHistPIDTPCKaonRec)->Fill(trk->GetTPCmomentum(), trk->GetTPCsignal()*trk->Charge()); // Reconstructed PIDKaon histo
    		return kSpKaon;
  	}
  	if( ( nsigmaPion   < nsigmaKaon ) && ( nsigmaPion < nsigmaProton ) && (nsigmaPion   < fNSigmaPID))
	{
		if(hman->GetPIDHistogram(kHistPIDTPCPionRec))
    			hman->GetPIDHistogram(kHistPIDTPCPionRec)->Fill(trk->GetTPCmomentum(), trk->GetTPCsignal()*trk->Charge()); // Reconstructed PIDPion histo
    		return kSpPion;
  	}
  	if( ( nsigmaProton < nsigmaKaon ) && ( nsigmaProton < nsigmaPion ) && (nsigmaProton < fNSigmaPID))
	{
		if(hman->GetPIDHistogram(kHistPIDTPCProtonRec))
    			hman->GetPIDHistogram(kHistPIDTPCProtonRec)->Fill(trk->GetTPCmomentum(), trk->GetTPCsignal()*trk->Charge()); // Reconstructed PIDProton histo
    		return kSpProton;
  	}
  	// else, return undefined
  	return kSpUndefined;

}

Long64_t AliSpectraBothPID::Merge(TCollection* list)
{
  // Merging interface.
  // Returns the number of merged objects (including this).

  Printf("Merging");

  if (!list)
    return 0;

  if (list->IsEmpty())
    return 1;

  TIterator* iter = list->MakeIterator();
  TObject* obj;

  // Actually, we don't do anything here...
  // collections of all histograms
  //  TList collections;

  Int_t count = 0;

  while ((obj = iter->Next())) {
    AliSpectraBothPID* entry = dynamic_cast<AliSpectraBothPID*> (obj);
    if (entry == 0) 
      continue;

    // TH1I * histo = entry->GetHistoCuts();      
    // collections.Add(histo);
    count++;
  }
  
  //  fHistoCuts->Merge(&collections);
  
  delete iter;
  Printf("OK");
  return count+1;
}
Commit	Line	Data
239a080a	1	#include "AliSpectraBothPID.h"
	2	#include "AliAODEvent.h"
	3	#include "TH1F.h"
	4	#include "TH2F.h"
	5	#include "TList.h"
	6	#include "AliAODTrack.h"
	7	#include "AliAODMCParticle.h"
	8	#include "AliPIDResponse.h"
	9	#include "AliAnalysisManager.h"
	10	#include "AliInputEventHandler.h"
	11	#include "AliSpectraBothTrackCuts.h"
	12
	13	ClassImp(AliSpectraBothPID)
	14
2d98dd91	15	AliSpectraBothPID::AliSpectraBothPID() : TNamed("PID", "PID object"), fPIDType(kNSigmacircleTPCTOF), fNSigmaPID(3), fPIDResponse(0),fshiftTPC(0),fshiftTOF(0)
239a080a	16	{
	17
	18	}
	19
	20	AliSpectraBothPID::AliSpectraBothPID(BothPIDType_t pidType) : TNamed("PID", "PID object"), fPIDType(pidType), fNSigmaPID(3), fPIDResponse(0), fshiftTPC(0),fshiftTOF(0)
	21	{
	22
	23
	24
	25	}
	26
	27
	28
	29	void AliSpectraBothPID::FillQAHistos(AliSpectraBothHistoManager * hman, AliVTrack * track, AliSpectraBothTrackCuts * trackCuts)
	30	{
	31
	32	// fill a bunch of QA histos
	33
	34
	35	// Get PID response object, if needed
	36	if(!fPIDResponse)
	37	{
	38	AliAnalysisManager *man = AliAnalysisManager::GetAnalysisManager();
	39	AliInputEventHandler* inputHandler = (AliInputEventHandler*)(man->GetInputEventHandler());
	40	fPIDResponse = inputHandler->GetPIDResponse();
	41	}
b3ea73e1	42	Bool_t ycut[3]={false,false,false};
	43	ycut[0]=trackCuts->CheckYCut ((BothParticleSpecies_t) kSpPion);
	44	ycut[1]=trackCuts->CheckYCut ((BothParticleSpecies_t) kSpKaon);
	45	ycut[2]=trackCuts->CheckYCut ((BothParticleSpecies_t) kSpProton);
	46
	47
239a080a	48	//Response
	49	AliVParticle inEvHMain = dynamic_cast<AliVParticle >(track);
	50
	51	hman->GetPIDHistogram(kHistPIDTPC)->Fill(track->GetTPCmomentum(), track->GetTPCsignal()*track->Charge()); // PID histo
	52	hman->GetPIDHistogram(kHistPIDTPCPion)->Fill(track->GetTPCmomentum(),fPIDResponse->GetTPCResponse().GetExpectedSignal(track->GetTPCmomentum(),AliPID::kPion)*track->Charge()); // Expected PIDPion histo
	53	hman->GetPIDHistogram(kHistPIDTPCKaon)->Fill(track->GetTPCmomentum(),fPIDResponse->GetTPCResponse().GetExpectedSignal(track->GetTPCmomentum(),AliPID::kKaon)*track->Charge()); // Expected PIDKaon histo
	54	hman->GetPIDHistogram(kHistPIDTPCProton)->Fill(track->GetTPCmomentum(),fPIDResponse->GetTPCResponse().GetExpectedSignal(track->GetTPCmomentum(),AliPID::kProton)*track->Charge()); // Expected PIDProton histo
	55
b3ea73e1	56	Double_t nsigmaTPCkProton = fPIDResponse->NumberOfSigmasTPC(inEvHMain, AliPID::kProton)+fshiftTPC;
	57	Double_t nsigmaTPCkKaon = fPIDResponse->NumberOfSigmasTPC(inEvHMain, AliPID::kKaon)+fshiftTPC;
	58	Double_t nsigmaTPCkPion = fPIDResponse->NumberOfSigmasTPC(inEvHMain, AliPID::kPion)+fshiftTPC;
239a080a	59	Double_t nsigmaTOFkProton=100.0,nsigmaTOFkKaon=100.0,nsigmaTOFkPion=100.0;
	60
	61	Double_t nsigmaTPCTOFkProton = TMath::Abs(nsigmaTPCkProton);
	62	Double_t nsigmaTPCTOFkKaon = TMath::Abs(nsigmaTPCkKaon);
	63	Double_t nsigmaTPCTOFkPion = TMath::Abs(nsigmaTPCkPion);
	64	if(fPIDType==kNSigmaTOF)
	65	{
	66	nsigmaTPCTOFkProton = 100.0;
	67	nsigmaTPCTOFkKaon = 100.0;
	68	nsigmaTPCTOFkPion =100.0;
	69
	70	}
0ab8c127	71	Float_t tmppt=TMath::Min(trackCuts->GetPtTOFMatchingPion(),TMath::Min(trackCuts->GetPtTOFMatchingKaon(),trackCuts->GetPtTOFMatchingProton()));
b39b6404	72	if(track->Pt()>=tmppt&&trackCuts->CheckTOFMatchingParticleType(kSpProton)&&trackCuts->CheckTOFMatchingParticleType(kSpKaon)&&trackCuts->CheckTOFMatchingParticleType(kSpPion))
239a080a	73	{
	74
	75	hman->GetPIDHistogram(kHistPIDTOF)->Fill(track->P(),(track->GetTOFsignal()/100)*track->Charge()); // PID histo
	76
b3ea73e1	77	nsigmaTOFkProton = fPIDResponse->NumberOfSigmasTOF(inEvHMain, AliPID::kProton)+fshiftTOF;
	78	nsigmaTOFkKaon = fPIDResponse->NumberOfSigmasTOF(inEvHMain, AliPID::kKaon)+fshiftTOF;
	79	nsigmaTOFkPion = fPIDResponse->NumberOfSigmasTOF(inEvHMain, AliPID::kPion)+fshiftTOF;
239a080a	80
b3ea73e1	81	//TOF
	82	if(ycut[0])
	83	{
	84	hman->GetPtHistogram(kHistNSigPionTOF)->Fill(track->P(),nsigmaTOFkPion );
	85	hman->GetPtHistogram(kHistNSigPionPtTOF)->Fill(track->Pt(),nsigmaTOFkPion );
	86	}
	87	if(ycut[1])
	88	{
	89	hman->GetPtHistogram(kHistNSigKaonTOF)->Fill(track->P(),nsigmaTOFkKaon);
	90	hman->GetPtHistogram(kHistNSigKaonPtTOF)->Fill(track->Pt(),nsigmaTOFkKaon );
	91	}
	92	if(ycut[2])
	93	{
	94	hman->GetPtHistogram(kHistNSigProtonTOF)->Fill(track->P(),nsigmaTOFkProton );
	95	hman->GetPtHistogram(kHistNSigProtonPtTOF)->Fill(track->Pt(),nsigmaTOFkProton );
	96	}
	97
239a080a	98
2d98dd91	99	nsigmaTPCTOFkProton = TMath::Sqrt((nsigmaTPCkProtonnsigmaTPCkProton+nsigmaTOFkProtonnsigmaTOFkProton)/2.0);
	100	nsigmaTPCTOFkKaon = TMath::Sqrt((nsigmaTPCkKaonnsigmaTPCkKaon+nsigmaTOFkKaonnsigmaTOFkKaon)/2.0);
	101	nsigmaTPCTOFkPion = TMath::Sqrt((nsigmaTPCkPionnsigmaTPCkPion+nsigmaTOFkPionnsigmaTOFkPion)/2.0);
239a080a	102
239a080a	103	}
239a080a	104
b3ea73e1	105	if(ycut[0])
	106	{
	107	hman->GetPtHistogram(kHistNSigPionTPC)->Fill(track->P(),nsigmaTPCkPion );
	108	hman->GetPtHistogram(kHistNSigPionPtTPC)->Fill(track->Pt(),nsigmaTPCkPion );
	109	hman->GetPtHistogram(kHistNSigPionTPCTOF)->Fill(track->P(),nsigmaTPCTOFkPion);
	110	hman->GetPtHistogram(kHistNSigPionPtTPCTOF)->Fill(track->Pt(),nsigmaTPCTOFkPion);
	111	}
	112	if(ycut[1])
	113	{
	114	hman->GetPtHistogram(kHistNSigKaonTPC)->Fill(track->P(),nsigmaTPCkKaon );
	115	hman->GetPtHistogram(kHistNSigKaonPtTPC)->Fill(track->Pt(),nsigmaTPCkKaon );
	116	hman->GetPtHistogram(kHistNSigKaonTPCTOF)->Fill(track->P(),nsigmaTPCTOFkKaon);
	117	hman->GetPtHistogram(kHistNSigKaonPtTPCTOF)->Fill(track->Pt(),nsigmaTPCTOFkKaon);
	118
	119	}
	120	if(ycut[2])
	121	{
	122	hman->GetPtHistogram(kHistNSigProtonTPC)->Fill(track->P(),nsigmaTPCkProton );
	123	hman->GetPtHistogram(kHistNSigProtonPtTPC)->Fill(track->Pt(),nsigmaTPCkProton );
	124	hman->GetPtHistogram(kHistNSigProtonTPCTOF)->Fill(track->P(),nsigmaTPCTOFkProton);
	125	hman->GetPtHistogram(kHistNSigProtonPtTPCTOF)->Fill(track->Pt(),nsigmaTPCTOFkProton);
	126	}
239a080a	127
	128	}
	129
	130	Int_t AliSpectraBothPID::GetParticleSpecie(AliAODMCParticle * part)
	131	{
	132	// return PID according to MC truth
	133	switch(TMath::Abs(part->PdgCode()))
	134	{
	135	case 2212:
	136	return kSpProton;
	137	break;
	138	case 321:
	139	return kSpKaon;
	140	break;
	141	case 211:
	142	return kSpPion;
	143	break;
	144	default:
	145	return kSpUndefined;
	146	}
	147	}
	148
	149	Int_t AliSpectraBothPID::GetParticleSpecie(TParticle* part)
	150	{
	151	// return PID according to MC truth
	152	switch(TMath::Abs(part->GetPdgCode()))
	153	{
	154	case 2212:
	155	return kSpProton;
	156	break;
	157	case 321:
	158	return kSpKaon;
	159	break;
	160	case 211:
	161	return kSpPion;
	162	break;
	163	default:
	164	return kSpUndefined;
	165	}
	166	}
	167
	168	Int_t AliSpectraBothPID::GetParticleSpecie(AliSpectraBothHistoManager * hman,AliVTrack * trk, AliSpectraBothTrackCuts * trackCuts, Bool_t* rec)
	169	{
	170	// return PID according to detectors
	171	// Get PID response object, if needed
	172
	173	// guess the particle based on the smaller nsigma
239a080a	174	rec[kSpPion]=false;
	175	rec[kSpKaon]=false;
	176	rec[kSpProton]=false;
	177
	178	if(!fPIDResponse)
	179	{
	180	AliAnalysisManager *man = AliAnalysisManager::GetAnalysisManager();
	181	AliInputEventHandler* inputHandler = (AliInputEventHandler*)(man->GetInputEventHandler());
	182	fPIDResponse = inputHandler->GetPIDResponse();
	183	}
	184
	185	if(!fPIDResponse)
	186	{
	187	AliFatal("Cannot get pid response");
	188	return 0;
	189	}
	190
	191
	192	// Compute nsigma for each hypthesis
	193	AliVParticle inEvHMain = dynamic_cast<AliVParticle >(trk);
	194
	195	// --- TPC
	196	Double_t nsigmaTPCkProton = TMath::Abs(fPIDResponse->NumberOfSigmasTPC(inEvHMain, AliPID::kProton)+fshiftTPC);
	197	Double_t nsigmaTPCkKaon = TMath::Abs(fPIDResponse->NumberOfSigmasTPC(inEvHMain, AliPID::kKaon)+fshiftTPC);
	198	Double_t nsigmaTPCkPion = TMath::Abs(fPIDResponse->NumberOfSigmasTPC(inEvHMain, AliPID::kPion)+fshiftTPC);
	199
	200
	201
0ab8c127	202	Double_t nsigmaTOFkProton=0.0;
	203	Double_t nsigmaTOFkKaon=0.0;
	204	Double_t nsigmaTOFkPion=0.0;
2d98dd91	205
2d98dd91	206
239a080a	207
	208	Double_t nsigmaTPCTOFkProton = TMath::Abs(nsigmaTPCkProton);
	209	Double_t nsigmaTPCTOFkKaon = TMath::Abs(nsigmaTPCkKaon);
	210	Double_t nsigmaTPCTOFkPion = TMath::Abs(nsigmaTPCkPion);
	211	if(fPIDType==kNSigmaTOF)
	212	{
0ab8c127	213	nsigmaTPCTOFkProton = 10000.0;
	214	nsigmaTPCTOFkKaon = 10000.0;
	215	nsigmaTPCTOFkPion =10000.0;
	216	nsigmaTOFkProton = 10000.0;
	217	nsigmaTOFkKaon = 10000.0;
	218	nsigmaTOFkPion =10000.0;
239a080a	219
	220	}
	221
	222
0ab8c127	223	if(trackCuts->GetUseTypeDependedTOFCut())
0ab8c127	224	{
b39b6404	225	if(trackCuts->CheckTOFMatchingParticleType(kSpPion)&&trk->Pt()>=trackCuts->GetPtTOFMatchingPion())
0ab8c127	226	{
	227	nsigmaTOFkPion = TMath::Abs(fPIDResponse->NumberOfSigmasTOF(inEvHMain, AliPID::kPion)+fshiftTOF);
	228	nsigmaTPCTOFkPion = TMath::Sqrt((nsigmaTPCkPionnsigmaTPCkPion+nsigmaTOFkPionnsigmaTOFkPion)/2.);
	229	}
b39b6404	230	else if (trk->Pt()>=trackCuts->GetPtTOFMatchingPion()) // tracks without proper flags
0ab8c127	231	{
	232	nsigmaTOFkPion = 10000.0;
	233	nsigmaTPCTOFkPion =10000.0;
	234
	235	}
	236
b39b6404	237	if(trackCuts->CheckTOFMatchingParticleType(kSpKaon)&&trk->Pt()>=trackCuts->GetPtTOFMatchingKaon())
0ab8c127	238	{
	239	nsigmaTOFkKaon = TMath::Abs(fPIDResponse->NumberOfSigmasTOF(inEvHMain, AliPID::kKaon)+fshiftTOF);
	240	nsigmaTPCTOFkKaon = TMath::Sqrt((nsigmaTPCkKaonnsigmaTPCkKaon+nsigmaTOFkKaonnsigmaTOFkKaon)/2.);
	241	}
b39b6404	242	else if (trk->Pt()>=trackCuts->GetPtTOFMatchingKaon()) // tracks without proper flags
0ab8c127	243	{
	244	nsigmaTOFkKaon = 10000.0;
	245	nsigmaTPCTOFkKaon =10000.0;
239a080a	246
0ab8c127	247	}
0ab8c127	248
b39b6404	249	if(trackCuts->CheckTOFMatchingParticleType(kSpProton)&&trk->Pt()>=trackCuts->GetPtTOFMatchingProton())
0ab8c127	250	{
	251	nsigmaTOFkProton = TMath::Abs(fPIDResponse->NumberOfSigmasTOF(inEvHMain, AliPID::kProton)+fshiftTOF);
	252	nsigmaTPCTOFkProton = TMath::Sqrt((nsigmaTPCkProtonnsigmaTPCkProton+nsigmaTOFkProtonnsigmaTOFkProton)/2.);
	253	}
b39b6404	254	else if (trk->Pt()>=trackCuts->GetPtTOFMatchingProton()) // tracks without proper flags
0ab8c127	255	{
	256	nsigmaTOFkProton = 10000.0;
	257	nsigmaTPCTOFkProton =10000.0;
	258
	259	}
	260
	261
	262	}
b39b6404	263	else if(trk->Pt()>=trackCuts->GetPtTOFMatching())
239a080a	264	{
	265	nsigmaTOFkProton = TMath::Abs(fPIDResponse->NumberOfSigmasTOF(inEvHMain, AliPID::kProton)+fshiftTOF);
	266	nsigmaTOFkKaon = TMath::Abs(fPIDResponse->NumberOfSigmasTOF(inEvHMain, AliPID::kKaon)+fshiftTOF);
	267	nsigmaTOFkPion = TMath::Abs(fPIDResponse->NumberOfSigmasTOF(inEvHMain, AliPID::kPion)+fshiftTOF);
	268	// the TOF info is used in combined
	269	nsigmaTPCTOFkProton = TMath::Sqrt((nsigmaTPCkProtonnsigmaTPCkProton+nsigmaTOFkProtonnsigmaTOFkProton)/2.);
	270	nsigmaTPCTOFkKaon = TMath::Sqrt((nsigmaTPCkKaonnsigmaTPCkKaon+nsigmaTOFkKaonnsigmaTOFkKaon)/2.);
	271	nsigmaTPCTOFkPion = TMath::Sqrt((nsigmaTPCkPionnsigmaTPCkPion+nsigmaTOFkPionnsigmaTOFkPion)/2.);
	272	}
	273	else
	274	{
	275	if (fPIDType==kNSigmaTOF)
	276	return kSpUndefined;
	277
	278	}
	279
	280	// select the nsigma to be used for the actual PID
0ab8c127	281	Double_t nsigmaPion=10000.0;
	282	Double_t nsigmaKaon=10000.0;
	283	Double_t nsigmaProton=10000.0;
239a080a	284
	285	switch (fPIDType)
	286	{
	287	case kNSigmaTPC:
	288	nsigmaProton = nsigmaTPCkProton;
	289	nsigmaKaon = nsigmaTPCkKaon ;
	290	nsigmaPion = nsigmaTPCkPion ;
	291	break;
	292	case kNSigmaTOF:
	293	nsigmaProton = nsigmaTOFkProton;
	294	nsigmaKaon = nsigmaTOFkKaon ;
	295	nsigmaPion = nsigmaTOFkPion ;
	296	break;
2d98dd91	297	case kNSigmacircleTPCTOF:
239a080a	298	nsigmaProton = nsigmaTPCTOFkProton;
	299	nsigmaKaon = nsigmaTPCTOFkKaon ;
	300	nsigmaPion = nsigmaTPCTOFkPion ;
	301	break;
2d98dd91	302	case kNSigmasquareTPCTOF:
	303	nsigmaProton = TMath::Max(nsigmaTPCkProton,nsigmaTOFkProton);
	304	nsigmaKaon = TMath::Max(nsigmaTPCkKaon,nsigmaTOFkKaon);
	305	nsigmaPion = TMath::Max(nsigmaTPCkPion,nsigmaTOFkPion) ;
	306	break;
	307
239a080a	308	}
	309
	310	if(nsigmaPion < fNSigmaPID)
	311	rec[kSpPion]=true;
	312	if(nsigmaKaon < fNSigmaPID)
	313	rec[kSpKaon]=true;
	314	if(nsigmaProton < fNSigmaPID)
	315	rec[kSpProton]=true;
	316
	317	if( ( nsigmaKaon==nsigmaPion ) && ( nsigmaKaon==nsigmaProton ))
	318	return kSpUndefined;//if is the default value for the three
	319	if( ( nsigmaKaon < nsigmaPion ) && ( nsigmaKaon < nsigmaProton ) && (nsigmaKaon < fNSigmaPID))
	320	{
be25efef	321	if(hman->GetPIDHistogram(kHistPIDTPCKaonRec))
be25efef	322	hman->GetPIDHistogram(kHistPIDTPCKaonRec)->Fill(trk->GetTPCmomentum(), trk->GetTPCsignal()*trk->Charge()); // Reconstructed PIDKaon histo
239a080a	323	return kSpKaon;
	324	}
	325	if( ( nsigmaPion < nsigmaKaon ) && ( nsigmaPion < nsigmaProton ) && (nsigmaPion < fNSigmaPID))
	326	{
be25efef	327	if(hman->GetPIDHistogram(kHistPIDTPCPionRec))
be25efef	328	hman->GetPIDHistogram(kHistPIDTPCPionRec)->Fill(trk->GetTPCmomentum(), trk->GetTPCsignal()*trk->Charge()); // Reconstructed PIDPion histo
239a080a	329	return kSpPion;
	330	}
	331	if( ( nsigmaProton < nsigmaKaon ) && ( nsigmaProton < nsigmaPion ) && (nsigmaProton < fNSigmaPID))
	332	{
be25efef	333	if(hman->GetPIDHistogram(kHistPIDTPCProtonRec))
be25efef	334	hman->GetPIDHistogram(kHistPIDTPCProtonRec)->Fill(trk->GetTPCmomentum(), trk->GetTPCsignal()*trk->Charge()); // Reconstructed PIDProton histo
239a080a	335	return kSpProton;
	336	}
	337	// else, return undefined
	338	return kSpUndefined;
	339
	340	}
	341
	342	Long64_t AliSpectraBothPID::Merge(TCollection* list)
	343	{
	344	// Merging interface.
	345	// Returns the number of merged objects (including this).
	346
	347	Printf("Merging");
	348
	349	if (!list)
	350	return 0;
	351
	352	if (list->IsEmpty())
	353	return 1;
	354
	355	TIterator* iter = list->MakeIterator();
	356	TObject* obj;
	357
	358	// Actually, we don't do anything here...
	359	// collections of all histograms
	360	// TList collections;
	361
	362	Int_t count = 0;
	363
	364	while ((obj = iter->Next())) {
	365	AliSpectraBothPID* entry = dynamic_cast<AliSpectraBothPID*> (obj);
	366	if (entry == 0)
	367	continue;
	368
	369	// TH1I * histo = entry->GetHistoCuts();
	370	// collections.Add(histo);
	371	count++;
	372	}
	373
	374	// fHistoCuts->Merge(&collections);
	375
	376	delete iter;
	377	Printf("OK");
	378	return count+1;
	379	}
	380