]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - PWG3/hfe/AliHFEextraCuts.cxx
Cleanup the code. Fix memory leak. Now inherit from AliAnalysisTaskSE (Antoine, Phili...
[u/mrichter/AliRoot.git] / PWG3 / hfe / AliHFEextraCuts.cxx
index aa3a6add6caae9de77774d6b2040dfebec0af915..ff24116f4a075569d5fc859e2a5c87d8c2c9ab59 100644 (file)
@@ -12,6 +12,9 @@
 * about the suitability of this software for any purpose. It is          *
 * provided "as is" without express or implied warranty.                  *
 **************************************************************************/
+
+/* $Id$ */
+
 //
 // Extra cuts implemented by the ALICE Heavy Flavour Electron Group
 // Cuts stored here:
@@ -22,6 +25,7 @@
 // Authors:
 //   Markus Fasel <M.Fasel@gsi.de>
 //
+#include <TBits.h>
 #include <TClass.h>
 #include <TH1F.h>
 #include <TH2F.h>
 #include <TString.h>
 #include <TMath.h>
 
+#include "AliAODTrack.h"
+#include "AliAODPid.h"
+#include "AliESDEvent.h"
+#include "AliESDVertex.h"
 #include "AliESDtrack.h"
 #include "AliLog.h"
 #include "AliMCParticle.h"
+#include "AliVEvent.h"
+#include "AliVTrack.h"
+#include "AliVParticle.h"
+#include "AliVertexerTracks.h"
+#include "AliVVertex.h"
 
 #include "AliHFEextraCuts.h"
 
@@ -40,11 +53,16 @@ ClassImp(AliHFEextraCuts)
 //______________________________________________________
 AliHFEextraCuts::AliHFEextraCuts(const Char_t *name, const Char_t *title):
   AliCFCutBase(name, title),
+  fEvent(NULL),
   fCutCorrelation(0),
   fRequirements(0),
+  fMinNClustersTPC(0),
   fClusterRatioTPC(0.),
   fMinTrackletsTRD(0),
   fPixelITS(0),
+  fTOFpid(kFALSE),
+  fTPCclusterDef(0),
+  fTPCclusterRatioDef(0),
   fCheck(kFALSE),
   fQAlist(0x0) ,
   fDebugLevel(0)
@@ -58,11 +76,16 @@ AliHFEextraCuts::AliHFEextraCuts(const Char_t *name, const Char_t *title):
 //______________________________________________________
 AliHFEextraCuts::AliHFEextraCuts(const AliHFEextraCuts &c):
   AliCFCutBase(c),
+  fEvent(c.fEvent),
   fCutCorrelation(c.fCutCorrelation),
   fRequirements(c.fRequirements),
+  fMinNClustersTPC(c.fMinNClustersTPC),
   fClusterRatioTPC(c.fClusterRatioTPC),
   fMinTrackletsTRD(c.fMinTrackletsTRD),
   fPixelITS(c.fPixelITS),
+  fTOFpid(c.fTOFpid),
+  fTPCclusterDef(c.fTPCclusterDef),
+  fTPCclusterRatioDef(c.fTPCclusterRatioDef),
   fCheck(c.fCheck),
   fQAlist(0x0),
   fDebugLevel(0)
@@ -75,7 +98,7 @@ AliHFEextraCuts::AliHFEextraCuts(const AliHFEextraCuts &c):
   if(IsQAOn()){
     fIsQAOn = kTRUE;
     fQAlist = dynamic_cast<TList *>(c.fQAlist->Clone());
-    fQAlist->SetOwner();
+    if(fQAlist) fQAlist->SetOwner();
   }
 }
 
