Update in application of weights for Bayesian PID (Jeremy)

[u/mrichter/AliRoot.git] / PWGHF / vertexingHF / AliVertexingHFUtils.cxx

diff --git a/PWGHF/vertexingHF/AliVertexingHFUtils.cxx b/PWGHF/vertexingHF/AliVertexingHFUtils.cxx
index 5441ccb331ee5caea74d0f15cfed579d6755fe70..94b48d926d9d6f73d990900b148408d97db33f07 100644 (file)
--- a/PWGHF/vertexingHF/AliVertexingHFUtils.cxx
+++ b/PWGHF/vertexingHF/AliVertexingHFUtils.cxx
@@ -22,6 +22,7 @@
 #include <TF1.h>
 #include "AliAODEvent.h"
 #include "AliAODMCHeader.h"
+#include "AliGenEventHeader.h"
 #include "AliAODMCParticle.h"
 #include "AliAODRecoDecayHF.h"
 #include "AliVertexingHFUtils.h"
@@ -63,6 +64,24 @@ AliVertexingHFUtils::AliVertexingHFUtils(Int_t k):
 }
 
 
+//______________________________________________________________________
+void AliVertexingHFUtils::ComputeSignificance(Double_t signal, Double_t  errsignal, Double_t  background, Double_t  errbackground, Double_t &significance,Double_t &errsignificance){
+  // calculate significance from S, B and errors
+
+
+  Double_t errSigSq=errsignal*errsignal;
+  Double_t errBkgSq=errbackground*errbackground;
+  Double_t sigPlusBkg=signal+background;
+  if(sigPlusBkg>0. && signal>0.){
+    significance =  signal/TMath::Sqrt(signal+background);
+    errsignificance = significance*TMath::Sqrt((errSigSq+errBkgSq)/(4.*sigPlusBkg*sigPlusBkg)+(background/sigPlusBkg)*errSigSq/signal/signal);
+  }else{
+    significance=0.;
+    errsignificance=0.;
+  }
+  return;
+
+}
 //______________________________________________________________________
 Double_t AliVertexingHFUtils::Pol(Double_t x, Int_t k){
   // compute chi from polynomial approximation
@@ -294,6 +313,49 @@ void AliVertexingHFUtils::AveragePt(Float_t& averagePt, Float_t& errorPt,Float_t
   delete hMassDpt;
   delete hMptSig;
 
+}
+//____________________________________________________________________________
+Double_t AliVertexingHFUtils::GetTrueImpactParameterDzero(AliAODMCHeader *mcHeader, TClonesArray* arrayMC, AliAODMCParticle *partD) {
+  // true impact parameter calculation for Dzero
+
+  if(!partD || !arrayMC || !mcHeader) return 99999.;
+  Int_t code=partD->GetPdgCode();
+  if(TMath::Abs(code)!=421) return 99999.;
+
+  Double_t vtxTrue[3];
+  mcHeader->GetVertex(vtxTrue);
+  Double_t origD[3];
+  partD->XvYvZv(origD);
+  Short_t charge=partD->Charge();
+  Double_t pXdauTrue[2],pYdauTrue[2],pZdauTrue[2];
+  for(Int_t iDau=0; iDau<2; iDau++){
+    pXdauTrue[iDau]=0.;
+    pYdauTrue[iDau]=0.;
+    pZdauTrue[iDau]=0.;
+  }
+
+  Int_t nDau=partD->GetNDaughters();
+  Int_t labelFirstDau = partD->GetDaughter(0); 
+  if(nDau==2){
+    for(Int_t iDau=0; iDau<2; iDau++){
+      Int_t ind = labelFirstDau+iDau;
+      AliAODMCParticle* part = dynamic_cast<AliAODMCParticle*>(arrayMC->At(ind));
+      if(!part){
+       printf("Daughter particle not found in MC array");
+       return 99999.;
+      } 
+      pXdauTrue[iDau]=part->Px();
+      pYdauTrue[iDau]=part->Py();
+      pZdauTrue[iDau]=part->Pz();
+    }
+  }else{
+    return 99999.;
+  }
+
+  Double_t d0dummy[2]={0.,0.};
+  AliAODRecoDecayHF aodDvsMC(vtxTrue,origD,2,charge,pXdauTrue,pYdauTrue,pZdauTrue,d0dummy);
+  return aodDvsMC.ImpParXY();
+
 }
 //____________________________________________________________________________
 Double_t AliVertexingHFUtils::GetTrueImpactParameterDplus(AliAODMCHeader *mcHeader, TClonesArray* arrayMC, AliAODMCParticle *partD) {
@@ -406,4 +468,18 @@ Double_t AliVertexingHFUtils::GetCorrectedNtracklets(TProfile* estimatorAvg, Dou
 
   return correctedNacc;
 }
-
+//______________________________________________________________________
+TString AliVertexingHFUtils::GetGenerator(Int_t label, AliAODMCHeader* header){
+   Int_t nsumpart=0;
+   TList *lh=header->GetCocktailHeaders();
+   Int_t nh=lh->GetEntries();
+   for(Int_t i=0;i<nh;i++){
+     AliGenEventHeader* gh=(AliGenEventHeader*)lh->At(i);
+     TString genname=gh->GetName();
+     Int_t npart=gh->NProduced();
+     if(label>=nsumpart && label<(nsumpart+npart)) return genname;
+     nsumpart+=npart;
+   }
+   TString empty="";
+   return empty;
+ }