-add nsigma electron PID corrections
authorjbook <jbook@f7af4fe6-9843-0410-8265-dc069ae4e863>
Thu, 4 Apr 2013 13:29:29 +0000 (13:29 +0000)
committerjbook <jbook@f7af4fe6-9843-0410-8265-dc069ae4e863>
Thu, 4 Apr 2013 13:29:29 +0000 (13:29 +0000)
-add Nacc and NTrk for nanoAODs

PWGDQ/dielectron/AliDielectronPID.cxx
PWGDQ/dielectron/AliDielectronPID.h
PWGDQ/dielectron/AliDielectronVarManager.h

index 4732c9b..e3fac40 100644 (file)
@@ -50,6 +50,8 @@ TGraph *AliDielectronPID::fgFitCorr=0x0;
 Double_t AliDielectronPID::fgCorr=0.0;
 Double_t AliDielectronPID::fgCorrdEdx=1.0;
 TF1 *AliDielectronPID::fgFunEtaCorr=0x0;
+TF1 *AliDielectronPID::fgFunCntrdCorr=0x0;
+TF1 *AliDielectronPID::fgFunWdthCorr=0x0;
 TGraph *AliDielectronPID::fgdEdxRunCorr=0x0;
 
 AliDielectronPID::AliDielectronPID() :
@@ -396,6 +398,8 @@ Bool_t AliDielectronPID::IsSelectedTPC(AliVTrack * const part, Int_t icut)
 
   if (fPartType[icut]==AliPID::kElectron){
     numberOfSigmas-=fgCorr;
+    numberOfSigmas-=GetCntrdCorr(part);
+    numberOfSigmas/=GetWdthCorr(part);
   }
   
   // test if we are supposed to use a function for the cut
@@ -648,3 +652,41 @@ Double_t AliDielectronPID::GetEtaCorr(const AliVTrack *track)
   if (!fgFunEtaCorr) return 1;
   return fgFunEtaCorr->Eval(track->Eta());
 }
+
+//______________________________________________
+void AliDielectronPID::SetCentroidCorrFunction(TF1 *fun, UInt_t varx, UInt_t vary, UInt_t varz)
+{
+  fun->GetHistogram()->GetXaxis()->SetUniqueID(varx);
+  fun->GetHistogram()->GetYaxis()->SetUniqueID(vary);
+  fun->GetHistogram()->GetZaxis()->SetUniqueID(varz);
+  fgFunCntrdCorr=fun;
+}
+//______________________________________________
+void AliDielectronPID::SetWidthCorrFunction(TF1 *fun, UInt_t varx, UInt_t vary, UInt_t varz)
+{
+  fun->GetHistogram()->GetXaxis()->SetUniqueID(varx);
+  fun->GetHistogram()->GetYaxis()->SetUniqueID(vary);
+  fun->GetHistogram()->GetZaxis()->SetUniqueID(varz);
+  fgFunWdthCorr=fun;
+}
+
+//______________________________________________
+Double_t AliDielectronPID::GetPIDCorr(const AliVTrack *track, TF1 *fun)
+{
+  //
+  // return correction value
+  //
+
+  //Fill only event and vparticle values (otherwise we end up in a circle)
+  Double_t values[AliDielectronVarManager::kNMaxValues];
+  AliDielectronVarManager::FillVarVParticle(track,values);
+
+  Int_t dim=fun->GetNdim();
+  Double_t var[3] = {0.,0.,0.};
+  if(dim>0) var[0] = values[fun->GetHistogram()->GetXaxis()->GetUniqueID()];
+  if(dim>1) var[1] = values[fun->GetHistogram()->GetYaxis()->GetUniqueID()];
+  if(dim>2) var[2] = values[fun->GetHistogram()->GetZaxis()->GetUniqueID()];
+  Double_t corr = fun->Eval(var[0],var[1],var[2]);
+  // printf(" %d-dim CORR value: %f (track %p) \n",dim,corr,track);
+  return corr;
+}
index 911664a..8a6bb0f 100644 (file)
@@ -79,9 +79,15 @@ public:
   static TGraph *GetCorrGraphdEdx()  { return fgdEdxRunCorr; }
 
   static void SetEtaCorrFunction(TF1 *fun) {fgFunEtaCorr=fun;}
+  static void SetCentroidCorrFunction(TF1 *fun, UInt_t varx, UInt_t vary=0, UInt_t varz=0);
+  static void SetWidthCorrFunction(TF1 *fun, UInt_t varx, UInt_t vary=0, UInt_t varz=0);
   static TF1* GetEtaCorrFunction() { return fgFunEtaCorr; }
+  static TF1* GetCentroidCorrFunction() { return fgFunCntrdCorr; }
+  static TF1* GetWidthCorrFunction() { return fgFunWdthCorr; }
 
   static Double_t GetEtaCorr(const AliVTrack *track);
+  static Double_t GetCntrdCorr(const AliVTrack *track) { return (fgFunCntrdCorr ? GetPIDCorr(track,fgFunCntrdCorr) : 0.0); }
+  static Double_t GetWdthCorr(const AliVTrack *track)  { return (fgFunWdthCorr  ? GetPIDCorr(track,fgFunWdthCorr)  : 1.0); }
 
 private:
   enum {kNmaxPID=30};
@@ -111,9 +117,11 @@ private:
   static Double_t fgCorrdEdx;     //!dEdx correction value for current run. Set if fgFitCorr is set and SetCorrVal(run)
                                   // was called
   static TF1    *fgFunEtaCorr;    //function for eta correction of electron sigma
+  static TF1    *fgFunCntrdCorr;  //function for correction of electron sigma (centroid)
+  static TF1    *fgFunWdthCorr;   //function for correction of electron sigma (width)
   static TGraph *fgdEdxRunCorr;   //run by run correction for dEdx
 