@@ -86,18 +109,23 @@ AliHFEextraCuts &AliHFEextraCuts::operator=(const AliHFEextraCuts &c){
   //
   if(this != &c){
     AliCFCutBase::operator=(c);
+    fEvent = c.fEvent;
     fCutCorrelation = c.fCutCorrelation;
     fRequirements = c.fRequirements;
     fClusterRatioTPC = c.fClusterRatioTPC;
+    fMinNClustersTPC = c.fMinNClustersTPC;
     fMinTrackletsTRD = c.fMinTrackletsTRD;
     fPixelITS = c.fPixelITS;
+    fTPCclusterDef = c.fTPCclusterDef;
+    fTPCclusterRatioDef = c.fTPCclusterRatioDef;
+    fTOFpid = c.fTOFpid;
     fCheck = c.fCheck;
     fDebugLevel = c.fDebugLevel;
     memcpy(fImpactParamCut, c.fImpactParamCut, sizeof(Float_t) * 4);
     if(IsQAOn()){
       fIsQAOn = kTRUE;
       fQAlist = dynamic_cast<TList *>(c.fQAlist->Clone());
-      fQAlist->SetOwner();
+      if(fQAlist) fQAlist->SetOwner();
     }else fQAlist = 0x0;
   }
   return *this;
@@ -108,10 +136,25 @@ AliHFEextraCuts::~AliHFEextraCuts(){
   //
   // Destructor
   //
-  if(fQAlist){
-    fQAlist->Delete();
-    delete fQAlist;
+  if(fQAlist) delete fQAlist;
+}
+
+//______________________________________________________
+void AliHFEextraCuts::SetRecEventInfo(const TObject *event){
+  //
+  // Set Virtual event an make a copy
+  //
+  if (!event) {
+    AliError("Pointer to AliVEvent !");
+    return;
+  }
+  TString className(event->ClassName());
+  if (! (className.CompareTo("AliESDEvent")==0 || className.CompareTo("AliAODEvent")==0)) {
+    AliError("argument must point to an AliESDEvent or AliAODEvent !");
+    return ;
   }
+  fEvent = (AliVEvent*) event;
+
 }
 
 //______________________________________________________
@@ -119,14 +162,16 @@ Bool_t AliHFEextraCuts::IsSelected(TObject *o){
   //
   // Steering function for the track selection
   //
-  if(TString(o->IsA()->GetName()).CompareTo("AliESDtrack") == 0){
-    return CheckESDCuts(dynamic_cast<AliESDtrack *>(o));
+  TString type = o->IsA()->GetName();
+  AliDebug(2, Form("Object type %s", type.Data()));
+  if(!type.CompareTo("AliESDtrack") || !type.CompareTo("AliAODTrack")){
+    return CheckRecCuts(dynamic_cast<AliVTrack *>(o));
   }
-  return CheckMCCuts(dynamic_cast<AliMCParticle *>(o));
+  return CheckMCCuts(dynamic_cast<AliVParticle *>(o));
 }
 
 //______________________________________________________
-Bool_t AliHFEextraCuts::CheckESDCuts(AliESDtrack *track){
+Bool_t AliHFEextraCuts::CheckRecCuts(AliVTrack *track){
   //
   // Checks cuts on reconstructed tracks
   // returns true if track is selected
@@ -135,19 +180,26 @@ Bool_t AliHFEextraCuts::CheckESDCuts(AliESDtrack *track){
   //
   AliDebug(1, "Called");
   ULong64_t survivedCut = 0;   // Bitmap for cuts which are passed by the track, later to be compared with fRequirements
-  if(IsQAOn()) FillQAhistosESD(track, kBeforeCuts);
+  if(IsQAOn()) FillQAhistosRec(track, kBeforeCuts);
   // Apply cuts
-  Float_t impactR, impactZ, ratioTPC;
-  track->GetImpactParameters(impactR, impactZ);
+  Float_t impactR, impactZ;
+  Double_t hfeimpactR, hfeimpactnsigmaR;
+  Double_t hfeimpactRcut, hfeimpactnsigmaRcut;
+  Bool_t tofstep = kTRUE;
+  GetImpactParameters(track, impactR, impactZ);
+  if(TESTBIT(fRequirements, kMinHFEImpactParamR) || TESTBIT(fRequirements, kMinHFEImpactParamNsigmaR)){
+    // Protection for PbPb
+    GetHFEImpactParameterCuts(track, hfeimpactRcut, hfeimpactnsigmaRcut);
+    GetHFEImpactParameters(track, hfeimpactR, hfeimpactnsigmaR);
+  }
+  UInt_t nclsTPC = GetTPCncls(track);
   // printf("Check TPC findable clusters: %d, found Clusters: %d\n", track->GetTPCNclsF(), track->GetTPCNcls());
-  ratioTPC = track->GetTPCNclsF() > 0. ? static_cast<Float_t>(track->GetTPCNcls())/static_cast<Float_t>(track->GetTPCNclsF()) : 1.;
+  Double_t ratioTPC = GetTPCclusterRatio(track);
   UChar_t trdTracklets;
-  trdTracklets = track->GetTRDntrackletsPID();
+  trdTracklets = GetTRDnTrackletsPID(track);
   UChar_t itsPixel = track->GetITSClusterMap();
-  Int_t det, status1, status2;
-  Float_t xloc, zloc;
-  track->GetITSModuleIndexInfo(0, det, status1, xloc, zloc);
-  track->GetITSModuleIndexInfo(1, det, status2, xloc, zloc);
+  Int_t status1 = GetITSstatus(track, 0);
+  Int_t status2 = GetITSstatus(track, 1);
   Bool_t statusL0 = CheckITSstatus(status1);
   Bool_t statusL1 = CheckITSstatus(status2);
   if(TESTBIT(fRequirements, kMinImpactParamR)){
@@ -166,6 +218,14 @@ Bool_t AliHFEextraCuts::CheckESDCuts(AliESDtrack *track){
     // cut on max. Impact Parameter in Z direction
     if(TMath::Abs(impactZ) <= fImpactParamCut[3]) SETBIT(survivedCut, kMaxImpactParamZ);
   }
+  if(TESTBIT(fRequirements, kMinHFEImpactParamR)){
+    // cut on min. HFE Impact Parameter in Radial direction
+    if(TMath::Abs(hfeimpactR) >= hfeimpactRcut) SETBIT(survivedCut, kMinHFEImpactParamR);
+  }
+  if(TESTBIT(fRequirements, kMinHFEImpactParamNsigmaR)){
+    // cut on max. HFE Impact Parameter n sigma in Radial direction
+    if(TMath::Abs(hfeimpactnsigmaR) >= hfeimpactnsigmaRcut) SETBIT(survivedCut, kMinHFEImpactParamNsigmaR);
+  }
   if(TESTBIT(fRequirements, kClusterRatioTPC)){
     // cut on min ratio of found TPC clusters vs findable TPC clusters
     if(ratioTPC >= fClusterRatioTPC) SETBIT(survivedCut, kClusterRatioTPC);
@@ -175,6 +235,11 @@ Bool_t AliHFEextraCuts::CheckESDCuts(AliESDtrack *track){
     AliDebug(1, Form("Min TRD cut: [%d|%d]\n", fMinTrackletsTRD, trdTracklets));
     if(trdTracklets >= fMinTrackletsTRD) SETBIT(survivedCut, kMinTrackletsTRD);
   }
+  if(TESTBIT(fRequirements, kMinNClustersTPC)){
+    // cut on minimum number of TRD tracklets
+    AliDebug(1, Form("Min TPC cut: [%d|%d]\n", fMinNClustersTPC, nclsTPC));
+    if(nclsTPC >= fMinNClustersTPC) SETBIT(survivedCut, kMinNClustersTPC);
+  }
   if(TESTBIT(fRequirements, kPixelITS)){
     // cut on ITS pixel layers
     AliDebug(1, "ITS cluster Map: ");
@@ -218,19 +283,27 @@ Bool_t AliHFEextraCuts::CheckESDCuts(AliESDtrack *track){
     }
     AliDebug(1, Form("Survived Cut: %s\n", TESTBIT(survivedCut, kPixelITS) ? "YES" : "NO"));
   }
-  if(fRequirements == survivedCut){
+  if(fTOFpid){
+    // cut on TOF pid
+    tofstep = kFALSE;
+    if(track->GetStatus() & AliESDtrack::kTOFpid) tofstep = kTRUE;
+  
+  }
+  if(fRequirements == survivedCut && tofstep){
     //
     // Track selected
     //
-    if(IsQAOn()) FillQAhistosESD(track, kAfterCuts);
+    AliDebug(2, "Track Survived cuts\n");
+    if(IsQAOn()) FillQAhistosRec(track, kAfterCuts);
     return kTRUE;
   }
+  AliDebug(2, "Track cut");
   if(IsQAOn()) FillCutCorrelation(survivedCut);
   return kFALSE;
 }
 
 //______________________________________________________
-Bool_t AliHFEextraCuts::CheckMCCuts(AliMCParticle */*track*/) const {
+Bool_t AliHFEextraCuts::CheckMCCuts(AliVParticle */*track*/) const {
   //
   // Checks cuts on Monte Carlo tracks
   // returns true if track is selected
@@ -241,33 +314,38 @@ Bool_t AliHFEextraCuts::CheckMCCuts(AliMCParticle */*track*/) const {
 }
 
 //______________________________________________________
-void AliHFEextraCuts::FillQAhistosESD(AliESDtrack *track, UInt_t when){
+void AliHFEextraCuts::FillQAhistosRec(AliVTrack *track, UInt_t when){
   //
   // Fill the QA histograms for ESD tracks
   // Function can be called before cuts or after cut application (second argument)
   //
-  TList *container = dynamic_cast<TList *>(fQAlist->At(when));
+  const Int_t kNhistos = 6;
   Float_t impactR, impactZ;
-  track->GetImpactParameters(impactR, impactZ);
-  (dynamic_cast<TH1F *>(container->At(0)))->Fill(impactR);
-  (dynamic_cast<TH1F *>(container->At(1)))->Fill(impactZ);
+  GetImpactParameters(track, impactR, impactZ);
+  TH1 *htmp = NULL;
+  if((htmp = dynamic_cast<TH1F *>(fQAlist->At(0 + when * kNhistos)))) htmp->Fill(impactR);
+  if((htmp = dynamic_cast<TH1F *>(fQAlist->At(1 + when * kNhistos)))) htmp->Fill(impactZ);
   // printf("TPC findable clusters: %d, found Clusters: %d\n", track->GetTPCNclsF(), track->GetTPCNcls());
-  (dynamic_cast<TH1F *>(container->At(2)))->Fill(track->GetTPCNclsF() > 0. ? static_cast<Float_t>(track->GetTPCNcls())/static_cast<Float_t>(track->GetTPCNclsF()) : 1.);
-  (dynamic_cast<TH1F *>(container->At(3)))->Fill(track->GetTRDntrackletsPID());
+  if((htmp = dynamic_cast<TH1F *>(fQAlist->At(2 + when * kNhistos)))) htmp->Fill(GetTPCclusterRatio(track));
+  if((htmp = dynamic_cast<TH1F *>(fQAlist->At(3 + when * kNhistos)))) htmp->Fill(GetTRDnTrackletsPID(track));
+  if((htmp = dynamic_cast<TH1F *>(fQAlist->At(4 + when * kNhistos)))) htmp->Fill(GetTPCncls(track));
   UChar_t itsPixel = track->GetITSClusterMap();
-  TH1 *pixelHist = dynamic_cast<TH1F *>(container->At(4));
-  Int_t firstEntry = pixelHist->GetXaxis()->GetFirst();
-  if(!((itsPixel & BIT(0)) || (itsPixel & BIT(1))))
-    pixelHist->Fill(firstEntry + 3);
-  else{
-    if(itsPixel & BIT(0)){
-      pixelHist->Fill(firstEntry);
-      if(itsPixel & BIT(1)) pixelHist->Fill(firstEntry + 2);
-      else pixelHist->Fill(firstEntry + 4);
-    }
-    if(itsPixel & BIT(1)){
-      pixelHist->Fill(firstEntry + 1);
-      if(!(itsPixel & BIT(0))) pixelHist->Fill(firstEntry + 5);
+  TH1 *pixelHist = dynamic_cast<TH1F *>(fQAlist->At(5 + when * kNhistos));
+  //Int_t firstEntry = pixelHist->GetXaxis()->GetFirst();
+  if(pixelHist){
+    Double_t firstEntry = 0.5;
+    if(!((itsPixel & BIT(0)) || (itsPixel & BIT(1))))
+      pixelHist->Fill(firstEntry + 3);
+    else{
+      if(itsPixel & BIT(0)){
+        pixelHist->Fill(firstEntry);
+        if(itsPixel & BIT(1)) pixelHist->Fill(firstEntry + 2);
+        else pixelHist->Fill(firstEntry + 4);
+      }
+      if(itsPixel & BIT(1)){
+        pixelHist->Fill(firstEntry + 1);
+        if(!(itsPixel & BIT(0))) pixelHist->Fill(firstEntry + 5);
+      }
     }
   }
 }
@@ -286,13 +364,16 @@ void AliHFEextraCuts::FillCutCorrelation(ULong64_t survivedCut){
   //
   // Fill cut correlation histograms for tracks that didn't pass cuts
   //
-  TH2 *correlation = dynamic_cast<TH2F *>(fQAlist->At(2));
-  for(Int_t icut = 0; icut < kNcuts; icut++){
-    if(!TESTBIT(fRequirements, icut)) continue;
-    for(Int_t jcut = icut; jcut < kNcuts; jcut++){
-      if(!TESTBIT(fRequirements, jcut)) continue;
-      if(TESTBIT(survivedCut, icut) && TESTBIT(survivedCut, jcut))
-       correlation->Fill(icut, jcut);
+  const Int_t kNhistos = 6;
+  TH2 *correlation = dynamic_cast<TH2F *>(fQAlist->At(2 * kNhistos));
+  if(correlation){
+    for(Int_t icut = 0; icut < kNcuts; icut++){
+      if(!TESTBIT(fRequirements, icut)) continue;
+      for(Int_t jcut = icut; jcut < kNcuts; jcut++){
+        if(!TESTBIT(fRequirements, jcut)) continue;
+        if(TESTBIT(survivedCut, icut) && TESTBIT(survivedCut, jcut))
+               correlation->Fill(icut, jcut);
+      }
     }
   }
 }
@@ -306,30 +387,36 @@ void AliHFEextraCuts::AddQAHistograms(TList *qaList){
   // Additionally a histogram with the cut correlation is created and stored
   // in the top directory
   //
-  TList *histos[2];
+
+  const Int_t kNhistos = 6;
   TH1 *histo1D = 0x0;
   TH2 *histo2D = 0x0;
-  histos[0] = new TList();
-  histos[0]->SetName("BeforeCut");
-  histos[0]->SetOwner();
-  histos[1] = new TList();
-  histos[1]->SetName("AfterCut");
-  histos[1]->SetOwner();
   TString cutstr[2] = {"before", "after"};
+
+  if(!fQAlist) fQAlist = new TList;  // for internal representation, not owner
   for(Int_t icond = 0; icond < 2; icond++){
-    histos[icond]->AddAt((histo1D = new TH1F(Form("impactParamR%s", cutstr[icond].Data()), "Radial Impact Parameter", 100, 0, 10)), 0);
+    qaList->AddAt((histo1D = new TH1F(Form("%s_impactParamR%s",GetName(),cutstr[icond].Data()), "Radial Impact Parameter", 100, 0, 10)), 0 + icond * kNhistos);
+    fQAlist->AddAt(histo1D, 0 + icond * kNhistos);
     histo1D->GetXaxis()->SetTitle("Impact Parameter");
     histo1D->GetYaxis()->SetTitle("Number of Tracks");
-    histos[icond]->AddAt((histo1D = new TH1F(Form("impactParamZ%s", cutstr[icond].Data()), "Z Impact Parameter", 200, 0, 20)), 1);
+    qaList->AddAt((histo1D = new TH1F(Form("%s_impactParamZ%s",GetName(),cutstr[icond].Data()), "Z Impact Parameter", 200, 0, 20)), 1 + icond * kNhistos);
+    fQAlist->AddAt(histo1D, 1 + icond * kNhistos);
     histo1D->GetXaxis()->SetTitle("Impact Parameter");
     histo1D->GetYaxis()->SetTitle("Number of Tracks");
-    histos[icond]->AddAt((histo1D = new TH1F(Form("tpcClr%s", cutstr[icond].Data()), "Cluster Ratio TPC", 10, 0, 1)), 2);
+    qaList->AddAt((histo1D = new TH1F(Form("%s_tpcClr%s",GetName(),cutstr[icond].Data()), "Cluster Ratio TPC", 100, 0, 1)), 2 + icond * kNhistos);
+    fQAlist->AddAt(histo1D, 2 + icond * kNhistos);
     histo1D->GetXaxis()->SetTitle("Cluster Ratio TPC");
     histo1D->GetYaxis()->SetTitle("Number of Tracks");
-    histos[icond]->AddAt((histo1D = new TH1F(Form("trdTracklets%s", cutstr[icond].Data()), "Number of TRD tracklets", 7, 0, 7)), 3);
+    qaList->AddAt((histo1D = new TH1F(Form("%s_trdTracklets%s",GetName(),cutstr[icond].Data()), "Number of TRD tracklets", 7, 0, 7)), 3 + icond * kNhistos);
+    fQAlist->AddAt(histo1D, 3 + icond * kNhistos);
     histo1D->GetXaxis()->SetTitle("Number of TRD Tracklets");
     histo1D->GetYaxis()->SetTitle("Number of Tracks");
-    histos[icond]->AddAt((histo1D = new TH1F(Form("itsPixel%s", cutstr[icond].Data()), "ITS Pixel Hits", 6, 0, 6)), 4);
+    qaList->AddAt((histo1D = new TH1F(Form("%s_tpcClusters%s",GetName(),cutstr[icond].Data()), "Number of TPC clusters", 161, 0, 160)), 4 + icond * kNhistos);
+    fQAlist->AddAt(histo1D, 4 + icond * kNhistos);
+    histo1D->GetXaxis()->SetTitle("Number of TPC clusters");
+    histo1D->GetYaxis()->SetTitle("Number of Tracks");
+    qaList->AddAt((histo1D = new TH1F(Form("%s_itsPixel%s",GetName(),cutstr[icond].Data()), "ITS Pixel Hits", 6, 0, 6)), 5 + icond * kNhistos);
+    fQAlist->AddAt(histo1D, 5 + icond * kNhistos);
     histo1D->GetXaxis()->SetTitle("ITS Pixel");
     histo1D->GetYaxis()->SetTitle("Number of Tracks");
     Int_t first = histo1D->GetXaxis()->GetFirst();
@@ -337,20 +424,15 @@ void AliHFEextraCuts::AddQAHistograms(TList *qaList){
     for(Int_t ilabel = 0; ilabel < 6; ilabel++)
       histo1D->GetXaxis()->SetBinLabel(first + ilabel, binNames[ilabel].Data());
   }
-  fQAlist = new TList();
-  fQAlist->SetOwner();
-  fQAlist->SetName("HFelectronExtraCuts");
-  fQAlist->AddAt(histos[0], 0);
-  fQAlist->AddAt(histos[1], 1);
   // Add cut correlation
-  fQAlist->AddAt((histo2D = new TH2F("cutcorrellation", "Cut Correlation", kNcuts, 0, kNcuts - 1, kNcuts, 0, kNcuts -1)), 2);
-  TString labels[kNcuts] = {"MinImpactParamR", "MaxImpactParamR", "MinImpactParamZ", "MaxImpactParamZ", "ClusterRatioTPC", "MinTrackletsTRD", "ITSpixel"};
+  qaList->AddAt((histo2D = new TH2F(Form("%s_cutcorrelation",GetName()), "Cut Correlation", kNcuts, 0, kNcuts - 1, kNcuts, 0, kNcuts -1)), 2 * kNhistos);
+  fQAlist->AddAt(histo2D, 2 * kNhistos);
+  TString labels[kNcuts] = {"MinImpactParamR", "MaxImpactParamR", "MinImpactParamZ", "MaxImpactParamZ", "ClusterRatioTPC", "MinTrackletsTRD", "ITSpixel", "kMinHFEImpactParamR", "kMinHFEImpactParamNsigmaR", "TPC Number of clusters"};
   Int_t firstx = histo2D->GetXaxis()->GetFirst(), firsty = histo2D->GetYaxis()->GetFirst();
   for(Int_t icut = 0; icut < kNcuts; icut++){
     histo2D->GetXaxis()->SetBinLabel(firstx + icut, labels[icut].Data());
     histo2D->GetYaxis()->SetBinLabel(firsty + icut, labels[icut].Data());
   }
-  qaList->AddLast(fQAlist);
 }
 
 //______________________________________________________
@@ -374,3 +456,193 @@ Bool_t AliHFEextraCuts::CheckITSstatus(Int_t itsStatus) const {
   }
   return status;
 }
+
+//______________________________________________________
+Int_t AliHFEextraCuts::GetTRDnTrackletsPID(AliVTrack *track){
+       //
+       // Get Number of TRD tracklets
+       //
+       Int_t nTracklets = 0;
+       if(!TString(track->IsA()->GetName()).CompareTo("AliESDtrack")){
+               AliESDtrack *esdtrack = dynamic_cast<AliESDtrack *>(track);
+               if(esdtrack) nTracklets = esdtrack->GetTRDntrackletsPID();
+       } else if(!TString(track->IsA()->GetName()).CompareTo("AliAODTrack")){
+               AliAODTrack *aodtrack = dynamic_cast<AliAODTrack *>(track);
+               AliAODPid *pidobject = NULL;
+    if(aodtrack) pidobject = aodtrack->GetDetPid();
+               // this is normally NOT the way to do this, but due to limitation in the
+               // AOD track it is not possible in a different way
+               if(pidobject){
+                       Float_t *trdmom = pidobject->GetTRDmomentum();
+                       for(Int_t ily = 0; ily < 6; ily++){
+                               if(trdmom[ily] > -1) nTracklets++;
+                       }
+               } else nTracklets = 6;  // No Cut possible
+       }
+       return nTracklets;
+}
+
+//______________________________________________________
+Int_t AliHFEextraCuts::GetITSstatus(AliVTrack *track, Int_t layer){
+       //
+       // Check ITS layer status
+       //
+       Int_t status = 0;
+       if(!TString(track->IsA()->GetName()).CompareTo("AliESDtrack")){
+               Int_t det;
+               Float_t xloc, zloc;
+               AliESDtrack *esdtrack = dynamic_cast<AliESDtrack *>(track);
+               if(esdtrack) esdtrack->GetITSModuleIndexInfo(layer, det, status, xloc, zloc);
+       }
+       return status;
+}
+
+//______________________________________________________
+UInt_t AliHFEextraCuts::GetTPCncls(AliVTrack *track){
+       //
+       // Get Number of findable clusters in the TPC
+       //
+       Int_t nClusters = 0; // in case no Information available consider all clusters findable
+       TString type = track->IsA()->GetName();
+       if(!type.CompareTo("AliESDtrack")){
+    AliESDtrack *esdtrack = dynamic_cast<AliESDtrack *>(track);
+    if(esdtrack){ // coverity
+      if(TESTBIT(fTPCclusterDef, kFoundIter1)){
+                   nClusters = esdtrack->GetTPCNclsIter1();
+        AliDebug(2, ("Using def kFoundIter1"));
+      } else if(TESTBIT(fTPCclusterDef, kCrossedRows)){
+        AliDebug(2, ("Using def kCrossedRows"));
+        nClusters = static_cast<UInt_t>(esdtrack->GetTPCClusterInfo(2,1));
+      } else{
+        AliDebug(2, ("Using def kFound"));
+                   nClusters = esdtrack->GetTPCNcls();
+      }
+    }
+  }
+       else if(!type.CompareTo("AliAODTrack")){
+               AliAODTrack *aodtrack = dynamic_cast<AliAODTrack *>(track);
+    if(aodtrack){
+           const TBits &tpcmap = aodtrack->GetTPCClusterMap();
+           for(UInt_t ibit = 0; ibit < tpcmap.GetNbits(); ibit++)
+                   if(tpcmap.TestBitNumber(ibit)) nClusters++;
+    }
+       }
+       return nClusters;
+}
+
+//______________________________________________________
+Double_t AliHFEextraCuts::GetTPCclusterRatio(AliVTrack *track){
+  //
+  // Get Ratio of found / findable clusters for different definitions
+  // Only implemented for ESD tracks
+  //
+  Double_t clusterRatio = 1.; // in case no Information available consider all clusters findable
+  TString type = track->IsA()->GetName();
+  if(!type.CompareTo("AliESDtrack")){
+    AliESDtrack *esdtrack = dynamic_cast<AliESDtrack *>(track);
+    if(esdtrack){ // coverity
+      if(TESTBIT(fTPCclusterRatioDef, kCROverFindable)){
+        AliDebug(2, "Using ratio def kCROverFindable");
+        clusterRatio = esdtrack->GetTPCNclsF() ? esdtrack->GetTPCClusterInfo(2,1)/static_cast<Double_t>(esdtrack->GetTPCNclsF()) : 1.; // crossed rows/findable
+      } else if(TESTBIT(fTPCclusterRatioDef, kFoundOverCR)){
+        AliDebug(2, "Using ratio def kFoundOverCR");
+        clusterRatio = esdtrack->GetTPCClusterInfo(2,0); // found/crossed rows
+      } else if(TESTBIT(fTPCclusterRatioDef, kFoundOverFindableIter1)){
+        AliDebug(2, "Using ratio def kFoundOverFindableIter1");
+        clusterRatio = esdtrack->GetTPCNclsFIter1() ? static_cast<Double_t>(esdtrack->GetTPCNclsIter1())/static_cast<Double_t>(esdtrack->GetTPCNclsFIter1()) : 1.; // crossed
+      } else {
+        AliDebug(2, "Using ratio def kFoundOverFindable");
+        clusterRatio = esdtrack->GetTPCNclsF() ? static_cast<Double_t>(esdtrack->GetTPCNcls())/static_cast<Double_t>(esdtrack->GetTPCNclsF()) : 1.; // found/findable
+      }
+    }
+  }
+  return clusterRatio;
+}
+
+//______________________________________________________
+void AliHFEextraCuts::GetImpactParameters(AliVTrack *track, Float_t &radial, Float_t &z){
+       //
+       // Get impact parameter
+       //
+       TString type = track->IsA()->GetName();
+       if(!type.CompareTo("AliESDtrack")){
+               AliESDtrack *esdtrack = dynamic_cast<AliESDtrack *>(track);
+    if(esdtrack) esdtrack->GetImpactParameters(radial, z);
+       }
+       else if(!type.CompareTo("AliAODTrack")){
+               AliAODTrack *aodtrack = dynamic_cast<AliAODTrack *>(track);
+    if(aodtrack){
+                 Double_t xyz[3];
+                 aodtrack->XYZAtDCA(xyz);
+                 z = xyz[2];
+                 radial = TMath::Sqrt(xyz[0]*xyz[0] + xyz[1]+xyz[1]);
+    }
+       }
+}
+
+//______________________________________________________
+void AliHFEextraCuts::GetHFEImpactParameters(AliVTrack *track, Double_t &dcaxy, Double_t &dcansigmaxy){
+       //
+       // Get HFE impact parameter (with recalculated primary vertex)
+       //
+       dcaxy=0;
+       dcansigmaxy=0;
+  if(!fEvent){
+    AliDebug(1, "No Input event available\n");
+    return;
+  }
+  const Double_t kBeampiperadius=3.;
+  TString type = track->IsA()->GetName();
+  Double_t dca[2]={-999.,-999.};
+  Double_t cov[3]={-999.,-999.,-999.};
+
+  // recalculate primary vertex
+  AliVertexerTracks vertexer(fEvent->GetMagneticField());
+  vertexer.SetITSMode();
+  vertexer.SetMinClusters(4);
+       Int_t skipped[2];
+  skipped[0] = track->GetID();
+  vertexer.SetSkipTracks(1,skipped);
+  AliVVertex *vtxESDSkip = vertexer.FindPrimaryVertex(fEvent);
+  vertexer.SetSkipTracks(1,skipped);
+  if(vtxESDSkip->GetNContributors()<2) return;
+
+  // Getting the DCA
+  // Propagation always done on a working copy to not disturb the track params of the original track
+  AliESDtrack *esdtrack = NULL;
+  if(!TString(track->IsA()->GetName()).CompareTo("AliESDtrack")){
+    // Case ESD track: take copy constructor
+    AliESDtrack *tmptrack = dynamic_cast<AliESDtrack *>(track);
+    if(tmptrack) esdtrack = new AliESDtrack(*tmptrack);
+  } else {
+    // Case AOD track: take different constructor
+    esdtrack = new AliESDtrack(track);
+  }
+  if(esdtrack && esdtrack->PropagateToDCA(vtxESDSkip, fEvent->GetMagneticField(), kBeampiperadius, dca, cov)){
+    // protection
+    dcaxy = dca[0];
+    if(cov[0]) dcansigmaxy = dcaxy/TMath::Sqrt(cov[0]);
+    if(!cov[0]) dcansigmaxy = -99.;
+  }
+  delete esdtrack;
+  delete vtxESDSkip;
+}
+
+
+//______________________________________________________
+void AliHFEextraCuts::GetHFEImpactParameterCuts(AliVTrack *track, Double_t &hfeimpactRcut, Double_t &hfeimpactnsigmaRcut){
+       //
+       // Get HFE impact parameter cut(pt dependent)
+       //
+  
+  TString type = track->IsA()->GetName();
+  if(!type.CompareTo("AliESDtrack")){
+    AliESDtrack *esdtrack = dynamic_cast<AliESDtrack *>(track);
+    if(!esdtrack) return;
+    Double_t pt = esdtrack->Pt();      
+    //hfeimpactRcut=0.0064+0.078*exp(-0.56*pt);  // used Carlo's old parameter 
+    hfeimpactRcut=0.011+0.077*exp(-0.65*pt); // used Carlo's new parameter
+    hfeimpactnsigmaRcut=3; // 3 sigma trail cut
+  }
+}
+