]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - PWG3/vertexingHF/AliRDHFCuts.cxx
Skip tracks without ITS
[u/mrichter/AliRoot.git] / PWG3 / vertexingHF / AliRDHFCuts.cxx
index 8e8f000a00ab8869cd726e92d129daa6b448b017..e11577d52b36ce4fe3e87c1cc425545dd5d20d7c 100644 (file)
 #include "AliAODEvent.h"
 #include "AliVVertex.h"
 #include "AliESDVertex.h"
+#include "AliLog.h"
 #include "AliAODVertex.h"
 #include "AliESDtrack.h"
 #include "AliAODTrack.h"
 #include "AliESDtrackCuts.h"
+#include "AliCentrality.h"
 #include "AliAODRecoDecayHF.h"
 #include "AliRDHFCuts.h"
 
@@ -41,6 +43,7 @@ AliAnalysisCuts(name,title),
 fMinVtxType(3),
 fMinVtxContr(1),
 fMaxVtxRedChi2(1e6),
+fMaxVtxZ(1e6),
 fMinSPDMultiplicity(0),
 fTriggerMask(0),
 fTrackCuts(0),
@@ -57,7 +60,14 @@ fIsUpperCut(0),
 fUsePID(kFALSE),
 fPidHF(0),
 fWhyRejection(0),
