PID updates for TRD
authormorsch <morsch@f7af4fe6-9843-0410-8265-dc069ae4e863>
Mon, 19 Sep 2011 09:48:40 +0000 (09:48 +0000)
committermorsch <morsch@f7af4fe6-9843-0410-8265-dc069ae4e863>
Mon, 19 Sep 2011 09:48:40 +0000 (09:48 +0000)
STEER/AOD/AliAODTrack.cxx
STEER/AOD/AliAODTrack.h
STEER/STEERBase/AliPIDResponse.cxx
STEER/STEERBase/AliTRDPIDResponse.cxx
STEER/STEERBase/AliVTrack.h

index a8e717c..4462c58 100644 (file)
@@ -639,7 +639,7 @@ Double_t  AliAODTrack::GetTRDslice(Int_t plane, Int_t slice) const {
     return -1.;
   }
 
-  Int_t ns=fDetPid->GetTRDnSlices();
+  Int_t ns=fDetPid->GetTRDnSlices()/kTRDnPlanes;
   if ((slice<-1) || (slice>=ns)) {
     return -1.;
   }
@@ -652,6 +652,27 @@ Double_t  AliAODTrack::GetTRDslice(Int_t plane, Int_t slice) const {
   return q/ns;
 }
 
+//______________________________________________________________________________
+UChar_t AliAODTrack::GetTRDntrackletsPID() const{
+  //
+  // return number of tracklets calculated from the slices
+  //
+  if(!fDetPid) return -1;
+
+  Int_t ntracklets = 0,                                           // Number of tracklets / track
+        nSlicesTracklet = fDetPid->GetTRDnSlices()/kTRDnPlanes,   // Number of slices per tracklet
+        nSlicesNonZero = 0;                                       // Number of slices containing a dE/dx measurement
+  for(Int_t ily = 0; ily < kTRDnPlanes; ily++){
+    // a tracklet is found if it has at least one slice containing a dE/dx measurement
+    nSlicesNonZero = 0;
+    for(Int_t islice = 0; islice < nSlicesTracklet; islice++){
+      if(fDetPid->GetTRDsignal()[nSlicesTracklet * ily + islice] > 0.01) nSlicesNonZero++;
+    }
+    if(nSlicesNonZero) ntracklets++;
+  }
+  return ntracklets;
+}
+
 //______________________________________________________________________________
 Double_t AliAODTrack::GetTRDmomentum(Int_t plane, Double_t */*sp*/) const
 {
index 6448ac9..b3c84cf 100644 (file)
@@ -262,6 +262,7 @@ class AliAODTrack : public AliVTrack {
   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;
+  UChar_t   GetTRDntrackletsPID() const;
   void      GetHMPIDpid(Double_t *p) const { if (fDetPid) fDetPid->GetHMPIDprobs(p); }
 
   
index 27da122..2274207 100644 (file)
@@ -579,7 +579,7 @@ void AliPIDResponse::SetRecoInfo()
   if (fRun>=114737&&fRun<=117223)      { fLHCperiod="LHC10B"; fMCperiodTPC="LHC10D1";  }
   else if (fRun>=118503&&fRun<=121040) { fLHCperiod="LHC10C"; fMCperiodTPC="LHC10D1";  }
   else if (fRun>=122195&&fRun<=126437) { fLHCperiod="LHC10D"; fMCperiodTPC="LHC10F6A"; }
-  else if (fRun>=127719&&fRun<=130850) { fLHCperiod="LHC10E"; fMCperiodTPC="LHC10F6A"; }
+  else if (fRun>=127710&&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) {
@@ -768,14 +768,16 @@ Bool_t AliPIDResponse::IdentifiedAsElectronTRD(const AliVTrack *vtrack, Double_t
   Double_t probs[AliPID::kSPECIES];
   ComputeTRDProbability(vtrack, AliPID::kSPECIES, probs);
 
-  Int_t ntracklets=0;
-  Double_t p = 0;
+  Int_t ntracklets = vtrack->GetTRDntrackletsPID();
+  // Take mean of the TRD momenta in the given tracklets
+  Float_t p = 0, trdmomenta[AliVTrack::kTRDnPlanes];
+  Int_t nmomenta = 0;
   for(Int_t iPl=0;iPl<AliVTrack::kTRDnPlanes;iPl++){
     if(vtrack->GetTRDmomentum(iPl) > 0.){
-      ntracklets++;
-      p = vtrack->GetTRDmomentum(iPl); 
+      trdmomenta[nmomenta++] = vtrack->GetTRDmomentum(iPl); 
     }
   }
+  p = TMath::Mean(nmomenta, trdmomenta);
 
   return fTRDResponse.IdentifiedAsElectron(ntracklets, probs, p, efficiencyLevel);
 }
index 2f40336..e9256a1 100644 (file)
@@ -233,7 +233,8 @@ Bool_t AliTRDPIDResponse::CookdEdx(Int_t nSlice, const Double_t * const in, Doub
     break;
   case kLQ1D: // 1D LQ 
     out[0]= 0.;
-    for(Int_t islice = 0; islice < nSlice; islice++) out[0] += in[islice] * fGainNormalisationFactor;
+    for(Int_t islice = 0; islice < nSlice; islice++) 
+      if(in[islice] > 0) out[0] += in[islice] * fGainNormalisationFactor;   // Protect against negative values for slices having no dE/dx information
     if(out[0] < 1e-6) return kFALSE;
     break;
   default:
@@ -293,6 +294,7 @@ const TVectorD* AliTRDPIDResponse::GetParams(Int_t ntracklets, Double_t level) c
       currentLower = (*effLevel)[0];
     }
   }  
+  AliDebug(2, Form("Taking params for %d tracklets and %f electron Efficiency\n", ntracklets, currentLower));
 
   return parameters;
 }
index 5155277..9a98123 100644 (file)
@@ -55,6 +55,7 @@ public:
   virtual UShort_t GetTPCNcls() const { return 0;}\r
   virtual UShort_t GetTPCNclsF() const { return 0;}\r
   virtual Double_t GetTRDslice(Int_t /*plane*/, Int_t /*slice*/) const { return -1.; }\r
+  virtual UChar_t  GetTRDntrackletsPID() const { return 0;}\r
   \r
   virtual Int_t GetEMCALcluster()     const {return -1;}\r
   virtual void SetEMCALcluster(Int_t)       {;}\r