AliAODPid.cxx.diff make getter const
authormorsch <morsch@f7af4fe6-9843-0410-8265-dc069ae4e863>
Tue, 21 Jun 2011 19:53:45 +0000 (19:53 +0000)
committermorsch <morsch@f7af4fe6-9843-0410-8265-dc069ae4e863>
Tue, 21 Jun 2011 19:53:45 +0000 (19:53 +0000)
AliAODPid.h.diff             make getter const
AliAODpidUtil.cxx.diff       clean up TRD pid calculation
AliAODpidUtil.h.diff         clean up TRD pid calculation
AliAODTrack.h.diff           Add interface for HMPIDpid
AliESDpid.cxx.diff           clean up TRD pid calculation
AliPIDResponse.cxx.diff      add TRD pid calculation, add HMPID probs
AliPIDResponse.h.diff        add TRD pid calculation, add HMPID probs
AliTRDPIDResponse.cxx.diff   updates from Markus Fasel
AliTRDPIDResponse.h.diff     updates from Markus Fasel
AliVTrack.h.diff             Add interface for HMPIDpid
J. Wiechula

STEER/AliAODPid.cxx
STEER/AliAODPid.h
STEER/AliAODTrack.h
STEER/AliAODpidUtil.cxx
STEER/AliAODpidUtil.h
STEER/AliESDpid.cxx
STEER/AliPIDResponse.cxx
STEER/AliPIDResponse.h
STEER/AliTRDPIDResponse.cxx
STEER/AliTRDPIDResponse.h
STEER/AliVTrack.h

index 1f46271..5070cb8 100644 (file)
@@ -172,7 +172,7 @@ void AliAODPid::SetHMPIDprobs(Double_t hmpPid[5])
   for(Int_t i = 0; i < 5; i++ ) fHMPIDprobs[i] =  hmpPid[i];
 }
 //______________________________________________________________________________