-fRemoveDaughtersFromPrimary(kFALSE)
+fRemoveDaughtersFromPrimary(kFALSE),
+fOptPileup(0),
+fMinContrPileup(3),
+fMinDzPileup(0.6),
+fUseCentrality(0),
+fMinCentrality(0.),
+fMaxCentrality(100.)
+
 {
   //
   // Default Constructor
@@ -69,6 +79,7 @@ AliRDHFCuts::AliRDHFCuts(const AliRDHFCuts &source) :
   fMinVtxType(source.fMinVtxType),
   fMinVtxContr(source.fMinVtxContr),
   fMaxVtxRedChi2(source.fMaxVtxRedChi2),
+  fMaxVtxZ(source.fMaxVtxZ),
   fMinSPDMultiplicity(source.fMinSPDMultiplicity),
   fTriggerMask(source.fTriggerMask),
   fTrackCuts(0),
@@ -85,7 +96,15 @@ AliRDHFCuts::AliRDHFCuts(const AliRDHFCuts &source) :
   fUsePID(source.fUsePID),
   fPidHF(0),
   fWhyRejection(source.fWhyRejection),
-  fRemoveDaughtersFromPrimary(source.fRemoveDaughtersFromPrimary)
+  fRemoveDaughtersFromPrimary(source.fRemoveDaughtersFromPrimary),
+  fOptPileup(source.fOptPileup),
+  fMinContrPileup(source.fMinContrPileup),
+  fMinDzPileup(source.fMinDzPileup),
+  fUseCentrality(source.fUseCentrality),
+  fMinCentrality(source.fMinCentrality),
+  fMaxCentrality(source.fMaxCentrality)
+
+
 {
   //
   // Copy constructor
@@ -113,6 +132,7 @@ AliRDHFCuts &AliRDHFCuts::operator=(const AliRDHFCuts &source)
   fMinVtxType=source.fMinVtxType;
   fMinVtxContr=source.fMinVtxContr;
   fMaxVtxRedChi2=source.fMaxVtxRedChi2;
+  fMaxVtxZ=source.fMaxVtxZ;
   fMinSPDMultiplicity=source.fMinSPDMultiplicity;
   fTriggerMask=source.fTriggerMask;
   fnPtBins=source.fnPtBins;
@@ -123,6 +143,12 @@ AliRDHFCuts &AliRDHFCuts::operator=(const AliRDHFCuts &source)
   SetPidHF(source.GetPidHF());
   fWhyRejection=source.fWhyRejection;
   fRemoveDaughtersFromPrimary=source.fRemoveDaughtersFromPrimary;
+  fOptPileup=source.fOptPileup;
+  fMinContrPileup=source.fMinContrPileup;
+  fMinDzPileup=source.fMinDzPileup;
+  fUseCentrality=source.fUseCentrality;
+  fMinCentrality=source.fMinCentrality;
+  fMaxCentrality=source.fMaxCentrality;
 
   if(source.GetTrackCuts()) AddTrackCuts(source.GetTrackCuts());
   if(source.fPtBinLimits) SetPtBins(source.fnPtBinLimits,source.fPtBinLimits);
@@ -152,14 +178,16 @@ AliRDHFCuts::~AliRDHFCuts() {
   }
 }
 //---------------------------------------------------------------------------
-Bool_t AliRDHFCuts::IsEventSelected(AliVEvent *event) const {
+Bool_t AliRDHFCuts::IsEventSelected(AliVEvent *event) {
   //
   // Event selection
   // 
   //if(fTriggerMask && event->GetTriggerMask()!=fTriggerMask) return kFALSE;
 
-  // multiplicity cuts no implemented yet
+  fWhyRejection=0;
 
+  // multiplicity cuts no implemented yet
+   
   const AliVVertex *vertex = event->GetPrimaryVertex();
 
   if(!vertex) return kFALSE;
@@ -170,6 +198,48 @@ Bool_t AliRDHFCuts::IsEventSelected(AliVEvent *event) const {
 
   if(vertex->GetNContributors()<fMinVtxContr) return kFALSE; 
 
+  if(TMath::Abs(vertex->GetZ())>fMaxVtxZ) return kFALSE;
+
+  // switch to settings for 1-pad cls in TPC
+  if(fPidHF) {
+    if(event->GetRunNumber()>121693 && event->GetRunNumber()<136851) 
+      fPidHF->SetOnePad(kTRUE);
+    if(event->GetRunNumber()>=136851 && event->GetRunNumber()<=139517) 
+      fPidHF->SetPbPb(kTRUE);
+  }
+
+  if(fOptPileup==kRejectPileupEvent){
+    Int_t cutc=(Int_t)fMinContrPileup;
+    Double_t cutz=(Double_t)fMinDzPileup;
+    if(event->IsPileupFromSPD(cutc,cutz,3.,2.,10.)) {
+      fWhyRejection=1;
+      return kFALSE;
+    }
+  }
+
+  //centrality selection
+  if (!(fUseCentrality==kCentOff)){  
+    if(fUseCentrality<kCentOff||fUseCentrality>=kCentInvalid){
+      AliWarning("Centrality estimator not valid");
+      fWhyRejection=3;
+      return kFALSE;
+    }else{
+      Float_t centvalue=GetCentrality((AliAODEvent*)event);
+      if (centvalue<0.){
+       if (fWhyRejection==3) return kFALSE;
+       else return kTRUE;
+      }
+      else{
+       
+       if (centvalue<fMinCentrality || centvalue>fMaxCentrality){
+         fWhyRejection=2; 
+         return kFALSE; 
+       }
+      }
+    }
+  }
+
+
   return kTRUE;
 }
 //---------------------------------------------------------------------------
@@ -214,6 +284,11 @@ Bool_t AliRDHFCuts::IsDaughterSelected(AliAODTrack *track,const AliESDVertex *pr
   esdTrack.RelateToVertex(primary,0.,3.); 
   if(!cuts->IsSelected(&esdTrack)) retval = kFALSE;
  
+  if(fOptPileup==kRejectTracksFromPileupVertex){
+    // to be implemented
+    // we need either to have here the AOD Event, 
+    // or to have the pileup vertex object
+  }
   return retval; 
 }
 //---------------------------------------------------------------------------
@@ -286,6 +361,19 @@ void AliRDHFCuts::SetVarsForOpt(Int_t nVars,Bool_t *forOpt) {
 
   return;
 }
+
+//---------------------------------------------------------------------------
+void AliRDHFCuts::SetUseCentrality(Int_t flag) {
+  //
+  // set centrality estimator  
+  //
+  fUseCentrality=flag;
+  if(fUseCentrality<kCentOff||fUseCentrality>=kCentInvalid) AliWarning("Centrality estimator not valid");
+  return;
+}
+
+
 //---------------------------------------------------------------------------
 void AliRDHFCuts::SetCuts(Int_t nVars,Int_t nPtBins,Float_t **cutsRD) {
   //
@@ -347,6 +435,17 @@ void AliRDHFCuts::PrintAll() const {
   printf("Min SPD mult %d\n",fMinSPDMultiplicity);
   printf("Use PID %d\n",(Int_t)fUsePID);
   printf("Remove daughters from vtx %d\n",(Int_t)fRemoveDaughtersFromPrimary);
+  printf("Pileup rejection: %s\n",(fOptPileup > 0) ? "Yes" : "No");
+  if(fOptPileup==1) printf(" -- Reject pileup event");
+  if(fOptPileup==2) printf(" -- Reject tracks from pileup vtx");
+  if(fUseCentrality>0) {
+    TString estimator="";
+    if(fUseCentrality==1) estimator = "V0";
+    if(fUseCentrality==2) estimator = "Tracks";
+    if(fUseCentrality==3) estimator = "Tracklets";
+    if(fUseCentrality==4) estimator = "SPD clusters outer"; 
+    printf("Centrality class considered: %.1f-%.1f, estimated with %s",fMinCentrality,fMaxCentrality,estimator.Data());
+  }
 
   if(fVarNames){
     cout<<"Array of variables"<<endl;
@@ -459,6 +558,33 @@ Float_t AliRDHFCuts::GetCutValue(Int_t iVar,Int_t iPtBin) const {
   return fCutsRD[GetGlobalIndex(iVar,iPtBin)];
 }
 //-------------------------------------------------------------------
+Float_t AliRDHFCuts::GetCentrality(AliAODEvent* aodEvent,AliRDHFCuts::ECentrality estimator) const {
+  //
+  // Get centrality percentile
+  //
+  AliAODHeader *header=aodEvent->GetHeader();
+  AliCentrality *centrality=header->GetCentralityP();
+  Float_t cent=-999.;
+  if(!centrality) return cent;
+  else{
+    if (estimator==kCentV0M) cent=(Float_t)(centrality->GetCentralityPercentile("V0M"));
+    else {
+      if (estimator==kCentTRK) cent=(Float_t)(centrality->GetCentralityPercentile("TRK"));
+      else{
+       if (estimator==kCentTKL) cent=(Float_t)(centrality->GetCentralityPercentile("TKL"));
+       else{
+         if (estimator==kCentCL1) cent=(Float_t)(centrality->GetCentralityPercentile("CL1"));
+         else {
+           AliWarning("Centrality estimator not valid");
+
+         }
+       }
+      }
+    } 
+  }
+  return cent;
+}
+//-------------------------------------------------------------------
 Bool_t AliRDHFCuts::CompareCuts(const AliRDHFCuts *obj) const {
   //
   // Compare two cuts objects
@@ -477,14 +603,15 @@ Bool_t AliRDHFCuts::CompareCuts(const AliRDHFCuts *obj) const {
   if(fUsePID!=obj->fUsePID) { printf("Use PID %d  %d\n",(Int_t)fUsePID,(Int_t)obj->fUsePID); areEqual=kFALSE;}
 
   if(fRemoveDaughtersFromPrimary!=obj->fRemoveDaughtersFromPrimary) {printf("Remove daughters from vtx %d  %d\n",(Int_t)fRemoveDaughtersFromPrimary,(Int_t)obj->fRemoveDaughtersFromPrimary); areEqual=kFALSE;}
+  if(fTrackCuts){
+    if(fTrackCuts->GetMinNClusterTPC()!=obj->fTrackCuts->GetMinNClusterTPC()) {printf("MinNClsTPC %d  %d\n",fTrackCuts->GetMinNClusterTPC(),obj->fTrackCuts->GetMinNClusterTPC()); areEqual=kFALSE;}
 
-  if(fTrackCuts->GetMinNClusterTPC()!=obj->fTrackCuts->GetMinNClusterTPC()) {printf("MinNClsTPC %d  %d\n",fTrackCuts->GetMinNClusterTPC(),obj->fTrackCuts->GetMinNClusterTPC()); areEqual=kFALSE;}
-
-  if(fTrackCuts->GetMinNClustersITS()!=obj->fTrackCuts->GetMinNClustersITS()) {printf("MinNClsITS %d  %d\n",fTrackCuts->GetMinNClustersITS(),obj->fTrackCuts->GetMinNClustersITS()); areEqual=kFALSE;}
+    if(fTrackCuts->GetMinNClustersITS()!=obj->fTrackCuts->GetMinNClustersITS()) {printf("MinNClsITS %d  %d\n",fTrackCuts->GetMinNClustersITS(),obj->fTrackCuts->GetMinNClustersITS()); areEqual=kFALSE;}
 
-  if(TMath::Abs(fTrackCuts->GetMaxChi2PerClusterTPC()-obj->fTrackCuts->GetMaxChi2PerClusterTPC())>1.e-10) {printf("MaxChi2ClsTPC %f  %f\n",fTrackCuts->GetMaxChi2PerClusterTPC(),obj->fTrackCuts->GetMaxChi2PerClusterTPC()); areEqual=kFALSE;}
+    if(TMath::Abs(fTrackCuts->GetMaxChi2PerClusterTPC()-obj->fTrackCuts->GetMaxChi2PerClusterTPC())>1.e-10) {printf("MaxChi2ClsTPC %f  %f\n",fTrackCuts->GetMaxChi2PerClusterTPC(),obj->fTrackCuts->GetMaxChi2PerClusterTPC()); areEqual=kFALSE;}
 
-  if(fTrackCuts->GetClusterRequirementITS(AliESDtrackCuts::kSPD)!=obj->fTrackCuts->GetClusterRequirementITS(AliESDtrackCuts::kSPD)) {printf("ClusterReq SPD %d  %d\n",fTrackCuts->GetClusterRequirementITS(AliESDtrackCuts::kSPD),obj->fTrackCuts->GetClusterRequirementITS(AliESDtrackCuts::kSPD)); areEqual=kFALSE;}
+    if(fTrackCuts->GetClusterRequirementITS(AliESDtrackCuts::kSPD)!=obj->fTrackCuts->GetClusterRequirementITS(AliESDtrackCuts::kSPD)) {printf("ClusterReq SPD %d  %d\n",fTrackCuts->GetClusterRequirementITS(AliESDtrackCuts::kSPD),obj->fTrackCuts->GetClusterRequirementITS(AliESDtrackCuts::kSPD)); areEqual=kFALSE;}
+  }
 
   if(fCutsRD) {
    for(Int_t iv=0;iv<fnVars;iv++) {
@@ -499,3 +626,39 @@ Bool_t AliRDHFCuts::CompareCuts(const AliRDHFCuts *obj) const {
 
   return areEqual;
 }
+//---------------------------------------------------------------------------
+void AliRDHFCuts::MakeTable() const {
+  //
+  // print cuts values in table format
+  // 
+
+       TString ptString = "pT range";
+       if(fVarNames && fPtBinLimits && fCutsRD){
+               TString firstLine(Form("*       %-15s",ptString.Data()));
+               for (Int_t ivar=0; ivar<fnVars; ivar++){
+                       firstLine+=Form("*    %-15s  ",fVarNames[ivar].Data());
+                       if (ivar == fnVars){
+                               firstLine+="*\n";
+                       }
+               }
+               Printf("%s",firstLine.Data());
+               
+               for (Int_t ipt=0; ipt<fnPtBins; ipt++){
+                       TString line;
+                       if (ipt==fnPtBins-1){
+                               line=Form("*  %5.1f < pt < inf    ",fPtBinLimits[ipt]);
+                       }
+                       else{
+                               line=Form("*  %5.1f < pt < %4.1f   ",fPtBinLimits[ipt],fPtBinLimits[ipt+1]);
+                       }
+                       for (Int_t ivar=0; ivar<fnVars; ivar++){
+                               line+=Form("*     %-15f ",fCutsRD[GetGlobalIndex(ivar,ipt)]);
+                       }
+                       Printf("%s",line.Data());
+               }
+
+       }
+
+
+  return;
+}