X-Git-Url: http://git.uio.no/git/?a=blobdiff_plain;f=PWG4%2FPartCorrBase%2FAliCaloPID.cxx;h=46915a006f0905874305210053765736d925dd0b;hb=cb41788f6b3084b23fa5616bf78c0ce4ebae9b63;hp=0a5eb46f5f177ba466dcfcb44b30effcdb3de2bd;hpb=477d6cee166e6352a42084ae97f00e0a425eead6;p=u%2Fmrichter%2FAliRoot.git diff --git a/PWG4/PartCorrBase/AliCaloPID.cxx b/PWG4/PartCorrBase/AliCaloPID.cxx index 0a5eb46f5f1..46915a006f0 100755 --- a/PWG4/PartCorrBase/AliCaloPID.cxx +++ b/PWG4/PartCorrBase/AliCaloPID.cxx @@ -15,9 +15,26 @@ /* $Id: AliCaloPID.cxx 21839 2007-10-29 13:49:42Z gustavo $ */ //_________________________________________________________________________ -// Class for track/cluster acceptance selection -// Selection in Central barrel, EMCAL and PHOS -// +// Class for PID selection with calorimeters +// The Output of the 2 main methods GetPdg is a PDG number identifying the cluster, +// being kPhoton, kElectron, kPi0 ... as defined in the header file +// - GetPdg(const TString calo, const Double_t * pid, const Float_t energy) +// Reads the PID weights array of the ESDs and depending on its magnitude identifies the particle +// - GetPdg(const TString calo,const TLorentzVector mom, const AliAODCaloCluster * cluster) +// Recalcultes PID, the bayesian or any new one to be implemented in the future +// Right now only the possibility to recalculate EMCAL with bayesian and simple PID. +// In order to recalculate Bayesian, it is necessary to load the EMCALUtils library +// and do SwitchOnBayesianRecalculation(). +// To change the PID parameters from Low to High like the ones by default, use the constructor +// AliCaloPID(flux) +// where flux is AliCaloPID::kLow or AliCaloPID::kHigh +// If it is necessary to change the parameters use the constructor +// AliCaloPID(AliEMCALPIDUtils *utils) and set the parameters before. +// - SetPIDBits: Simple PID, depending on the thresholds fDispCut fTOFCut and even the +// result of the PID bayesian a different PID bit is set. +// +// All these methods can be called in the analysis you are interested. +// //*-- Author: Gustavo Conesa (LNF-INFN) ////////////////////////////////////////////////////////////////////////////// @@ -31,7 +48,7 @@ #include "AliCaloPID.h" #include "AliAODCaloCluster.h" #include "AliAODPWG4Particle.h" - +#include "AliEMCALPIDUtils.h" ClassImp(AliCaloPID) @@ -44,7 +61,8 @@ fPHOSPhotonWeight(0.), fPHOSPi0Weight(0.), fPHOSElectronWeight(0.), fPHOSChargeWeight(0.) , fPHOSNeutralWeight(0.), fPHOSWeightFormula(0), fPHOSPhotonWeightFormula(0x0), fPHOSPi0WeightFormula(0x0), -fDispCut(0.),fTOFCut(0.), fDebug(-1) +fDispCut(0.),fTOFCut(0.), fDebug(-1), +fRecalculateBayesian(kFALSE), fParticleFlux(kLow), fEMCALPIDUtils(new AliEMCALPIDUtils) { //Ctor @@ -52,6 +70,42 @@ fDispCut(0.),fTOFCut(0.), fDebug(-1) InitParameters(); } +//________________________________________________ +AliCaloPID::AliCaloPID(const Int_t flux) : +TObject(), fEMCALPhotonWeight(0.), fEMCALPi0Weight(0.), +fEMCALElectronWeight(0.), fEMCALChargeWeight(0.), +fEMCALNeutralWeight(0.), +fPHOSPhotonWeight(0.), fPHOSPi0Weight(0.), +fPHOSElectronWeight(0.), fPHOSChargeWeight(0.) , +fPHOSNeutralWeight(0.), fPHOSWeightFormula(0), +fPHOSPhotonWeightFormula(0x0), fPHOSPi0WeightFormula(0x0), +fDispCut(0.),fTOFCut(0.), fDebug(-1), +fRecalculateBayesian(kTRUE), fParticleFlux(flux), fEMCALPIDUtils(new AliEMCALPIDUtils) +{ + //Ctor + + //Initialize parameters + InitParameters(); +} + +//________________________________________________ +AliCaloPID::AliCaloPID(const TTask * emcalpid) : +TObject(), fEMCALPhotonWeight(0.), fEMCALPi0Weight(0.), +fEMCALElectronWeight(0.), fEMCALChargeWeight(0.), +fEMCALNeutralWeight(0.), +fPHOSPhotonWeight(0.), fPHOSPi0Weight(0.), +fPHOSElectronWeight(0.), fPHOSChargeWeight(0.) , +fPHOSNeutralWeight(0.), fPHOSWeightFormula(0), +fPHOSPhotonWeightFormula(0x0), fPHOSPi0WeightFormula(0x0), +fDispCut(0.),fTOFCut(0.), fDebug(-1), +fRecalculateBayesian(kTRUE), fParticleFlux(-1), fEMCALPIDUtils( (AliEMCALPIDUtils*) emcalpid) +{ + //Ctor + + //Initialize parameters + InitParameters(); +} + //____________________________________________________________________________ AliCaloPID::AliCaloPID(const AliCaloPID & pid) : TObject(pid), fEMCALPhotonWeight(pid.fEMCALPhotonWeight), @@ -68,7 +122,8 @@ fPHOSWeightFormula(pid.fPHOSWeightFormula), fPHOSPhotonWeightFormula(pid.fPHOSPhotonWeightFormula), fPHOSPi0WeightFormula(pid.fPHOSPi0WeightFormula), fDispCut(pid.fDispCut),fTOFCut(pid.fTOFCut), -fDebug(pid.fDebug) +fDebug(pid.fDebug), fRecalculateBayesian(pid.fRecalculateBayesian), +fParticleFlux(pid.fParticleFlux), fEMCALPIDUtils(pid.fEMCALPIDUtils) { // cpy ctor @@ -101,6 +156,10 @@ AliCaloPID & AliCaloPID::operator = (const AliCaloPID & pid) fTOFCut = pid.fTOFCut; fDebug = pid.fDebug; + fRecalculateBayesian = pid.fRecalculateBayesian; + fParticleFlux = pid.fParticleFlux; + fEMCALPIDUtils = pid.fEMCALPIDUtils; + return *this; } @@ -111,7 +170,7 @@ AliCaloPID::~AliCaloPID() { if(fPHOSPhotonWeightFormula) delete fPHOSPhotonWeightFormula ; if(fPHOSPi0WeightFormula) delete fPHOSPi0WeightFormula ; - + if(fEMCALPIDUtils) delete fEMCALPIDUtils ; } @@ -120,14 +179,14 @@ void AliCaloPID::InitParameters() { //Initialize the parameters of the PID. - fEMCALPhotonWeight = 0.8 ; + fEMCALPhotonWeight = 0.5 ; fEMCALPi0Weight = 0.5 ; - fEMCALElectronWeight = 0.8 ; + fEMCALElectronWeight = 0.5 ; fEMCALChargeWeight = 0.5 ; fEMCALNeutralWeight = 0.5 ; - fPHOSPhotonWeight = 0.75 ; - fPHOSPi0Weight = 0.8 ; + fPHOSPhotonWeight = 0.5 ; + fPHOSPi0Weight = 0.5 ; fPHOSElectronWeight = 0.5 ; fPHOSChargeWeight = 0.5 ; fPHOSNeutralWeight = 0.5 ; @@ -142,6 +201,17 @@ void AliCaloPID::InitParameters() fDispCut = 1.5; fTOFCut = 5.e-9; fDebug = -1; + + if(fRecalculateBayesian){ + if(fParticleFlux == kLow){ + printf("AliCaloPID::Init() - SetLOWFluxParam\n"); + fEMCALPIDUtils->SetLowFluxParam() ; + } + else if (fParticleFlux == kHigh){ + printf("AliCaloPID::Init() - SetHIGHFluxParam\n"); + fEMCALPIDUtils->SetHighFluxParam() ; + } + } } //_______________________________________________________________ @@ -153,11 +223,11 @@ Int_t AliCaloPID::GetPdg(const TString calo, const Double_t * pid, const Float_t abort(); } - Float_t wPh = fPHOSPhotonWeight ; + Float_t wPh = fPHOSPhotonWeight ; Float_t wPi0 = fPHOSPi0Weight ; - Float_t wE = fPHOSElectronWeight ; - Float_t wCh = fPHOSChargeWeight ; - Float_t wNe = fPHOSNeutralWeight ; + Float_t wE = fPHOSElectronWeight ; + Float_t wCh = fPHOSChargeWeight ; + Float_t wNe = fPHOSNeutralWeight ; if(calo == "PHOS" && fPHOSWeightFormula){ @@ -202,9 +272,11 @@ Int_t AliCaloPID::GetPdg(const TString calo, const Double_t * pid, const Float_t pdg = kNeutralUnknown ; } else{//EMCAL - if(pid[AliAODCluster::kPhoton] > wPh) pdg = kPhoton ; + + if(pid[AliAODCluster::kPhoton]+pid[AliAODCluster::kElectron] > wPh) pdg = kPhoton ; //temporal sollution until track matching for electrons is considered + //if(pid[AliAODCluster::kPhoton] > wPh) pdg = kPhoton ; else if(pid[AliAODCluster::kPi0] > wPi0) pdg = kPi0 ; - else if(pid[AliAODCluster::kElectron] > wE) pdg = kElectron ; + //else if(pid[AliAODCluster::kElectron] > wE) pdg = kElectron ; else if(chargedHadronWeight + neutralHadronWeight > wCh) pdg = kChargedHadron ; else if(neutralHadronWeight + chargedHadronWeight > wNe) pdg = kNeutralHadron ; else pdg = kNeutralUnknown ; @@ -220,13 +292,37 @@ Int_t AliCaloPID::GetPdg(const TString calo, const Double_t * pid, const Float_t //_______________________________________________________________ Int_t AliCaloPID::GetPdg(const TString calo,const TLorentzVector mom, const AliAODCaloCluster * cluster) const { //Recalculated PID with all parameters - if(fDebug > 0)printf("AliCaloPID::GetPdg: Calorimeter %s, E %3.2f, l0 %3.2f, l1 %3.2f, disp %3.2f, tof %1.11f, distCPV %3.2f, distToBC %1.1f, NMax %d\n", - calo.Data(),mom.E(),cluster->GetM02(),cluster->GetM20(),cluster->GetDispersion(),cluster->GetTOF(), - cluster->GetEmcCpvDistance(), cluster->GetDistToBadChannel(),cluster->GetNExMax()); - + Float_t lambda0 = cluster->GetM02(); + Float_t energy = mom.E(); + + if(fDebug > 0) printf("AliCaloPID::GetPdg: Calorimeter %s, E %3.2f, l0 %3.2f, l1 %3.2f, disp %3.2f, tof %1.11f, distCPV %3.2f, distToBC %1.1f, NMax %d\n", + calo.Data(),energy,lambda0,cluster->GetM20(),cluster->GetDispersion(),cluster->GetTOF(), + cluster->GetEmcCpvDistance(), cluster->GetDistToBadChannel(),cluster->GetNExMax()); + if(calo == "EMCAL") { - if(cluster->GetM02()< 0.25) return kPhoton ; - else return kNeutralHadron ; + //Recalculate Bayesian + if(fRecalculateBayesian){ + if(fDebug > 0) { + const Double_t *pid0 = cluster->PID(); + printf("AliCaloPID::GetPdg: BEFORE calo %s, ph %0.2f, pi0 %0.2f, el %0.2f, conv el %0.2f, hadrons: pion %0.2f, kaon %0.2f, proton %0.2f , neutron %0.2f, kaon %0.2f \n", + calo.Data(),pid0[AliAODCluster::kPhoton], pid0[AliAODCluster::kPi0], + pid0[AliAODCluster::kElectron], pid0[AliAODCluster::kEleCon], + pid0[AliAODCluster::kPion], pid0[AliAODCluster::kKaon], pid0[AliAODCluster::kProton], + pid0[AliAODCluster::kNeutron], pid0[AliAODCluster::kKaon0]); + } + + fEMCALPIDUtils->ComputePID(energy, lambda0); + Double_t pid[AliPID::kSPECIESN]; + for(Int_t i = 0; i < AliPID::kSPECIESN; i++) pid[i] = fEMCALPIDUtils->GetPIDFinal(i); + return GetPdg(calo, pid, energy); + + + } + + // If no use of bayesian, simple PID + if(lambda0 < 0.25) return kPhoton ; + //else return kNeutralHadron ; + else return kPi0 ; } // if(calo == "PHOS") { @@ -272,8 +368,10 @@ TString AliCaloPID::GetPIDParametersList() { parList+=onePar ; if(fPHOSWeightFormula){ - parList+="PHOS Photon Weight Formula: "+(fPHOSPhotonWeightFormula->GetExpFormula("p")); - parList+="PHOS Pi0 Weight Formula: "+(fPHOSPi0WeightFormula->GetExpFormula("p")); + sprintf(onePar,"PHOS Photon Weight Formula: %s\n",(fPHOSPhotonWeightFormula->GetExpFormula("p")).Data()) ; + parList+=onePar; + sprintf(onePar,"PHOS Pi0 Weight Formula: %s\n",(fPHOSPi0WeightFormula->GetExpFormula("p")).Data()) ; + parList+=onePar; } return parList; @@ -306,7 +404,8 @@ void AliCaloPID::Print(const Option_t * opt) const printf("TOF cut = %e\n",fTOFCut); printf("Dispersion cut = %2.2f\n",fDispCut); printf("Debug level = %d\n",fDebug); - + printf("Recalculate Bayesian? = %d\n",fRecalculateBayesian); + if(fRecalculateBayesian) printf("Particle Flux? = %d\n",fParticleFlux); printf(" \n"); }