[PWGHF,Trigger] Add track quality selection in HF filter

EventFiltering/PWGHF/HFFilter.cxx

@@ -1,4 +1,4 @@
 // Copyright 2019-2020 CERN and copyright holders of ALICE O2.
 // See https://alice-o2.web.cern.ch/copyright for details of the copyright holders.
 // All rights not expressly granted are reserved.
 //
@@ -76,7 +76,7 @@
   Configurable<LabeledArray<float>> nSigmaPidCuts{"nSigmaPidCuts", {cutsNsigma[0], 4, 8, labelsRowsNsigma, labelsColumnsNsigma}, "Nsigma cuts for ITS/TPC/TOF PID (except for V0 and cascades)"};
   // min and max pts for tracks and bachelors (except for V0 and cascades)
   Configurable<LabeledArray<float>> ptCuts{"ptCuts", {cutsPt[0], 2, 10, labelsRowsCutsPt, labelsColumnsCutsPt}, "minimum and maximum pT for bachelor tracks (except for V0 and cascades)"};
-
+  Configurable<LabeledArray<float>> trackQaulityCuts{"trackQaulityCuts", {cutsTrackQuality[0], 2, 7, labelsColumnsPtThresholdsForFemto, labelsColumnsTrackQuality}, "Track quality cuts for proton and deuteron)"};
   // parameters for high-pT triggers
   Configurable<LabeledArray<float>> ptThresholds{"ptThresholds", {cutsHighPtThresholds[0], 1, 2, labelsEmpty, labelsColumnsHighPtThresholds}, "pT treshold for high pT charm hadron candidates for kHighPt triggers in GeV/c"};
 
@@ -136,8 +136,8 @@
 
   // CCDB configuration
   o2::ccdb::CcdbApi ccdbApi;
   Service<o2::ccdb::BasicCCDBManager> ccdb;
   Configurable<std::string> url{"ccdb-url", "http://alice-ccdb.cern.ch", "url of the ccdb repository"};
   int currentRun{0}; // needed to detect if the run changed and trigger update of calibrations etc.
 
   // TPC PID calibrations
@@ -167,7 +167,7 @@
   o2::vertexing::DCAFitterN<3> dfBtoDstar;    // fitter for Beauty Hadron to D* vertex (3-prong vertex fitter)
   o2::vertexing::DCAFitterN<2> dfStrangeness; // fitter for V0s and cascades (2-prong vertex fitter)
 
   HistogramRegistry registry{"registry"};
   std::shared_ptr<TH1> hProcessedEvents;
 
   // QA histos
@@ -215,6 +215,7 @@
     helper.setCutsBtoJPsi(cutsBtoHadrons.cutsBtoJPsiX);
     helper.setNsigmaProtonCutsForFemto(std::array{nSigmaPidCuts->get(0u, 3u), nSigmaPidCuts->get(1u, 3u), nSigmaPidCuts->get(2u, 3u), nSigmaPidCuts->get(3u, 3u)});
     helper.setNsigmaDeuteronCutsForFemto(std::array{nSigmaPidCuts->get(0u, 6u), nSigmaPidCuts->get(1u, 6u), nSigmaPidCuts->get(2u, 6u), nSigmaPidCuts->get(3u, 6u)});
+    helper.setDeuteronTrackSelectionForFemto(trackQaulityCuts->get(1u, 0u), trackQaulityCuts->get(1u, 1u), trackQaulityCuts->get(1u, 2u), trackQaulityCuts->get(1u, 3u), trackQaulityCuts->get(1u, 4u), trackQaulityCuts->get(1u, 5u), trackQaulityCuts->get(1u, 6u));
     helper.setNsigmaProtonCutsForCharmBaryons(nSigmaPidCuts->get(0u, 0u), nSigmaPidCuts->get(1u, 0u));
     helper.setNsigmaPionKaonCutsForDzero(nSigmaPidCuts->get(0u, 1u), nSigmaPidCuts->get(1u, 1u));
     helper.setNsigmaKaonCutsFor3Prongs(nSigmaPidCuts->get(0u, 2u), nSigmaPidCuts->get(1u, 2u));
