SV: Relax TPC cuts and opt. debug output

Signed-off-by: Felix Schlepper <[email protected]>

Author: f3sch

Detectors/GlobalTrackingWorkflow/src/SecondaryVertexingSpec.cxx

@@ -96,6 +96,7 @@ void SecondaryVertexingSpec::init(InitContext& ic)
   mVertexer.setEnable3BodyDecays(mEnable3BodyVertices);
   mVertexer.setNThreads(ic.options().get<int>("threads"));
   mVertexer.setUseMC(mUseMC);
+  mVertexer.setDebug(ic.options().get<std::string>("debug"));
   if (mEnableStrangenessTracking) {
     mStrTracker.setCorrType(o2::base::PropagatorImpl<float>::MatCorrType::USEMatCorrLUT);
     mStrTracker.setConfigParams(&o2::strangeness_tracking::StrangenessTrackingParamConfig::Instance());
@@ -131,6 +132,7 @@ void SecondaryVertexingSpec::run(ProcessingContext& pc)
 
 void SecondaryVertexingSpec::endOfStream(EndOfStreamContext& ec)
 {
+  mVertexer.finalize();
   LOGF(info, "Secondary vertexing total timing: Cpu: %.3e Real: %.3e s in %d slots, nThreads = %d",
        mTimer.CpuTime(), mTimer.RealTime(), mTimer.Counter() - 1, mVertexer.getNThreads());
 }
@@ -243,6 +245,7 @@ DataProcessorSpec getSecondaryVertexingSpec(GTrackID::mask_t src, bool enableCas
 {
   std::vector<OutputSpec> outputs;
   Options opts{
+    {"debug", VariantType::String, "", {"Request debug output"}},
     {"material-lut-path", VariantType::String, "", {"Path of the material LUT file"}},
     {"threads", VariantType::Int, 1, {"Number of threads"}}};
   auto dataRequest = std::make_shared<DataRequest>();

Detectors/Vertexing/CMakeLists.txt

@@ -9,7 +9,8 @@
 # granted to it by virtue of its status as an Intergovernmental Organization
 # or submit itself to any jurisdiction.
 
-# add_compile_options(-O0 -g -fPIC)
+# add_compile_options(-O0 -g -fPIC -fno-omit-frame-pointer)
+add_compile_definitions(SVERTEXER_DEBUG=1)
 
 add_subdirectory(StrangenessTracking)
 

Detectors/Vertexing/include/DetectorsVertexing/SVertexer.h

@@ -32,11 +32,16 @@
 #include "DetectorsVertexing/SVertexHypothesis.h"
 #include "StrangenessTracking/StrangenessTracker.h"
 #include "DataFormatsTPC/TrackTPC.h"
-#include <numeric>
-#include <algorithm>
 #include "GPUO2InterfaceRefit.h"
 #include "TPCFastTransform.h"
 #include "DataFormatsTPC/PIDResponse.h"
+#include "CommonUtils/TreeStreamRedirector.h"
+#include "CommonUtils/EnumFlags.h"
+
+#include <numeric>
+#include <algorithm>
+#include <span>
+#include <memory>
 
 namespace o2
 {
@@ -104,6 +109,7 @@ class SVertexer
     bool hasTPC = false;
     int8_t nITSclu = -1;
     bool compatibleProton = false; // dE/dx compatibility with proton hypothesis (FIXME: use better, uint8_t compat mask?)
+    float trk2beam = 0.;
     bool hasITS() const
     {
       return nITSclu > 0;
@@ -117,6 +123,7 @@ class SVertexer
   void setEnableCascades(bool v) { mEnableCascades = v; }
   void setEnable3BodyDecays(bool v) { mEnable3BodyDecays = v; }
   void init();
+  void finalize();
   void process(const o2::globaltracking::RecoContainer& recoTracks, o2::framework::ProcessingContext& pc);
   void produceOutput(o2::framework::ProcessingContext& pc);
   int getNV0s() const { return mNV0s; }
@@ -150,17 +157,47 @@ class SVertexer
   std::array<size_t, 3> getNFitterCalls() const;
   void setSources(GIndex::mask_t src) { mSrc = src; }
 
+  void setDebug(const std::string& s) { mDebugFlags.set(s); }
+
  private:
   template <class TVI, class TCI, class T3I, class TR>
-  void extractPVReferences(const TVI& v0s, TR& vtx2V0Refs, const TCI& cascades, TR& vtx2CascRefs, const T3I& vtxs3, TR& vtx2body3Refs);
-  bool checkV0(const TrackCand& seed0, const TrackCand& seed1, int iP, int iN, int ithread);
+  void
+    extractPVReferences(const TVI& v0s, TR& vtx2V0Refs, const TCI& cascades, TR& vtx2CascRefs, const T3I& vtxs3, TR& vtx2body3Refs);
+
+  enum class CheckV0Status : uint8_t {
+    None = 0,
+    OK,
+    RejTPCTgl,
+    RejTPC2Beam,
+    RejTPCDCA2,
+    RejFitter,
+    RejMinR2ToMeanVertex,
+    RejCausality,
+    RejProp,
+    RejPt2,
+    RejTgl,
+    RejHypo,
+    RejCascade,
+    RejDCA2,
+    RejNoCand,
+    RejNoCandNo3Body,
+  };
+  CheckV0Status checkV0(const TrackCand& seed0, const TrackCand& seed1, int iP, int iN, int ithread);
   int checkCascades(const V0Index& v0Idx, const V0& v0, float rv0, std::array<float, 3> pV0, float p2V0, int avoidTrackID, int posneg, VBracket v0vlist, int ithread);
   int check3bodyDecays(const V0Index& v0Idx, const V0& v0, float rv0, std::array<float, 3> pV0, float p2V0, int avoidTrackID, int posneg, VBracket v0vlist, int ithread);
   void setupThreads();
   void buildT2V(const o2::globaltracking::RecoContainer& recoTracks);
   void updateTimeDependentParams();
   bool acceptTrack(const GIndex gid, const o2::track::TrackParCov& trc) const;
-  bool processTPCTrack(const o2::tpc::TrackTPC& trTPC, GIndex gid, int vtxid);
+
+  enum class ProcessTPCTrackStatus : uint8_t {
+    Constrained,
+    OK,
+    RejMaxX,
+    RejCorr,
+    Rej2Beam,
+  };
+  ProcessTPCTrackStatus processTPCTrack(const o2::tpc::TrackTPC& trTPC, GIndex gid, int vtxid);
   float correctTPCTrack(TrackCand& trc, const o2::tpc::TrackTPC& tTPC, float tmus, float tmusErr) const;
 
   uint64_t getPairIdx(GIndex id1, GIndex id2) const
@@ -224,6 +261,57 @@ class SVertexer
   bool mEnableCascades = true;
   bool mEnable3BodyDecays = false;
   bool mUseMC = false;
+
+  /// Debug information
+  enum class SVDebug : uint8_t {
+    TPCCandPool, // check which TPC tracks do not go into the pool
+    Unassigned,  // check unassigned tracks
+    CheckV0,     // checkV0 status
+  };
+  utils::EnumFlags<SVDebug> mDebugFlags;
+  std::unique_ptr<utils::TreeStreamRedirector> mDebugStream;
+
+#ifdef SVERTEXER_DEBUG
+  template <bool log = true, typename... Args>
+    requires(sizeof...(Args) % 2 == 0 && sizeof...(Args) > 0)
+  void debugTree(SVDebug f, const char* tname, Args&&... args)
+  {
+    static struct LogLogThrottler {
+      size_t evCount{0};
+      size_t nextLog{1};
+      bool needToLog()
+      {
+        if (++evCount > nextLog) {
+          nextLog *= 2;
+          return true;
+        }
+        return false;
+      }
+    } logger;
+    if (mDebugFlags[f]) {
+      LOG_IF(info, logger.needToLog()) << "debugTree: " << tname << " dumped entries " << logger.evCount;
+      auto& stream = (*mDebugStream) << tname;
+      [&]<std::size_t... I>(std::index_sequence<I...>) {
+        const char* names[] = {std::get<2 * I>(std::forward_as_tuple(args...))...};
+        auto values = std::forward_as_tuple(std::get<2 * I + 1>(std::forward_as_tuple(args...))...);
+        ((stream << names[I] << [](auto&& val) -> decltype(auto) {
+           if constexpr (std::is_enum_v<std::decay_t<decltype(val)>>) {
+             return static_cast<std::underlying_type_t<std::decay_t<decltype(val)>>>(val);
+           } else {
+             return val;
+           }
+         }(std::get<I>(values))),
+         ...);
+      }(std::make_index_sequence<sizeof...(Args) / 2>{});
+      stream << "\n";
+    }
+  }
+#else // turn this into a no-op
+  template <typename... Args>
+  void debugTree(SVDebug f, const char* tname, Args&&... args)
+  {
+  }
+#endif
 };
 
 } // namespace vertexing

Detectors/Vertexing/include/DetectorsVertexing/SVertexerParams.h

@@ -93,20 +93,19 @@ struct SVertexerParams : public o2::conf::ConfigurableParamHelper<SVertexerParam
 
   // Cuts for TPC only tracks
   bool mExcludeTPCtracks = false;            // don't loop over TPC tracks if true (if loaded, dEdx info is used instead)
+  bool mTPCUseCollinearFit = true;          // use for TPC-only prongs collinear fit
   float mTPCTrackMaxX = -1.;                 // don't use TPC standalone tracks with X exceeding this;
   float mTPCTrack2Beam = 21.f;               // straight line for TPC track back to beamline
   bool mTPCTrackPhotonTune = true;           // use TPC-only photon tuning
-  int mTPCTrackMinNClusters = -1;            // minimum number of clusters
-  float mTPCTrackXY2Radius = 90.f;           // check for realistic conversion point, e.g., in material
-  float mTPCTrackD2R = 4.f;                  // check for maximal distance of pairs and their radii
-  float mTPCTrackDR = 90.f;                  // check if preliminary conversion point can lie somewhere reasonable
   float minTPCdEdx = 250;                    // starting from this dEdx value, tracks with p > minMomTPCdEdx are always accepted
   float minMomTPCdEdx = 0.8;                 // minimum p for tracks with dEdx > mMinTPCdEdx to be accepted
-  float maxV0TglAbsDiff = 0.3;               ///< max absolute difference in Tgl for V0 for photons only
+  float maxV0TglAbsDiff = 0.2;               ///< max absolute difference in Tgl for V0 for photons only
+  float maxV0BeamAbsDiff = 5.;               ///< max absolute difference in beamline extrapolation for V0 for photons only
   float mTPCTrackMaxChi2 = 4.;               ///< max DCA from prongs to vertex for photon TPC-only track only
   float mTPCTrackMaxDZIni = 8.;              ///< don't consider as a seed (circles intersection) if Z distance exceeds this, for photon TPC-only track only
   float mTPCTrackMaxDXYIni = 8.;             ///< don't consider as a seed (circles intersection) if XY distance exceeds this, for photon TPC-only track only
   float mTPCTrackMaxDCAXY2ToMeanVertex = 2.; ///< max DCA^2 of V0 from beam line (mean vertex) for prompt V0 candidates, for photon TPC-only track only
+  float maxTPCV0ToProngsRDiff = 90.;         ///< V0 radius cannot be lower than this ammount wrt minR of TPC-only contributors
 
   uint8_t mITSSAminNclu = 6;             // global requirement of at least this many ITS clusters if no TPC info present (N.B.: affects all secondary vertexing)
   uint8_t mITSSAminNcluCascades = 6;     // require at least this many ITS clusters if no TPC info present for cascade finding.
@@ -142,6 +141,9 @@ struct SVertexerParams : public o2::conf::ConfigurableParamHelper<SVertexerParam
   float pidCutsHe4L3body[SVertex3Hypothesis::NPIDParams] = {0.0025, 14, 0.07, 0.5}; // He4L -> He3 p pi-
   float pidCutsHe5L3body[SVertex3Hypothesis::NPIDParams] = {0.0025, 14, 0.07, 0.5}; // He5L -> He4 p pi-
 
+  // convienince getters
+  float getMaxV0ToProngsRDiff(bool isTPConly) const { return (isTPConly) ? maxTPCV0ToProngsRDiff : maxV0ToProngsRDiff; }
+
   O2ParamDef(SVertexerParams, "svertexer");
 };
 } // namespace vertexing