From e96b9916b66fc7ad958625dfea892bdf248b962a Mon Sep 17 00:00:00 2001 From: morsch Date: Wed, 10 Aug 2011 20:50:05 +0000 Subject: [PATCH] EMCAL PID Updates M. Weber --- STEER/STEERBase/AliEMCALPIDResponse.cxx | 50 +++++++++++- STEER/STEERBase/AliEMCALPIDResponse.h | 12 +-- STEER/STEERBase/AliPIDResponse.cxx | 104 ++++++++++++++++++++++-- STEER/STEERBase/AliPIDResponse.h | 8 ++ STEER/STEERBase/AliVTrack.h | 3 + 5 files changed, 163 insertions(+), 14 deletions(-) diff --git a/STEER/STEERBase/AliEMCALPIDResponse.cxx b/STEER/STEERBase/AliEMCALPIDResponse.cxx index 71745931177..6fce1a7fc81 100644 --- a/STEER/STEERBase/AliEMCALPIDResponse.cxx +++ b/STEER/STEERBase/AliEMCALPIDResponse.cxx @@ -81,6 +81,50 @@ fNorm(0) SetPtBoundary(); SetParametrizations(); } + +AliEMCALPIDResponse::AliEMCALPIDResponse(const AliEMCALPIDResponse &other): + fNorm(other.fNorm) +{ + // + // The copy constructor + // + for(Int_t i = 0; i < fNptBins; i++) + { + fPtCutMin[i] = 0.0; + for(Int_t j = 0; j < 2*AliPID::kSPECIES; j++) + { + fMeanEoP[j][i] = 0.0; + fSigmaEoP[j][i] = 0.0; + fProbLow[j][i] = 0.0; + fProbHigh[j][i] = 0.0; + } + } + + fPtCutMin[fNptBins] = 0.0; + SetPtBoundary(); + SetParametrizations(); +} + + +AliEMCALPIDResponse & AliEMCALPIDResponse::operator=( const AliEMCALPIDResponse& other) +{ + // + // The assignment operator + // + fNorm = other.fNorm; + for(Int_t i = 0; i < fNptBins; i++) + { + fPtCutMin[i] = 0.0; + for(Int_t j = 0; j < 2*AliPID::kSPECIES; j++) + { + fMeanEoP[j][i] = 0.0; + fSigmaEoP[j][i] = 0.0; + fProbLow[j][i] = 0.0; + fProbHigh[j][i] = 0.0; + } + } +} + //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ void AliEMCALPIDResponse::SetPtBoundary(){ // @@ -346,8 +390,8 @@ Double_t AliEMCALPIDResponse::ComputeEMCALProbability(Float_t pt, Float_t eop, I Double_t fRange = 5.0; // hardcoded Double_t nsigma = 0.0; Double_t sigma = 0.0; - Double_t proba = -1.; - + Double_t proba = 999.; + // Check the charge if( charge != -1 && charge != 1){ @@ -358,7 +402,7 @@ Double_t AliEMCALPIDResponse::ComputeEMCALProbability(Float_t pt, Float_t eop, I // default value (will be returned, if pt below threshold) for (Int_t species = 0; species < AliPID::kSPECIES; species++) { - pEMCAL[species] = -1; + pEMCAL[species] = 999.; } if( GetPtBin(pt) > -1 ){ diff --git a/STEER/STEERBase/AliEMCALPIDResponse.h b/STEER/STEERBase/AliEMCALPIDResponse.h index 92aecf74333..e77128e4e7b 100644 --- a/STEER/STEERBase/AliEMCALPIDResponse.h +++ b/STEER/STEERBase/AliEMCALPIDResponse.h @@ -15,13 +15,17 @@ ////////////////////////////////////////////////////////////////////////// #include "AliPID.h" +class TF1; class AliEMCALPIDResponse { public : - AliEMCALPIDResponse(); //ctor - virtual ~AliEMCALPIDResponse() {;} //dtor + AliEMCALPIDResponse(); //ctor + AliEMCALPIDResponse( const AliEMCALPIDResponse& other); //copy ructor + AliEMCALPIDResponse &operator=( const AliEMCALPIDResponse& other); //assignment operator + virtual ~AliEMCALPIDResponse() {;} //dtor + // Getters Int_t GetPtBin(Float_t pt) const; @@ -47,8 +51,6 @@ public : protected: private: - AliEMCALPIDResponse( AliEMCALPIDResponse& r); //dummy copy ructor - AliEMCALPIDResponse &operator=( AliEMCALPIDResponse& r); //dummy assignment operator TF1 *fNorm; // Gauss function for normalizing NON electron probabilities @@ -63,7 +65,7 @@ private: Float_t fProbHigh[2*AliPID::kSPECIES][fNptBins]; // probability above E/p threshold for NON electrons (charge dependent) - ClassDef(AliEMCALPIDResponse,0) + ClassDef(AliEMCALPIDResponse, 1) }; #endif // #ifdef AliEMCALPIDResponse_cxx diff --git a/STEER/STEERBase/AliPIDResponse.cxx b/STEER/STEERBase/AliPIDResponse.cxx index 612e020d005..4b3d4e6c39d 100644 --- a/STEER/STEERBase/AliPIDResponse.cxx +++ b/STEER/STEERBase/AliPIDResponse.cxx @@ -47,6 +47,7 @@ fITSResponse(isMC), fTPCResponse(), fTRDResponse(), fTOFResponse(), +fEMCALResponse(), fRange(5.), fITSPIDmethod(kITSTruncMean), fIsMC(isMC), @@ -65,7 +66,8 @@ fResolutionCorrection(0x0), fTRDPIDParams(0x0), fTRDPIDReference(0x0), fTOFTimeZeroType(kBest_T0), -fTOFres(100.) +fTOFres(100.), +fCurrentEvent(0x0) { // // default ctor @@ -95,6 +97,7 @@ fITSResponse(other.fITSResponse), fTPCResponse(other.fTPCResponse), fTRDResponse(other.fTRDResponse), fTOFResponse(other.fTOFResponse), +fEMCALResponse(other.fEMCALResponse), fRange(other.fRange), fITSPIDmethod(other.fITSPIDmethod), fIsMC(other.fIsMC), @@ -113,7 +116,8 @@ fResolutionCorrection(0x0), fTRDPIDParams(0x0), fTRDPIDReference(0x0), fTOFTimeZeroType(AliPIDResponse::kBest_T0), -fTOFres(100.) +fTOFres(100.), +fCurrentEvent(0x0) { // // copy ctor @@ -134,6 +138,7 @@ AliPIDResponse& AliPIDResponse::operator=(const AliPIDResponse &other) fTPCResponse=other.fTPCResponse; fTRDResponse=other.fTRDResponse; fTOFResponse=other.fTOFResponse; + fEMCALResponse=other.fEMCALResponse; fRange=other.fRange; fITSPIDmethod=other.fITSPIDmethod; fOADBPath=other.fOADBPath; @@ -154,6 +159,7 @@ AliPIDResponse& AliPIDResponse::operator=(const AliPIDResponse &other) memset(fTRDslicesForPID,0,sizeof(UInt_t)*2); fTOFTimeZeroType=AliPIDResponse::kBest_T0; fTOFres=100.; + fCurrentEvent=other.fCurrentEvent; } return *this; } @@ -171,13 +177,55 @@ Float_t AliPIDResponse::NumberOfSigmas(EDetCode detCode, const AliVParticle *tra case kDetTOF: return NumberOfSigmasTOF(track, type); break; // case kDetTRD: return ComputeTRDProbability(track, type); break; // case kDetPHOS: return ComputePHOSProbability(track, type); break; -// case kDetEMCAL: return ComputeEMCALProbability(track, type); break; +// case kDetEMCAL: return NumberOfSigmasEMCAL(track, type); break; // case kDetHMPID: return ComputeHMPIDProbability(track, type); break; default: return -999.; } } +//______________________________________________________________________________ +Float_t AliPIDResponse::NumberOfSigmasEMCAL(const AliVTrack *track, AliPID::EParticleType type) const { + + AliVCluster *matchedClus = NULL; + + Double_t mom = -1.; + Double_t pt = -1.; + Double_t EovP = -1.; + Double_t fClsE = -1.; + + Int_t nMatchClus = -1; + Int_t charge = 0; + + // Track matching + nMatchClus = track->GetEMCALcluster(); + if(nMatchClus > -1){ + + mom = track->P(); + pt = track->Pt(); + charge = track->Charge(); + + matchedClus = (AliVCluster*)fCurrentEvent->GetCaloCluster(nMatchClus); + + if(matchedClus){ + + // matched cluster is EMCAL + if(matchedClus->IsEMCAL()){ + + fClsE = matchedClus->E(); + EovP = fClsE/mom; + + + // NSigma value really meaningful only for electrons! + return fEMCALResponse.GetNumberOfSigmas(pt,EovP,type,charge); + } + } + } + + return -999; + +} + //______________________________________________________________________________ AliPIDResponse::EDetPidStatus AliPIDResponse::ComputePIDProbability (EDetCode detCode, const AliVTrack *track, Int_t nSpecies, Double_t p[]) const { @@ -366,15 +414,57 @@ AliPIDResponse::EDetPidStatus AliPIDResponse::ComputeTRDProbability (const AliV return kDetPidOk; } //______________________________________________________________________________ -AliPIDResponse::EDetPidStatus AliPIDResponse::ComputeEMCALProbability(const AliVTrack */*track*/, Int_t nSpecies, Double_t p[]) const +AliPIDResponse::EDetPidStatus AliPIDResponse::ComputeEMCALProbability (const AliVTrack *track, Int_t nSpecies, Double_t p[]) const { // // Compute PID response for the EMCAL // - // set flat distribution (no decision) + AliVCluster *matchedClus = NULL; + + Double_t mom = -1.; + Double_t pt = -1.; + Double_t EovP = -1.; + Double_t fClsE = -1.; + + Int_t nMatchClus = -1; + Int_t charge = 0; + + // Track matching + nMatchClus = track->GetEMCALcluster(); + + if(nMatchClus > -1){ + + mom = track->P(); + pt = track->Pt(); + charge = track->Charge(); + + matchedClus = (AliVCluster*)fCurrentEvent->GetCaloCluster(nMatchClus); + + if(matchedClus){ + + // matched cluster is EMCAL + if(matchedClus->IsEMCAL()){ + + fClsE = matchedClus->E(); + EovP = fClsE/mom; + + + // compute the probabilities + if( 999 != fEMCALResponse.ComputeEMCALProbability(pt,EovP,charge,p)){ + + // in case everything is OK + return kDetPidOk; + + } + } + } + } + + // in all other cases set flat distribution (no decision) for (Int_t j=0; jGetRunNumber(); if (fRun!=fOldRun){ diff --git a/STEER/STEERBase/AliPIDResponse.h b/STEER/STEERBase/AliPIDResponse.h index d2ba62e5f32..e34ae654427 100644 --- a/STEER/STEERBase/AliPIDResponse.h +++ b/STEER/STEERBase/AliPIDResponse.h @@ -15,6 +15,7 @@ #include "AliTPCPIDResponse.h" #include "AliTRDPIDResponse.h" #include "AliTOFPIDResponse.h" +#include "AliEMCALPIDResponse.h" #include "AliVParticle.h" #include "AliVTrack.h" @@ -53,11 +54,13 @@ public: AliTPCPIDResponse &GetTPCResponse() {return fTPCResponse;} AliTOFPIDResponse &GetTOFResponse() {return fTOFResponse;} AliTRDPIDResponse &GetTRDResponse() {return fTRDResponse;} + AliEMCALPIDResponse &GetEMCALResponse() {return fEMCALResponse;} Float_t NumberOfSigmas(EDetCode detCode, const AliVParticle *track, AliPID::EParticleType type) const; virtual Float_t NumberOfSigmasITS(const AliVParticle *track, AliPID::EParticleType type) const; virtual Float_t NumberOfSigmasTPC(const AliVParticle *track, AliPID::EParticleType type) const; + virtual Float_t NumberOfSigmasEMCAL(const AliVTrack *track, AliPID::EParticleType type) const; virtual Float_t NumberOfSigmasTOF(const AliVParticle *track, AliPID::EParticleType type) const = 0; virtual Bool_t IdentifiedAsElectronTRD(const AliVTrack *track, Double_t efficiencyLevel) const; @@ -80,6 +83,8 @@ public: void InitialiseEvent(AliVEvent *event, Int_t pass); void SetCurrentFile(const char* file) { fCurrentFile=file; } + AliVEvent * GetCurrentEvent() const {return fCurrentEvent;} + // User settings for the MC period and reco pass void SetMCperiod(const char *mcPeriod) {fMCperiodUser=mcPeriod;} void SetRecoPass(Int_t recoPass) {fRecoPassUser=recoPass;} @@ -92,6 +97,7 @@ protected: AliTPCPIDResponse fTPCResponse; //PID response function of the TPC AliTRDPIDResponse fTRDResponse; //PID response function of the TRD AliTOFPIDResponse fTOFResponse; //PID response function of the TOF + AliEMCALPIDResponse fEMCALResponse; //PID response function of the EMCAL Float_t fRange; // nSigma max in likelihood ITSPIDmethod fITSPIDmethod; // 0 = trunc mean; 1 = likelihood @@ -120,6 +126,8 @@ private: Int_t fTOFTimeZeroType; //! default start time type for tof (ESD) Float_t fTOFres; //! TOF resolution + + AliVEvent *fCurrentEvent; //! event currently being processed void ExecNewRun(); diff --git a/STEER/STEERBase/AliVTrack.h b/STEER/STEERBase/AliVTrack.h index ad597aaebeb..f565043853c 100644 --- a/STEER/STEERBase/AliVTrack.h +++ b/STEER/STEERBase/AliVTrack.h @@ -55,6 +55,9 @@ public: virtual UShort_t GetTPCNclsF() const { return 0;} virtual Double_t GetTRDslice(Int_t /*plane*/, Int_t /*slice*/) const { return -1.; } + virtual Int_t GetEMCALcluster() const {return -1;} + virtual Int_t GetPHOScluster() const {return -1;} + //pid info virtual Double_t GetITSsignal() const {return 0.;} virtual Double_t GetTPCsignal() const {return 0.;} -- 2.43.0