@@ -242,7 +243,7 @@
 
     hProcessedEvents = registry.add<TH1>("fProcessedEvents", "HF - event filtered;;counts", HistType::kTH1D, {{kNtriggersHF + 2, -0.5, +kNtriggersHF + 1.5}});
     for (auto iBin = 0; iBin < kNtriggersHF + 2; ++iBin) {
       if (iBin < 2)
         hProcessedEvents->GetXaxis()->SetBinLabel(iBin + 1, eventTitles[iBin].data());
       else
         hProcessedEvents->GetXaxis()->SetBinLabel(iBin + 1, hfTriggerNames[iBin - 2].data());
@@ -389,7 +390,7 @@
     for (const auto& collision : collisions) {
 
       bool keepEvent[kNtriggersHF]{false};
       if (!collision.sel8() || std::fabs(collision.posZ()) > 11.f) { // safety margin for Zvtx
         tags(keepEvent[kHighPt2P], keepEvent[kHighPt3P], keepEvent[kBeauty3P], keepEvent[kBeauty4P], keepEvent[kFemto2P], keepEvent[kFemto3P], keepEvent[kDoubleCharm2P], keepEvent[kDoubleCharm3P], keepEvent[kDoubleCharmMix], keepEvent[kV0Charm2P], keepEvent[kV0Charm3P], keepEvent[kCharmBarToXiBach], keepEvent[kSigmaCPPK], keepEvent[kSigmaC0K0], keepEvent[kPhotonCharm2P], keepEvent[kPhotonCharm3P], keepEvent[kSingleCharm2P], keepEvent[kSingleCharm3P], keepEvent[kSingleNonPromptCharm2P], keepEvent[kSingleNonPromptCharm3P], keepEvent[kCharmBarToXi2Bach], keepEvent[kPrCharm2P], keepEvent[kBtoJPsiKa], keepEvent[kBtoJPsiKstar], keepEvent[kBtoJPsiPhi], keepEvent[kBtoJPsiPrKa], keepEvent[kBtoJPsiPi]);
         continue;
       }
@@ -468,7 +469,7 @@
         std::vector<float> scores{};
         if (preselD0) {
           scores.insert(scores.end(), cand2Prong.mlProbSkimD0ToKPi().begin(), cand2Prong.mlProbSkimD0ToKPi().end());
           if (scores.size() != 3) {
             scores.resize(3);
             scores[0] = 2.;
             scores[1] = -1.;
@@ -1187,7 +1188,7 @@
             continue;
           }
 
           if (scores[iCharmPart].size() != 3) {
             scores[iCharmPart].resize(3);
             scores[iCharmPart][0] = 2.;
             scores[iCharmPart][1] = -1.;
@@ -1414,7 +1415,7 @@
               // check the candidate SigmaC++ charge
               std::array<int, 4> chargesSc = {trackFirst.sign(), trackSecond.sign(), trackThird.sign(), trackSoftPi.sign()};
               int chargeSc = std::accumulate(chargesSc.begin(), chargesSc.end(), 0); // SIGNED electric charge of SigmaC candidate
               if (std::abs(chargeSc) != 2) {
                 continue;
               }
 
@@ -1736,7 +1737,7 @@
           std::array<float, 2> dcaInfo;
           o2::track::TrackParCov trackParCasc;
           o2::track::TrackParCov trackParCascTrack;
           if (requireStrangenessTrackingAny < 2) { // needed for at least one of the two
             trackParCasc = o2::track::TrackParCov(cascCand.vtx, cascCand.mom, cascCand.cov, cascCand.sign, true);
             trackParCasc.setPID(o2::track::PID::XiMinus);
             trackParCasc.setAbsCharge(1); // to be sure

EventFiltering/PWGHF/HFFilterHelpers.h

@@ -313,6 +313,11 @@ constexpr float cutsNsigma[4][8] = {
 static const std::vector<std::string> labelsColumnsNsigma = {"PrFromLc", "PiKaFromDZero", "KaFrom3Prong", "PrForFemto", "PiKaFromCharmBaryon", "SoftKaonFromXicResoToSigmaC", "DeForFemto", "KaPrFromBeautyToJPsi"};
 static const std::vector<std::string> labelsRowsNsigma = {"TPC", "TOF", "Comb", "ITS"};
 
+// track cut
+constexpr float cutsTrackQuality[2][7] = {{0., 0., 0., 999., 999., 0., 0.},
+                                          {90, 80, 0.83, 160., 1., 5., 0.}};
+static const std::vector<std::string> labelsColumnsTrackQuality = {"minTpcCluster", "minTpcRow", "minTpcCrossedOverFound", "maxTpcShared", "maxTpcFracShared", "minItsCluster", "minItsIbCluster"};
+
 // high pt
 constexpr float cutsHighPtThresholds[1][2] = {{8., 8.}}; // 2-prongs, 3-prongs
 static const std::vector<std::string> labelsColumnsHighPtThresholds = {"2Prongs", "3Prongs"};
@@ -489,6 +494,18 @@ class HfFilterHelper
 
   void setNsigmaProtonCutsForFemto(std::array<float, 4> nSigmaCuts) { mNSigmaPrCutsForFemto = nSigmaCuts; }
   void setNsigmaDeuteronCutsForFemto(std::array<float, 4> nSigmaCuts) { mNSigmaDeCutsForFemto = nSigmaCuts; }
+
+  void setDeuteronTrackSelectionForFemto(float minTpcCluster, float minTpcRow, float minTpcCrossedOverFound, float maxTpcShared, float maxTpcFracShared, float minItsCluster, float minItsIbCluster)
+  {
+    mMinTpcCluster = minTpcCluster;
+    mMinTpcRow = minTpcRow;
+    mMinTpcCrossedOverFound = minTpcCrossedOverFound;
+    mMaxTpcShared = maxTpcShared;
+    mMaxTpcFracShared = maxTpcFracShared;
+    mMinItsCluster = minItsCluster;
+    mMinItsIbCluster = minItsIbCluster;
+  }
+
   void setNsigmaProtonCutsForCharmBaryons(float nSigmaTpc, float nSigmaTof)
   {
     mNSigmaTpcPrCutForCharmBaryons = nSigmaTpc;
@@ -754,7 +771,13 @@ class HfFilterHelper
   std::array<float, 2> mCosPaMinXiBach{-2.f, -2.f};                               // minimum cosine of pointing angle for XiBachelor candidates
   std::array<o2::framework::LabeledArray<double>, kNBeautyParticles> mCutsBhad{}; // selections for B-hadron candidates (DeltaMass, CPA, DecayLength, ImpactParameterProduct)
   o2::framework::LabeledArray<double> mCutsBhadToJPsi{};                          // selections for B->JPsi candidates (PtMinMu, DeltaMass, CPA, DecayLength)
-
+  float mMinTpcCluster{90.};                                                      // Minimum number of TPC clusters required on a track
+  float mMinTpcRow{80.};                                                          // Minimum number of TPC rows (pad rows) traversed by the track
+  float mMinTpcCrossedOverFound{0.83};                                            // Minimum ratio of crossed TPC rows over findable clusters
+  float mMaxTpcShared{160.};                                                      // Maximum allowed number of shared TPC clusters between tracks
+  float mMaxTpcFracShared{1.};                                                    // Maximum allowed fraction of shared TPC clusters relative to total clusters
+  float mMinItsCluster{1.};                                                       // Minimum required number of ITS clusters
+  float mMinItsIbCluster{1.};                                                     // Minimum required number of ITS clusters for IB
   // PID recalibrations
   int mTpcPidCalibrationOption{0};                          // Option for TPC PID calibration (0 -> AO2D, 1 -> postcalibrations, 2 -> alternative bethe bloch parametrisation)
   std::array<TH3F*, 8> mHistMapPiPrKaDe{};                  // Map for TPC PID postcalibrations for pions, kaon, protons and deuterons
@@ -941,9 +964,9 @@ inline bool HfFilterHelper::isSelectedTrack4Femto(const T1& track, const T2& tra
   }
 
   float NSigma = std::sqrt(NSigmaTPC * NSigmaTPC + NSigmaTOF * NSigmaTOF);
-
+  float momentum = track.p();
   if (trackSpecies == kProtonForFemto) {
-    if (pt <= ptThresholdPidStrategy) {
+    if (momentum <= ptThresholdPidStrategy) {
       if (NSigma > nSigmaCuts[2]) {
         return false;
       }
@@ -955,9 +978,32 @@ inline bool HfFilterHelper::isSelectedTrack4Femto(const T1& track, const T2& tra
   }
   // For deuterons: Determine whether to apply TOF based on pt threshold
   if (trackSpecies == kDeuteronForFemto) {
+
+    if (track.tpcNClsFound() < mMinTpcCluster) {
+      return false;
+    }
+    if (track.tpcNClsCrossedRows() < mMinTpcRow) {
+      return false;
+    }
+    if (track.tpcCrossedRowsOverFindableCls() < mMinTpcCrossedOverFound) {
+      return false;
+    }
+    if (track.tpcNClsShared() > mMaxTpcShared) {
+      return false;
+    }
+    if (track.tpcFractionSharedCls() > mMaxTpcFracShared) {
+      return false;
+    }
+    if (track.itsNCls() < mMinItsCluster) {
+      return false;
+    }
+    if (track.itsNClsInnerBarrel() < mMinItsIbCluster) {
+      return false;
+    }
+
     // Apply different PID strategy in different pt range
     // one side selection only
-    if (pt <= ptThresholdPidStrategy) {
+    if (momentum <= ptThresholdPidStrategy) {
       if (std::fabs(NSigmaTPC) > nSigmaCuts[0] || NSigmaITS < -nSigmaCuts[3]) { // Use TPC and ITS below the threshold, NSigmaITS for deuteron with a lower limit
         return false;
       }
@@ -972,9 +1018,9 @@ inline bool HfFilterHelper::isSelectedTrack4Femto(const T1& track, const T2& tra
     hTPCPID->Fill(track.p(), NSigmaTPC);
     if ((!forceTof || track.hasTOF())) {
       if (trackSpecies == kProtonForFemto)
-        hTOFPID->Fill(track.p(), NSigmaTOF);
-      else if (trackSpecies == kDeuteronForFemto && pt > ptThresholdPidStrategy)
-        hTOFPID->Fill(track.p(), NSigmaTOF);
+        hTOFPID->Fill(momentum, NSigmaTOF);
+      else if (trackSpecies == kDeuteronForFemto && momentum > ptThresholdPidStrategy)
+        hTOFPID->Fill(momentum, NSigmaTOF);
     }
   }