-void AliAODPid::GetHMPIDprobs(Double_t *p) 
+void AliAODPid::GetHMPIDprobs(Double_t *p) const
 {
   //
   // Set the HMPID PID probablities that are read from ESD
index 3b86cb5..2465c3f 100644 (file)
@@ -46,7 +46,7 @@ class AliAODPid : public TObject {
   Float_t*  GetTRDmomentum()           {return  fTRDmomentum;}
   Double_t  GetTOFsignal()       const {return  fTOFesdsignal;}
   Double_t  GetHMPIDsignal()     const {return  fHMPIDsignal;}
-  void  GetHMPIDprobs(Double_t *p);            
+  void      GetHMPIDprobs(Double_t *p) const;
 
   void      GetIntegratedTimes(Double_t timeint[5])  const; 
   void      GetEMCALPosition  (Double_t emcalpos[3]) const;
index 262914a..9cf8a98 100644 (file)
@@ -231,9 +231,10 @@ class AliAODTrack : public AliVTrack {
   UShort_t  GetTPCsignalN()      const { return fDetPid?fDetPid->GetTPCsignalN():0;    }
   Double_t  GetTPCmomentum()     const { return fDetPid?fDetPid->GetTPCmomentum():0.;  }
   Double_t  GetTOFsignal()       const { return fDetPid?fDetPid->GetTOFsignal():0.; }
-  void GetIntegratedTimes(Double_t *times) const {if (fDetPid) fDetPid->GetIntegratedTimes(times); }
+  void      GetIntegratedTimes(Double_t *times) const {if (fDetPid) fDetPid->GetIntegratedTimes(times); }
   Double_t  GetTRDslice(Int_t plane, Int_t slice) const;
-  Double_t GetTRDmomentum(Int_t plane, Double_t */*sp*/=0x0) const;
+  Double_t  GetTRDmomentum(Int_t plane, Double_t */*sp*/=0x0) const;
+  void      GetHMPIDpid(Double_t *p) const { if (fDetPid) fDetPid->GetHMPIDprobs(p); }
 
   
   AliAODPid    *GetDetPid() const { return fDetPid; }
index 0dddc91..a005b67 100644 (file)
@@ -223,20 +223,7 @@ void AliAODpidUtil::MakeTOFPID(const AliAODTrack *track, Double_t *p) const
 //_________________________________________________________________________
 void AliAODpidUtil::MakeTRDPID(const AliAODTrack *track,Double_t *p) const
 {
-  
-  // Method to recalculate the TRD PID probabilities
-  if ((track->GetStatus()&AliESDtrack::kTRDout )==0)   return;
-
-  AliAODPid *pidObj = track->GetDetPid();
-  Float_t *mom=pidObj->GetTRDmomentum();
-  Int_t ntracklets=0;
-  for(Int_t iPl=0;iPl<6;iPl++){
-   if(mom[iPl]>0.) ntracklets++;
-  }
-   if(ntracklets<4) return;
-
-  Double_t* dedx=pidObj->GetTRDsignal();
-  Bool_t norm=kTRUE;
-  fTRDResponse.GetResponse(pidObj->GetTRDnSlices(),dedx,mom,p,norm);
+  ComputeTRDProbability(track, AliPID::kSPECIES, p); 
   return;
 }
+
index 0ce178f..451ed6c 100644 (file)
@@ -25,7 +25,7 @@ class AliVParticle;
 class AliAODpidUtil : public AliPIDResponse  {
 public:
   //TODO: isMC???
-  AliAODpidUtil(Bool_t isMC = kFALSE): AliPIDResponse(isMC), fRange(5.) {;}
+  AliAODpidUtil(Bool_t isMC = kFALSE): AliPIDResponse(isMC) {;}
   virtual ~AliAODpidUtil() {;}
 
   Int_t MakePID(const AliAODTrack *track,Double_t *p) const;
@@ -38,9 +38,8 @@ public:
   virtual Float_t NumberOfSigmasTOF(const AliVParticle *vtrack, AliPID::EParticleType type) const;
   
 private:
-  Float_t           fRange;          // nSigma max in likelihood
   
-  ClassDef(AliAODpidUtil,2)  // PID calculation class
+  ClassDef(AliAODpidUtil,3)  // PID calculation class
 };
 
 inline Float_t AliAODpidUtil::NumberOfSigmasTOF(const AliVParticle *vtrack, AliPID::EParticleType type) const {
index 220c211..589d963 100644 (file)
@@ -246,16 +246,8 @@ void AliESDpid::MakeTRDPID(AliESDtrack *track) const
   //
   // Method to recalculate the TRD PID probabilities
   //
-  if((track->GetStatus()&AliESDtrack::kTRDout)==0) return;
-  Double_t prob[AliPID::kSPECIES]; Float_t mom[6];
-  Double_t dedx[48];  // Allocate space for the maximum number of TRD slices
-  for(Int_t ilayer = 0; ilayer < 6; ilayer++){
-    mom[ilayer] = track->GetTRDmomentum(ilayer);
-    for(Int_t islice = 0; islice < track->GetNumberOfTRDslices(); islice++){
-      dedx[ilayer*track->GetNumberOfTRDslices()+islice] = track->GetTRDslice(ilayer, islice);
-    }
-  }
-  fTRDResponse.GetResponse(track->GetNumberOfTRDslices(), dedx, mom, prob);
+  Double_t prob[AliPID::kSPECIES];
+  ComputeTRDProbability(track, AliPID::kSPECIES, prob);
   track->SetTRDpid(prob);
 }
 //_________________________________________________________________________
@@ -302,6 +294,9 @@ void AliESDpid::CombinePID(AliESDtrack *track) const
 }
 //_________________________________________________________________________
 Bool_t AliESDpid::CheckTOFMatching(AliESDtrack *track) const{
+  //
+  // Check pid matching of TOF with TPC as reference
+  //
     Bool_t status = kFALSE;
     
     Double_t exptimes[5];
index 3e068cd..f64b012 100644 (file)
@@ -34,6 +34,8 @@
 #include <AliVTrack.h>
 #include <AliLog.h>
 #include <AliPID.h>
+#include <AliOADBContainer.h>
+#include <AliTRDPIDReference.h>
 
 #include "AliPIDResponse.h"
 
@@ -53,12 +55,15 @@ fBeamType("PP"),
 fLHCperiod(),
 fMCperiodTPC(),
 fMCperiodUser(),
+fCurrentFile(),
 fRecoPass(0),
 fRecoPassUser(-1),
 fRun(0),
 fOldRun(0),
 fArrPidResponseMaster(0x0),
 fResolutionCorrection(0x0),
+fTRDPIDParams(0x0),
+fTRDPIDReference(0x0),
 fTOFTimeZeroType(kBest_T0),
 fTOFres(100.)
 {
@@ -68,6 +73,8 @@ fTOFres(100.)
   AliLog::SetClassDebugLevel("AliPIDResponse",10);
   AliLog::SetClassDebugLevel("AliESDpid",10);
   AliLog::SetClassDebugLevel("AliAODpidUtil",10);
+
+  memset(fTRDslicesForPID,0,sizeof(UInt_t)*2);
 }
 
 //______________________________________________________________________________
@@ -77,6 +84,8 @@ AliPIDResponse::~AliPIDResponse()
   // dtor
   //
   delete fArrPidResponseMaster;
+  delete fTRDPIDParams;
+  delete fTRDPIDReference;
 }
 
 //______________________________________________________________________________
@@ -94,18 +103,22 @@ fBeamType("PP"),
 fLHCperiod(),
 fMCperiodTPC(),
 fMCperiodUser(other.fMCperiodUser),
+fCurrentFile(),
 fRecoPass(0),
 fRecoPassUser(other.fRecoPassUser),
 fRun(0),
 fOldRun(0),
 fArrPidResponseMaster(0x0),
 fResolutionCorrection(0x0),
+fTRDPIDParams(0x0),
+fTRDPIDReference(0x0),
 fTOFTimeZeroType(AliPIDResponse::kBest_T0),
 fTOFres(100.)
 {
   //
   // copy ctor
   //
+  memset(fTRDslicesForPID,0,sizeof(UInt_t)*2);
 }
 
 //______________________________________________________________________________
@@ -129,12 +142,16 @@ AliPIDResponse& AliPIDResponse::operator=(const AliPIDResponse &other)
     fLHCperiod="";
     fMCperiodTPC="";
     fMCperiodUser=other.fMCperiodUser;
+    fCurrentFile="";
     fRecoPass=0;
     fRecoPassUser=other.fRecoPassUser;
     fRun=0;
     fOldRun=0;
     fArrPidResponseMaster=0x0;
     fResolutionCorrection=0x0;
+    fTRDPIDParams=0x0;
+    fTRDPIDReference=0x0;
+    memset(fTRDslicesForPID,0,sizeof(UInt_t)*2);
     fTOFTimeZeroType=AliPIDResponse::kBest_T0;
     fTOFres=100.;
   }
@@ -325,7 +342,7 @@ AliPIDResponse::EDetPidStatus AliPIDResponse::ComputeTOFProbability  (const AliV
   return kDetPidOk;
 }
 //______________________________________________________________________________
-AliPIDResponse::EDetPidStatus AliPIDResponse::ComputeTRDProbability  (const AliVTrack */*track*/, Int_t nSpecies, Double_t p[]) const
+AliPIDResponse::EDetPidStatus AliPIDResponse::ComputeTRDProbability  (const AliVTrack *track, Int_t nSpecies, Double_t p[]) const
 {
   //
   // Compute PID response for the
@@ -333,7 +350,20 @@ AliPIDResponse::EDetPidStatus AliPIDResponse::ComputeTRDProbability  (const AliV
 
   // set flat distribution (no decision)
   for (Int_t j=0; j<nSpecies; j++) p[j]=1./nSpecies;
-  return kDetNoSignal;
+  if((track->GetStatus()&AliVTrack::kTRDout)==0) return kDetNoSignal;
+
+  Float_t mom[6];
+  Double_t dedx[48];  // Allocate space for the maximum number of TRD slices
+  Int_t nslices = fTRDslicesForPID[1] - fTRDslicesForPID[0] + 1;
+  AliDebug(1, Form("First Slice: %d, Last Slice: %d, Number of slices: %d",  fTRDslicesForPID[0], fTRDslicesForPID[1], nslices));
+  for(UInt_t ilayer = 0; ilayer < 6; ilayer++){
+    mom[ilayer] = track->GetTRDmomentum(ilayer);
+    for(UInt_t islice = fTRDslicesForPID[0]; islice <= fTRDslicesForPID[1]; islice++){
+      dedx[ilayer*nslices+islice-fTRDslicesForPID[0]] = track->GetTRDslice(ilayer, islice);
+    }
+  }
+  fTRDResponse.GetResponse(nslices, dedx, mom, p);
+  return kDetPidOk;
 }
 //______________________________________________________________________________
 AliPIDResponse::EDetPidStatus AliPIDResponse::ComputeEMCALProbability(const AliVTrack */*track*/, Int_t nSpecies, Double_t p[]) const
@@ -358,7 +388,7 @@ AliPIDResponse::EDetPidStatus AliPIDResponse::ComputePHOSProbability (const AliV
   return kDetNoSignal;
 }
 //______________________________________________________________________________
-AliPIDResponse::EDetPidStatus AliPIDResponse::ComputeHMPIDProbability(const AliVTrack */*track*/, Int_t nSpecies, Double_t p[]) const
+AliPIDResponse::EDetPidStatus AliPIDResponse::ComputeHMPIDProbability(const AliVTrack *track, Int_t nSpecies, Double_t p[]) const
 {
   //
   // Compute PID response for the HMPID
@@ -366,7 +396,11 @@ AliPIDResponse::EDetPidStatus AliPIDResponse::ComputeHMPIDProbability(const AliV
 
   // set flat distribution (no decision)
   for (Int_t j=0; j<nSpecies; j++) p[j]=1./nSpecies;
-  return kDetNoSignal;
+  if((track->GetStatus()&AliVTrack::kHMPIDpid)==0) return kDetNoSignal;
+
+  track->GetHMPIDpid(p);
+  
+  return kDetPidOk;
 }
 
 //______________________________________________________________________________
@@ -445,6 +479,7 @@ void AliPIDResponse::SetRecoInfo()
     
   fBeamType="PP";
   
+  TPRegexp reg(".*(LHC11[a-z]+[0-9]+[a-z]*)/.*");
   //find the period by run number (UGLY, but not stored in ESD and AOD... )
   if (fRun>=114737&&fRun<=117223)      { fLHCperiod="LHC10B"; fMCperiodTPC="LHC10D1";  }
   else if (fRun>=118503&&fRun<=121040) { fLHCperiod="LHC10C"; fMCperiodTPC="LHC10D1";  }
@@ -452,7 +487,12 @@ void AliPIDResponse::SetRecoInfo()
   else if (fRun>=127719&&fRun<=130850) { fLHCperiod="LHC10E"; fMCperiodTPC="LHC10F6A"; }
   else if (fRun>=133004&&fRun<=135029) { fLHCperiod="LHC10F"; fMCperiodTPC="LHC10F6A"; }
   else if (fRun>=135654&&fRun<=136377) { fLHCperiod="LHC10G"; fMCperiodTPC="LHC10F6A"; }
-  else if (fRun>=136851&&fRun<=139517) { fLHCperiod="LHC10H"; fMCperiodTPC="LHC10H8"; fBeamType="PBPB"; }
+  else if (fRun>=136851&&fRun<=139517) {
+    fLHCperiod="LHC10H";
+    fMCperiodTPC="LHC10H8";
+    if (reg.MatchB(fCurrentFile)) fMCperiodTPC="LHC11A10";
+    fBeamType="PBPB";
+  }
   else if (fRun>=139699) { fLHCperiod="LHC11A"; fMCperiodTPC="LHC10F6A"; }
   
 }
@@ -481,11 +521,11 @@ void AliPIDResponse::SetTPCPidResponseMaster()
   
   TString fileName(Form("%s/COMMON/PID/data/TPCPIDResponse.root", fOADBPath.Data()));
   
-  TFile f(fileName.Data());
-  if (f.IsOpen() && !f.IsZombie()){
-    fArrPidResponseMaster=dynamic_cast<TObjArray*>(f.Get("TPCPIDResponse"));
-    f.Close();
+  TFile *f=TFile::Open(fileName.Data());
+  if (f && f->IsOpen() && !f->IsZombie()){
+    fArrPidResponseMaster=dynamic_cast<TObjArray*>(f->Get("TPCPIDResponse"));
   }
+  delete f;
   
   if (!fArrPidResponseMaster){
     AliFatal(Form("Could not retrieve the TPC pid response from: %s",fileName.Data()));
@@ -593,3 +633,52 @@ void AliPIDResponse::SetTPCParametrisation()
   if (fResolutionCorrection) AliInfo(Form("Setting multiplicity correction function: %s",fResolutionCorrection->GetName()));
 }
 
+//______________________________________________________________________________
+void AliPIDResponse::SetTRDPidResponseMaster()
+{
+  //
+  // Load the TRD pid params and references from the OADB
+  //
+  if(fTRDPIDParams) return;
+  AliOADBContainer contParams; 
+  contParams.InitFromFile(Form("%s/COMMON/PID/data/TRDPIDParams.root", fOADBPath.Data()), "AliTRDPIDParams");
+  fTRDPIDParams = (TObjArray *)contParams.GetObject(fRun);
+
+  AliOADBContainer contRefs;
+  contRefs.InitFromFile(Form("%s/COMMON/PID/dReferencesLQ1D.root", fOADBPath.Data()), "AliTRDPIDReference");
+  fTRDPIDReference = (AliTRDPIDReference *)contRefs.GetObject(fRun);
+}
+
+//______________________________________________________________________________
+void AliPIDResponse::InitializeTRDResponse(){
+  //
+  // Set PID Params and references to the TRD PID response
+  // 
+  fTRDResponse.SetPIDParams(fTRDPIDParams);
+  fTRDResponse.Load(fTRDPIDReference);
+  if(fLHCperiod == "LHC10b" || fLHCperiod == "LHC10c" || fLHCperiod == "LHC10d" || fLHCperiod == "LHC10e"){
+    fTRDslicesForPID[0] = 0;
+    fTRDslicesForPID[1] = 7;
+  }
+}
+
+//_________________________________________________________________________
+Bool_t AliPIDResponse::IdentifiedAsElectronTRD(const AliVTrack *vtrack, Double_t efficiencyLevel) const {
+  //
+  // Check whether track is identified as electron under a given electron efficiency hypothesis
+  //
+  Double_t probs[AliPID::kSPECIES];
+  ComputeTRDProbability(vtrack, AliPID::kSPECIES, probs);
+
+  Int_t ntracklets=0;
+  Double_t p = 0;
+  for(Int_t iPl=0;iPl<AliVTrack::kTRDnPlanes;iPl++){
+    if(vtrack->GetTRDmomentum(iPl) > 0.){
+      ntracklets++;
+      p = vtrack->GetTRDmomentum(iPl); 
+    }
+  }
+
+  return fTRDResponse.IdentifiedAsElectron(ntracklets, probs, p, efficiencyLevel);
+}
+
index c769da0..d2ba62e 100644 (file)
@@ -59,6 +59,7 @@ public:
   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 NumberOfSigmasTOF(const AliVParticle *track, AliPID::EParticleType type) const = 0;
+  virtual Bool_t IdentifiedAsElectronTRD(const AliVTrack *track, Double_t efficiencyLevel) const;
 
   EDetPidStatus ComputePIDProbability  (EDetCode detCode, const AliVTrack *track, Int_t nSpecies, Double_t p[]) const;
   
@@ -73,9 +74,11 @@ public:
 
   void SetITSPIDmethod(ITSPIDmethod pmeth) { fITSPIDmethod = pmeth; }
   virtual void SetTOFResponse(AliVEvent */*event*/,EStartTimeType_t /*option*/) {;}
-
+  void SetTRDslicesForPID(UInt_t slice1, UInt_t slice2) {fTRDslicesForPID[0]=slice1;fTRDslicesForPID[1]=slice2;}
+  
   void SetOADBPath(const char* path) {fOADBPath=path;}
   void InitialiseEvent(AliVEvent *event, Int_t pass);
+  void SetCurrentFile(const char* file) { fCurrentFile=file; }
 
   // User settings for the MC period and reco pass
   void SetMCperiod(const char *mcPeriod) {fMCperiodUser=mcPeriod;}
@@ -102,6 +105,7 @@ private:
   TString fLHCperiod;                  //! LHC period
   TString fMCperiodTPC;                //! corresponding MC period to use for the TPC splines
   TString fMCperiodUser;               //  MC prodution requested by the user
+  TString fCurrentFile;                //! name of currently processed file
   Int_t   fRecoPass;                   //! reconstruction pass
   Int_t   fRecoPassUser;               //  reconstruction pass explicitly set by the user
   Int_t   fRun;                        //! current run number
@@ -110,6 +114,10 @@ private:
   TObjArray *fArrPidResponseMaster;    //!  TPC pid splines
   TF1       *fResolutionCorrection;    //! TPC resolution correction
 
+  TObjArray *fTRDPIDParams;             //! TRD PID Params
+  AliTRDPIDReference *fTRDPIDReference; //! TRD PID References
+  UInt_t fTRDslicesForPID[2];           //! TRD PID slices
+
   Int_t   fTOFTimeZeroType;            //! default start time type for tof (ESD)
   Float_t fTOFres;                     //! TOF resolution
   
@@ -127,6 +135,10 @@ private:
   void SetTPCParametrisation();
   Double_t GetTPCMultiplicityBin(const AliVEvent * const event);
 
+  //TRD
+  void SetTRDPidResponseMaster();
+  void InitializeTRDResponse();
+
   //TOF
   
   //
index 3021cb9..c10b5cc 100644 (file)
 //  Markus Fasel <M.Fasel@gsi.de>
 //  Anton Andronic <A.Andronic@gsi.de>
 //
-#include <TClass.h>
 #include <TAxis.h>
+#include <TClass.h>
+#include <TDirectory.h>
 #include <TFile.h>
 #include <TH1.h>
 #include <TKey.h>
+#include <TMath.h>
 #include <TObjArray.h>
-#include <TString.h>
 #include <TROOT.h> 
+#include <TString.h>
 #include <TSystem.h>
-#include <TDirectory.h>
+#include <TVectorT.h>
 
 #include "AliLog.h"
 
@@ -44,6 +46,7 @@ ClassImp(AliTRDPIDResponse)
 AliTRDPIDResponse::AliTRDPIDResponse():
   TObject()
   ,fkPIDReference(NULL)
+  ,fkPIDParams(NULL)
   ,fGainNormalisationFactor(1.)
   ,fPIDmethod(kLQ1D)
 {
@@ -56,6 +59,7 @@ AliTRDPIDResponse::AliTRDPIDResponse():
 AliTRDPIDResponse::AliTRDPIDResponse(const AliTRDPIDResponse &ref):
   TObject(ref)
   ,fkPIDReference(ref.fkPIDReference)
+  ,fkPIDParams(ref.fkPIDParams)
   ,fGainNormalisationFactor(ref.fGainNormalisationFactor)
   ,fPIDmethod(ref.fPIDmethod)
 {
@@ -75,6 +79,7 @@ AliTRDPIDResponse &AliTRDPIDResponse::operator=(const AliTRDPIDResponse &ref){
   TObject::operator=(ref);
   fGainNormalisationFactor = ref.fGainNormalisationFactor;
   fkPIDReference = ref.fkPIDReference;
+  fkPIDParams = ref.fkPIDParams;
   fPIDmethod = ref.fPIDmethod;
   
   return *this;
@@ -237,3 +242,56 @@ Bool_t AliTRDPIDResponse::CookdEdx(Int_t nSlice, const Double_t * const in, Doub
   return kTRUE;
 }
 
+//____________________________________________________________
+Bool_t AliTRDPIDResponse::IdentifiedAsElectron(Int_t nTracklets, const Double_t *like, Double_t p, Double_t level) const {
+  //
+  // Check whether particle is identified as electron assuming a certain electron efficiency level
+  // Only electron and pion hypothesis is taken into account
+  //
+  // Inputs:
+  //         Number of tracklets
+  //         Likelihood values
+  //         Momentum
+  //         Electron efficiency level
+  //
+  // If the function fails when the params are not accessible, the function returns true
+  //
+  if(!fkPIDParams){
+    AliError("No PID Param object available");
+    return kTRUE;
+  } 
+  Double_t probEle = like[AliPID::kElectron]/(like[AliPID::kElectron] + like[AliPID::kPion]);
+  const TVectorD *params = GetParams(nTracklets, level);
+  if(!params){
+    AliError("No Params found for the given configuration");
+    return kTRUE;
+  }
+  Double_t threshold = 1. - (*params)[0] - (*params)[1] * p - (*params)[2] * TMath::Exp(-(*params)[3] * p);
+  if(probEle > TMath::Max(TMath::Min(threshold, 0.99), 0.2)) return kTRUE; // truncate the threshold upperwards to 0.999 and lowerwards to 0.2 and exclude unphysical values
+  return kFALSE;
+}
+
+//____________________________________________________________
+const TVectorD* AliTRDPIDResponse::GetParams(Int_t ntracklets, Double_t level) const {
+  //
+  // returns the threshold for a given number of tracklets and a given efficiency level
+  //tby definition the lower of step is given.
+  //
+  if(ntracklets > 6 || ntracklets) return NULL;
+  TObjArray * entry = dynamic_cast<TObjArray *>(fkPIDParams->At(ntracklets - 1));
+  if(!entry) return NULL;
+  
+  TObjArray*cut = NULL;
+  TVectorF *effLevel = NULL; const TVectorD *parameters = NULL;
+  Float_t currentLower = 0.;
+  TIter cutIter(entry);
+  while((cut = dynamic_cast<TObjArray *>(cutIter()))){
+    effLevel = static_cast<TVectorF *>(cut->At(0));
+    if(effLevel[0] > currentLower && effLevel[0] <= level){
+      // New Lower entry found
+      parameters = static_cast<const TVectorD *>(cut->At(1));
+    }
+  }  
+
+  return parameters;
+}
index 31b0ff1..cb462a4 100644 (file)
@@ -28,6 +28,8 @@
 #include "AliPID.h"
 #endif
 
+#include <TVectorDfwd.h>
+
 class TObjArray;
 class AliVTrack;
 class AliTRDPIDReference;
@@ -64,17 +66,23 @@ class AliTRDPIDResponse : public TObject {
     void      SetOwner();
     void      SetPIDmethod(ETRDPIDMethod m) {fPIDmethod=m;}
     void      SetGainNormalisationFactor(Double_t gainFactor) { fGainNormalisationFactor = gainFactor; }
+    void      SetPIDParams(const TObjArray * params) { fkPIDParams = params; }
 
     Bool_t    Load(const Char_t *filename = NULL, const Char_t *refName = "RefTRDLQ1D");
     Bool_t    Load(const AliTRDPIDReference *ref) { fkPIDReference = ref; return kTRUE; }
   
+    Bool_t    IdentifiedAsElectron(Int_t nTracklets, const Double_t *like, Double_t p, Double_t level) const;
+  
   private:
     Bool_t    CookdEdx(Int_t nSlice, const Double_t * const in, Double_t *out) const;
     Double_t  GetProbabilitySingleLayer(Int_t species, Double_t plocal, Double_t dEdx) const;
+    const TVectorD * GetParams(Int_t ntracklets, Double_t level) const;
     
     const AliTRDPIDReference *fkPIDReference;   // PID References
-    Double_t  fGainNormalisationFactor;  // Gain normalisation factor
-    ETRDPIDMethod   fPIDmethod;   // PID method selector  
+    const TObjArray *fkPIDParams;               // PID Params
+    Double_t  fGainNormalisationFactor;         // Gain normalisation factor
+    ETRDPIDMethod   fPIDmethod;                 // PID method selector
+      
   
   ClassDef(AliTRDPIDResponse, 3)    // Tool for TRD PID
 };
index 80636f6..87735d6 100644 (file)
@@ -61,9 +61,10 @@ public:
   virtual UShort_t  GetTPCsignalN()      const {return 0 ;}
   virtual Double_t  GetTPCmomentum()     const {return 0.;}
   virtual Double_t  GetTOFsignal()       const {return 0.;}
-  virtual void GetIntegratedTimes(Double_t */*times*/) const { return; }
-  virtual Double_t GetTRDmomentum(Int_t /*plane*/, Double_t */*sp*/=0x0) const {return 0.;}
-  
+  virtual void      GetIntegratedTimes(Double_t */*times*/) const { return; }
+  virtual Double_t  GetTRDmomentum(Int_t /*plane*/, Double_t */*sp*/=0x0) const {return 0.;}
+  virtual void      GetHMPIDpid(Double_t */*p*/) const {;}
+    
   virtual ULong_t  GetStatus() const = 0;
   virtual Bool_t   GetXYZ(Double_t *p) const = 0;
   virtual Double_t GetBz() const;