-
+  static Double_t GetPIDCorr(const AliVTrack *track, TF1 *fun);
   
   Bool_t IsSelectedITS(AliVTrack * const part, Int_t icut);
   Bool_t IsSelectedTPC(AliVTrack * const part, Int_t icut);
@@ -125,7 +133,7 @@ private:
   AliDielectronPID(const AliDielectronPID &c);
   AliDielectronPID &operator=(const AliDielectronPID &c);
 
-  ClassDef(AliDielectronPID,4)         // Dielectron PID
+  ClassDef(AliDielectronPID,5)         // Dielectron PID
 };
 
 #endif
index fcfcf86..e449566 100644 (file)
@@ -381,6 +381,7 @@ public:
   virtual ~AliDielectronVarManager();
   static void Fill(const TObject* particle, Double_t * const values);
   static void FillVarMCParticle2(const AliVParticle *p1, const AliVParticle *p2, Double_t * const values);
+  static void FillVarVParticle(const AliVParticle *particle,         Double_t * const values);
 
   static void InitESDpid(Int_t type=0);
   static void InitAODpidUtil(Int_t type=0);
@@ -417,7 +418,7 @@ private:
 
   static const char* fgkParticleNames[kNMaxValues][3];  //variable names
 
-  static void FillVarVParticle(const AliVParticle *particle,         Double_t * const values);
+
   static void FillVarESDtrack(const AliESDtrack *particle,           Double_t * const values);
   static void FillVarAODTrack(const AliAODTrack *particle,           Double_t * const values);
   static void FillVarMCParticle(const AliMCParticle *particle,       Double_t * const values);
@@ -689,8 +690,7 @@ inline void AliDielectronVarManager::FillVarESDtrack(const AliESDtrack *particle
   // nsigma to Electron band
   // TODO: for the moment we set the bethe bloch parameters manually
   //       this should be changed in future!
-  
-  values[AliDielectronVarManager::kTPCnSigmaEle]=fgPIDResponse->NumberOfSigmasTPC(particle,AliPID::kElectron)-AliDielectronPID::GetCorrVal();
+  values[AliDielectronVarManager::kTPCnSigmaEle]=(fgPIDResponse->NumberOfSigmasTPC(particle,AliPID::kElectron)-AliDielectronPID::GetCorrVal()-AliDielectronPID::GetCntrdCorr(particle)) / AliDielectronPID::GetWdthCorr(particle);
   values[AliDielectronVarManager::kTPCnSigmaPio]=fgPIDResponse->NumberOfSigmasTPC(particle,AliPID::kPion);
   values[AliDielectronVarManager::kTPCnSigmaMuo]=fgPIDResponse->NumberOfSigmasTPC(particle,AliPID::kMuon);
   values[AliDielectronVarManager::kTPCnSigmaKao]=fgPIDResponse->NumberOfSigmasTPC(particle,AliPID::kKaon);
@@ -836,7 +836,7 @@ inline void AliDielectronVarManager::FillVarAODTrack(const AliAODTrack *particle
     values[AliDielectronVarManager::kTPCsignalN]     = tpcSignalN;
     values[AliDielectronVarManager::kTPCsignalNfrac] = tpcNcls>0?tpcSignalN/tpcNcls:0;
     values[AliDielectronVarManager::kTPCclsDiff]     = tpcSignalN-tpcNcls;
-    Double_t tpcNsigmaEle=fgPIDResponse->NumberOfSigmasTPC(particle,AliPID::kElectron);
+    Double_t tpcNsigmaEle=(fgPIDResponse->NumberOfSigmasTPC(particle,AliPID::kElectron)-AliDielectronPID::GetCntrdCorr(particle))/AliDielectronPID::GetWdthCorr(particle);
     Double_t tpcNsigmaPio=fgPIDResponse->NumberOfSigmasTPC(particle,AliPID::kPion);
     Double_t tpcNsigmaMuo=fgPIDResponse->NumberOfSigmasTPC(particle,AliPID::kMuon);
     Double_t tpcNsigmaKao=fgPIDResponse->NumberOfSigmasTPC(particle,AliPID::kKaon);
@@ -891,9 +891,9 @@ inline void AliDielectronVarManager::FillVarAODTrack(const AliAODTrack *particle
     values[AliDielectronVarManager::kTOFnSigmaKao]=tofNsigmaKao;
     values[AliDielectronVarManager::kTOFnSigmaPro]=tofNsigmaPro;
 
-       values[AliDielectronVarManager::kTOFmismProb] = fgPIDResponse->GetTOFMismatchProbability(particle);
+    values[AliDielectronVarManager::kTOFmismProb] = fgPIDResponse->GetTOFMismatchProbability(particle);
   
-       pid->SetTPCsignal(origdEdx);
+    pid->SetTPCsignal(origdEdx);
   }
 
   //EMCAL PID information
@@ -1803,6 +1803,10 @@ inline void AliDielectronVarManager::FillVarAODEvent(const AliAODEvent *event, D
     values[AliDielectronVarManager::kCentrality] = header->GetCentrality();
   // nanoAODs (w/o AliEventPlane branch) should have the tpc event plane angle stored in the header
   if(!header->GetEventplaneP()) {
+
+    values[AliDielectronVarManager::kNTrk] = header->GetRefMultiplicity();    // overwritten datamembers in "our" nanoAODs
+    values[AliDielectronVarManager::kNacc] = header->GetRefMultiplicityPos(); // overwritten datamembers in "our" nanoAODs
+
     // TPC
     values[AliDielectronVarManager::kTPCrpH2uc]  = header->GetEventplane();
     values[AliDielectronVarManager::kTPCmagH2uc] = header->GetEventplaneMag();