Merge branch 'master' of https://git.cern.ch/reps/AliRoot

[u/mrichter/AliRoot.git] / STEER / AOD / AliAODpidUtil.cxx

diff --git a/STEER/AOD/AliAODpidUtil.cxx b/STEER/AOD/AliAODpidUtil.cxx
index 03885f2cf54c5e81f256f3c542cdbef7524db97a..f7c41b06cceef8463d5484a08140972f5f10543f 100644 (file)
--- a/STEER/AOD/AliAODpidUtil.cxx
+++ b/STEER/AOD/AliAODpidUtil.cxx
@@ -33,6 +33,7 @@
 #include "AliESDtrack.h"
 #include "AliAODMCHeader.h"
 #include "AliAODMCParticle.h"
+#include "AliTOFPIDParams.h"
 
 #include <AliDetectorPID.h>
 
@@ -44,8 +45,6 @@ Float_t AliAODpidUtil::GetTPCsignalTunedOnData(const AliVTrack *t) const {
     Float_t dedx = track->GetTPCsignalTunedOnData();
     if(dedx > 0) return dedx;
 
-    Double_t mom = t->GetTPCmomentum();
-
     dedx = t->GetTPCsignal();
     track->SetTPCsignalTunedOnData(dedx);
 
@@ -56,47 +55,53 @@ Float_t AliAODpidUtil::GetTPCsignalTunedOnData(const AliVTrack *t) const {
     
     AliAODMCHeader *mcHeader = dynamic_cast<AliAODMCHeader*>(track->GetAODEvent()->GetList()->FindObject(AliAODMCHeader::StdBranchName()));
     if (mcHeader) {
+        
        TClonesArray *mcArray = (TClonesArray*)track->GetAODEvent()->GetList()->FindObject(AliAODMCParticle::StdBranchName());
        
        Bool_t kGood = kTRUE;
-       
-       Int_t iS = TMath::Abs(((AliAODMCParticle*)mcArray->At(TMath::Abs(t->GetLabel())))->GetPdgCode());
-       if(iS==AliPID::ParticleCode(AliPID::kElectron)){
+
+       if ( mcArray->At(TMath::Abs(t->GetLabel())) != NULL ) {  // protects against label-0 tracks e.g. the initial proton for Phythia events
+         Int_t iS = TMath::Abs(((AliAODMCParticle*)mcArray->At(TMath::Abs(t->GetLabel())))->GetPdgCode());
+         if(iS==AliPID::ParticleCode(AliPID::kElectron)){
            type = AliPID::kElectron;
-       }
-       else if(iS==AliPID::ParticleCode(AliPID::kMuon)){
+         }
+         else if(iS==AliPID::ParticleCode(AliPID::kMuon)){
            type = AliPID::kMuon;
-       }
-       else if(iS==AliPID::ParticleCode(AliPID::kPion)){
+         }
+         else if(iS==AliPID::ParticleCode(AliPID::kPion)){
            type = AliPID::kPion;
-       }
-       else if(iS==AliPID::ParticleCode(AliPID::kKaon)){
+         }
+         else if(iS==AliPID::ParticleCode(AliPID::kKaon)){
            type = AliPID::kKaon;
-       }
-       else if(iS==AliPID::ParticleCode(AliPID::kProton)){
+         }
+         else if(iS==AliPID::ParticleCode(AliPID::kProton)){
            type = AliPID::kProton;
-       }
-       else if(iS==AliPID::ParticleCode(AliPID::kDeuteron)){ // d
+         }
+         else if(iS==AliPID::ParticleCode(AliPID::kDeuteron)){ // d
            type = AliPID::kDeuteron;
-       }
-       else if(iS==AliPID::ParticleCode(AliPID::kTriton)){ // t
+         }
+         else if(iS==AliPID::ParticleCode(AliPID::kTriton)){ // t
            type = AliPID::kTriton;
-       }
-       else if(iS==AliPID::ParticleCode(AliPID::kHe3)){ // 3He
+         }
+         else if(iS==AliPID::ParticleCode(AliPID::kHe3)){ // 3He
            type = AliPID::kHe3;
-       }
-       else if(iS==AliPID::ParticleCode(AliPID::kAlpha)){ // 4He
+         }
+         else if(iS==AliPID::ParticleCode(AliPID::kAlpha)){ // 4He
            type = AliPID::kAlpha;
-       }
+         }
        else
            kGood = kFALSE;
+       } else kGood = kFALSE;
 
        if(kGood){
-           Double_t bethe=fTPCResponse.GetExpectedSignal(mom,type);
-           Double_t sigma=fTPCResponse.GetExpectedSigma(mom,t->GetTPCsignalN(),type);
-           dedx = gRandom->Gaus(bethe,sigma);
-
-           if(iS == AliPID::ParticleCode(AliPID::kHe3) || iS == AliPID::ParticleCode(AliPID::kAlpha)) dedx *= 5;
+           //TODO maybe introduce different dEdxSources?
+        Double_t bethe = fTPCResponse.GetExpectedSignal(track, type, AliTPCPIDResponse::kdEdxDefault, this->UseTPCEtaCorrection(),
+                                                        this->UseTPCMultiplicityCorrection());
+        Double_t sigma = fTPCResponse.GetExpectedSigma(track, type, AliTPCPIDResponse::kdEdxDefault, this->UseTPCEtaCorrection(),
+                                                       this->UseTPCMultiplicityCorrection());
+        dedx = gRandom->Gaus(bethe,sigma);
+        
+//         if(iS == AliPID::ParticleCode(AliPID::kHe3) || iS == AliPID::ParticleCode(AliPID::kAlpha)) dedx *= 5;
        }
 
     }
@@ -105,24 +110,77 @@ Float_t AliAODpidUtil::GetTPCsignalTunedOnData(const AliVTrack *t) const {
     return dedx;
 }
 //_________________________________________________________________________
-Float_t AliAODpidUtil::NumberOfSigmasTOF(const AliVParticle *vtrack, AliPID::EParticleType type) const
+Float_t AliAODpidUtil::GetTOFsignalTunedOnData(const AliVTrack *t) const {
+    AliAODTrack *track = (AliAODTrack *) t;
+    Double_t tofSignal = track->GetTOFsignalTunedOnData();
+
+    if(tofSignal <  99999) return (Float_t)tofSignal; // it has been already set
+
+    // read additional mismatch fraction
+    Float_t addmism = GetTOFPIDParams()->GetTOFadditionalMismForMC();
+    if(addmism > 1.){
+      Float_t centr = GetCurrentCentrality();
+      if(centr > 50) addmism *= 0.1667;
+      else if(centr > 20) addmism *= 0.33;
+    }
+
+    AliAODPid *pidObj = track->GetDetPid();
+    tofSignal = pidObj->GetTOFsignal() + fTOFResponse.GetTailRandomValue(t->Pt(),t->Eta(),pidObj->GetTOFsignal(),addmism);
+    track->SetTOFsignalTunedOnData(tofSignal);
+    return (Float_t)tofSignal;
+}
+
+//_________________________________________________________________________
+Float_t AliAODpidUtil::GetSignalDeltaTOFold(const AliVParticle *vtrack, AliPID::EParticleType type, Bool_t ratio/*=kFALSE*/) const
 {
   //
   // Number of sigma implementation for the TOF
   //
   
   AliAODTrack *track=(AliAODTrack*)vtrack;
+  AliAODPid *pidObj = track->GetDetPid();
+  if (!pidObj) return -9999.;
+  Double_t tofTime = 99999;
+  if (fTuneMConData && ((fTuneMConDataMask & kDetTOF) == kDetTOF) ) tofTime = (Double_t)this->GetTOFsignalTunedOnData((AliVTrack*)vtrack);
+  else tofTime=pidObj->GetTOFsignal();
+  const Double_t expTime=fTOFResponse.GetExpectedSignal((AliVTrack*)vtrack,type);
+  Double_t sigmaTOFPid[AliPID::kSPECIES];
+  pidObj->GetTOFpidResolution(sigmaTOFPid);
+  AliAODEvent *event=(AliAODEvent*)track->GetAODEvent();
+  if (event) {  // protection if the user didn't call GetTrack, which sets the internal pointer
+    AliTOFHeader* tofH=(AliTOFHeader*)event->GetTOFHeader();
+    if (tofH && (TMath::Abs(sigmaTOFPid[0]) <= 1.E-16) ) { // new AOD
+        tofTime -= fTOFResponse.GetStartTime(vtrack->P());
+    }
+  } else {
+    AliError("pointer to AliAODEvent not found, please call GetTrack to set it");
+    return -9999.;
+  }
 
-  // look for cached value first
-  if (track->GetDetectorPID()){
-    return track->GetDetectorPID()->GetNumberOfSigmas(kTOF, type);
-  }  
+  Double_t delta=-9999.;
+
+  if (!ratio) delta=tofTime-expTime;
+  else if (expTime>1.e-20) delta=tofTime/expTime;
   
+  return delta;
+}
+
+//_________________________________________________________________________
+Float_t AliAODpidUtil::GetNumberOfSigmasTOFold(const AliVParticle *vtrack, AliPID::EParticleType type) const
+{
+  //
+  // Number of sigma implementation for the TOF
+  //
+  
+  AliAODTrack *track=(AliAODTrack*)vtrack;
+
   Bool_t oldAod=kTRUE;
-  Double_t sigTOF;
+  Double_t sigTOF=0.;
   AliAODPid *pidObj = track->GetDetPid();
   if (!pidObj) return -999.;
-  Double_t tofTime=pidObj->GetTOFsignal();
+  Double_t tofTime = 99999;
+  if (fTuneMConData && ((fTuneMConDataMask & kDetTOF) == kDetTOF) ) tofTime = (Double_t)this->GetTOFsignalTunedOnData((AliVTrack*)vtrack);
+  else tofTime=pidObj->GetTOFsignal();
   Double_t expTime=fTOFResponse.GetExpectedSignal((AliVTrack*)vtrack,type);
   Double_t sigmaTOFPid[AliPID::kSPECIES];
   pidObj->GetTOFpidResolution(sigmaTOFPid);