ITS: debug dump

Signed-off-by: Felix Schlepper <felix.schlepper@cern.ch>

Author: f3sch

Detectors/ITSMFT/ITS/tracking/include/ITStracking/Definitions.h

@@ -18,6 +18,8 @@
 #include <type_traits>
 #include <cstdint>
 
+#define ITS_ALLOW_DEBUG_TREES
+
 namespace o2::its
 {
 
@@ -26,10 +28,14 @@ enum class TrackletMode {
   Layer1Layer2 = 2
 };
 
-template <bool IsConst, typename T>
-using maybe_const = typename std::conditional<IsConst, const T, T>::type;
+#ifdef ITS_ALLOW_DEBUG_TREES
+enum DebugTreeFlags : uint32_t {
+  DebugTrackletCandidates = 1u << 0,
+  DebugTracklets = 1u << 1,
+  DebugCells = 1u << 2,
+};
 
-// simple implemnetion of logging with exp. backoff
+// Simple implementation of logging with exponential backoff.
 struct LogLogThrottler {
   uint64_t evCount{0};
   uint64_t nextLog{1};
@@ -50,6 +56,7 @@ struct LogLogThrottler {
     return false;
   }
 };
+#endif
 
 struct TimingStats {
   std::uint64_t calls = 0;
@@ -63,6 +70,9 @@ struct TimingStats {
   double averageTimeMs() const { return calls ? totalTimeMs / static_cast<double>(calls) : 0.; }
 };
 
+template <bool IsConst, typename T>
+using maybe_const = typename std::conditional<IsConst, const T, T>::type;
+
 } // namespace o2::its
 
 #endif

Detectors/ITSMFT/ITS/tracking/include/ITStracking/TimeFrame.h

@@ -14,6 +14,9 @@
 #define TRACKINGITSU_INCLUDE_TIMEFRAME_H_
 
 #include <array>
+#include <cstdint>
+#include <memory>
+#include <string>
 #include <vector>
 #include <utility>
 #include <algorithm>
@@ -26,6 +29,7 @@
 #include "ITStracking/Cell.h"
 #include "ITStracking/Cluster.h"
 #include "ITStracking/Configuration.h"
+#include "ITStracking/Definitions.h"
 #include "ITStracking/ClusterLines.h"
 #include "ITStracking/Tracklet.h"
 #include "ITStracking/IndexTableUtils.h"
@@ -35,11 +39,17 @@
 #include "ITStracking/TrackingTopology.h"
 #include "SimulationDataFormat/MCCompLabel.h"
 #include "SimulationDataFormat/MCTruthContainer.h"
-
 #include "DetectorsBase/Propagator.h"
 
 namespace o2
 {
+#ifdef ITS_ALLOW_DEBUG_TREES
+namespace utils
+{
+class TreeStreamRedirector;
+}
+#endif
+
 namespace gpu
 {
 class GPUChainITS;
@@ -292,6 +302,20 @@ struct TimeFrame {
   virtual bool isGPU() const noexcept { return false; }
   virtual const char* getName() const noexcept { return "CPU"; }
 
+#ifdef ITS_ALLOW_DEBUG_TREES
+  std::string mDBGTreeFileName{"dbg_its.root"};
+  void setDebugTreeFileName(std::string n) { mDBGTreeFileName = std::move(n); }
+  static void destroyDebugTree(o2::utils::TreeStreamRedirector*);
+  using DebugTreePtr = std::unique_ptr<o2::utils::TreeStreamRedirector, void (*)(o2::utils::TreeStreamRedirector*)>;
+  DebugTreePtr mDBGOut{nullptr, &TimeFrame::destroyDebugTree};
+  void closeDebugTree();
+  auto getDBGOut() { return mDBGOut.get(); }
+  uint32_t mDBGFlags{0};
+  bool isDebugFlag(DebugTreeFlags flag) const noexcept { return (mDBGFlags & static_cast<uint32_t>(flag)) != 0; }
+  void setDebugFlags(uint32_t flags) noexcept { mDBGFlags = flags; }
+  auto getDebugFlags() const noexcept { return mDBGFlags; }
+#endif
+
  protected:
   void prepareClusters(const TrackingParameters& trkParam, const int maxLayers = NLayers);
   float mBz = 5.;

Detectors/ITSMFT/ITS/tracking/include/ITStracking/TrackerTraits.h

@@ -81,6 +81,13 @@ class TrackerTraits
   virtual int getTFNumberOfCells() const { return mTimeFrame->getNumberOfCells(); }
 
  private:
+#ifdef ITS_ALLOW_DEBUG_TREES
+  void dumpTracklet(int iteration, int transitionId, int fromLayer, int toLayer, const Tracklet& tracklet);
+  void dumpTrackletCandidate(int iteration, int transitionId, int fromLayer, int toLayer, int pivotROF, int targetROF,
+                             const Cluster& currentCluster, const Cluster& nextCluster, const Vertex& pv,
+                             float sigmaZ, float phiCut, float deltaZ, bool accepted);
+#endif
+
   std::shared_ptr<BoundedMemoryResource> mMemoryPool;
   std::shared_ptr<tbb::task_arena> mTaskArena;
 

Detectors/ITSMFT/ITS/tracking/src/IOUtils.cxx

@@ -23,12 +23,6 @@
 #include "ITStracking/TrackingConfigParam.h"
 #include "ITSMFTReconstruction/ChipMappingITS.h"
 
-namespace
-{
-constexpr int PrimaryVertexLayerId{-1};
-constexpr int EventLabelsSeparator{-1};
-} // namespace
-
 using namespace o2::its;
 
 /// convert compact clusters to 3D spacepoints

Detectors/ITSMFT/ITS/tracking/src/TimeFrame.cxx

@@ -15,6 +15,10 @@
 
 #include <numeric>
 
+#include "ITStracking/Definitions.h"
+#ifdef ITS_ALLOW_DEBUG_TREES
+#include "CommonUtils/TreeStreamRedirector.h"
+#endif
 #include "Framework/Logger.h"
 #include "ITStracking/TimeFrame.h"
 #include "ITStracking/MathUtils.h"
@@ -37,12 +41,28 @@ struct ClusterHelper {
 
 namespace o2::its
 {
-
 constexpr float DefClusErrorRow = o2::itsmft::SegmentationAlpide::PitchRow * 0.5;
 constexpr float DefClusErrorCol = o2::itsmft::SegmentationAlpide::PitchCol * 0.5;
 constexpr float DefClusError2Row = DefClusErrorRow * DefClusErrorRow;
 constexpr float DefClusError2Col = DefClusErrorCol * DefClusErrorCol;
 
+#ifdef ITS_ALLOW_DEBUG_TREES
+template <int NLayers>
+void TimeFrame<NLayers>::destroyDebugTree(o2::utils::TreeStreamRedirector* stream)
+{
+  delete stream;
+}
+
+template <int NLayers>
+void TimeFrame<NLayers>::closeDebugTree()
+{
+  if (mDBGOut) {
+    mDBGOut->Close();
+    mDBGOut.reset();
+  }
+}
+#endif
+
 template <int NLayers>
 void TimeFrame<NLayers>::addPrimaryVertex(const Vertex& vert)
 {
@@ -265,6 +285,12 @@ void TimeFrame<NLayers>::initTrackerTopologies(gsl::span<const TrackingParameter
 template <int NLayers>
 void TimeFrame<NLayers>::initialise(const TrackingParameters& trkParam, const int maxLayers, const int iteration)
 {
+#ifdef ITS_ALLOW_DEBUG_TREES
+  if (mDBGFlags != 0 && !mDBGOut) {
+    mDBGOut.reset(new o2::utils::TreeStreamRedirector(mDBGTreeFileName.c_str(), "recreate"));
+  }
+#endif
+
   mTrackingTopologyView = iteration != constants::UnusedIndex ? mTrackerTopologies[iteration].getView() : (maxLayers == 3 ? mVertexingTopology.getView() : mDefaultTrackingTopology.getView());
 
   if (trkParam.PassFlags[IterationStep::FirstPass]) {

Detectors/ITSMFT/ITS/tracking/src/TrackerTraits.cxx

@@ -21,6 +21,10 @@
 #include <oneapi/tbb/blocked_range.h>
 #include <oneapi/tbb/enumerable_thread_specific.h>
 
+#include "ITStracking/Definitions.h"
+#ifdef ITS_ALLOW_DEBUG_TREES
+#include "CommonUtils/TreeStreamRedirector.h"
+#endif
 #include "DetectorsBase/Propagator.h"
 #include "GPUCommonMath.h"
 #include "ITStracking/BoundedAllocator.h"
@@ -149,9 +153,19 @@ void TrackerTraits<NLayers>::computeLayerTracklets(const int iteration, int iVer
                   continue;
                 }
                 const float deltaZ = o2::gpu::CAMath::Abs((tanLambda * (nextCluster.radius - currentCluster.radius)) + currentCluster.zCoordinate - nextCluster.zCoordinate);
+                const bool accepted = deltaZ / sigmaZ < mTrkParams[iteration].NSigmaCut &&
+                                      math_utils::isPhiDifferenceBelow(currentCluster.phi, nextCluster.phi, phiCut);
+
+#ifdef ITS_ALLOW_DEBUG_TREES
+                if constexpr (decltype(Tag)::value == PassMode::OnePass::value) {
+                  if (mTimeFrame->isDebugFlag(DebugTrackletCandidates)) {
+                    dumpTrackletCandidate(iteration, transitionId, transition.fromLayer, transition.toLayer, pivotROF, targetROF,
+                                          currentCluster, nextCluster, pv, sigmaZ, phiCut, deltaZ, accepted);
+                  }
+                }
+#endif
 
-                if (deltaZ / sigmaZ < mTrkParams[iteration].NSigmaCut &&
-                    math_utils::isPhiDifferenceBelow(currentCluster.phi, nextCluster.phi, phiCut)) {
+                if (accepted) {
                   const float phi{o2::gpu::CAMath::ATan2(currentCluster.yCoordinate - nextCluster.yCoordinate, currentCluster.xCoordinate - nextCluster.xCoordinate)};
                   const float tanL = (currentCluster.zCoordinate - nextCluster.zCoordinate) / (currentCluster.radius - nextCluster.radius);
                   if constexpr (decltype(Tag)::value == PassMode::OnePass::value) {
@@ -212,6 +226,14 @@ void TrackerTraits<NLayers>::computeLayerTracklets(const int iteration, int iVer
       std::sort(trkl.begin(), trkl.end());
       trkl.erase(std::unique(trkl.begin(), trkl.end()), trkl.end());
       trkl.shrink_to_fit();
+#ifdef ITS_ALLOW_DEBUG_TREES
+      if (mTimeFrame->isDebugFlag(DebugTracklets)) {
+        const auto& transition = topology.getTransition(transitionId);
+        for (const auto& tkl : trkl) {
+          dumpTracklet(iteration, transitionId, transition.fromLayer, transition.toLayer, tkl);
+        }
+      }
+#endif
       auto& lut{mTimeFrame->getTrackletsLookupTable()[transitionId]};
       if (!trkl.empty()) {
         for (const auto& tkl : trkl) {
@@ -903,18 +925,129 @@ void TrackerTraits<NLayers>::setBz(float bz)
 template <int NLayers>
 void TrackerTraits<NLayers>::setNThreads(int n, std::shared_ptr<tbb::task_arena>& arena)
 {
-#if defined(OPTIMISATION_OUTPUT)
-  mTaskArena = std::make_shared<tbb::task_arena>(1);
-#else
+#ifdef ITS_ALLOW_DEBUG_TREES
+  if (mTimeFrame->getDebugFlags() != 0) {
+    mTaskArena = std::make_shared<tbb::task_arena>(1);
+    LOGP(info, "Forcing ITS tracker to 1 thread for debug tree output.");
+    return;
+  }
+#endif
   if (arena == nullptr) {
     mTaskArena = std::make_shared<tbb::task_arena>(std::abs(n));
     LOGP(info, "Setting tracker with {} threads.", n);
   } else {
     mTaskArena = arena;
   }
-#endif
 }
 
+#ifdef ITS_ALLOW_DEBUG_TREES
+template <int NLayers>
+void TrackerTraits<NLayers>::dumpTracklet(int iteration, int transitionId, int fromLayer, int toLayer, const Tracklet& tracklet)
+{
+  auto* dbgOut = mTimeFrame->getDBGOut();
+  if (dbgOut == nullptr || !mTimeFrame->isDebugFlag(DebugTracklets)) {
+    return;
+  }
+
+  static LogLogThrottler logger;
+  if (logger.needToLog(iteration, transitionId)) {
+    LOGP(info, "debugTree: LayerTracklets:{}:{} dumped entries {}", iteration, transitionId, logger.evCount);
+  }
+
+  MCCompLabel label;
+  if (mTimeFrame->hasMCinformation()) {
+    const int currentId = mTimeFrame->getClusters()[fromLayer][tracklet.firstClusterIndex].clusterId;
+    const int nextId = mTimeFrame->getClusters()[toLayer][tracklet.secondClusterIndex].clusterId;
+    for (const auto& lab1 : mTimeFrame->getClusterLabels(fromLayer, currentId)) {
+      for (const auto& lab2 : mTimeFrame->getClusterLabels(toLayer, nextId)) {
+        if (lab1 == lab2 && lab1.isValid()) {
+          label = lab1;
+          break;
+        }
+      }
+      if (label.isValid()) {
+        break;
+      }
+    }
+  }
+
+  const float tglTracklet = tracklet.tanLambda;
+  const float etaTracklet = std::asinh(tglTracklet);
+  const bool accepted = true;
+
+  (*dbgOut) << "tracklets"
+            << "iter=" << iteration
+            << "transition=" << transitionId
+            << "fromLayer=" << fromLayer
+            << "toLayer=" << toLayer
+            << "lbl=" << label
+            << "tglTracklet=" << tglTracklet
+            << "etaTracklet=" << etaTracklet
+            << "acc=" << accepted
+            << "\n";
+}
+
+template <int NLayers>
+void TrackerTraits<NLayers>::dumpTrackletCandidate(int iteration, int transitionId, int fromLayer, int toLayer, int pivotROF, int targetROF,
+                                                   const Cluster& currentCluster, const Cluster& nextCluster, const Vertex& pv,
+                                                   float sigmaZ, float phiCut, float deltaZ, bool accepted)
+{
+  auto* dbgOut = mTimeFrame->getDBGOut();
+  if (dbgOut == nullptr || !mTimeFrame->isDebugFlag(DebugTrackletCandidates)) {
+    return;
+  }
+
+  static LogLogThrottler logger;
+  if (logger.needToLog(iteration, transitionId)) {
+    LOGP(info, "debugTree: LayerTrackletCandidates:{}:{} dumped entries {}", iteration, transitionId, logger.evCount);
+  }
+
+  MCCompLabel label;
+  if (mTimeFrame->hasMCinformation()) {
+    const int currentId = currentCluster.clusterId;
+    const int nextId = nextCluster.clusterId;
+    for (const auto& lab1 : mTimeFrame->getClusterLabels(fromLayer, currentId)) {
+      for (const auto& lab2 : mTimeFrame->getClusterLabels(toLayer, nextId)) {
+        if (lab1 == lab2 && lab1.isValid()) {
+          label = lab1;
+          break;
+        }
+      }
+      if (label.isValid()) {
+        break;
+      }
+    }
+  }
+
+  const float tglTracklet = (currentCluster.zCoordinate - nextCluster.zCoordinate) / (currentCluster.radius - nextCluster.radius);
+  const float etaTracklet = std::asinh(tglTracklet);
+  const float deltaPhiRaw = o2::gpu::CAMath::Abs(currentCluster.phi - nextCluster.phi);
+  const float deltaPhi = o2::gpu::CAMath::Min(deltaPhiRaw, o2::gpu::CAMath::Abs(deltaPhiRaw - o2::constants::math::TwoPI));
+  const float zTrackletAtPvR = currentCluster.zCoordinate + tglTracklet * (pv.getR() - currentCluster.radius);
+  const float deltaZTrackletPv = zTrackletAtPvR - pv.getZ();
+
+  (*dbgOut) << "trackletCandidates"
+            << "iter=" << iteration
+            << "transition=" << transitionId
+            << "fromLayer=" << fromLayer
+            << "toLayer=" << toLayer
+            << "pivotROF=" << pivotROF
+            << "targetROF=" << targetROF
+            << "lbl=" << label
+            << "pvZ=" << pv.getZ()
+            << "sigZ=" << sigmaZ
+            << "delZ=" << deltaZ
+            << "sigPhi=" << phiCut
+            << "delPhi=" << deltaPhi
+            << "tglTracklet=" << tglTracklet
+            << "etaTracklet=" << etaTracklet
+            << "zTrackletAtPvR=" << zTrackletAtPvR
+            << "deltaZTrackletPv=" << deltaZTrackletPv
+            << "acc=" << accepted
+            << "\n";
+}
+#endif
+
 template class TrackerTraits<7>;
 template void TrackerTraits<7>::processNeighbours<CellSeed>(int, int, int, const bounded_vector<CellSeed>&, const bounded_vector<int>&, const bounded_vector<int>&, bounded_vector<TrackSeed<7>>&, bounded_vector<int>&, bounded_vector<int>&);
 template void TrackerTraits<7>::processNeighbours<TrackSeed<7>>(int, int, int, const bounded_vector<TrackSeed<7>>&, const bounded_vector<int>&, const bounded_vector<int>&, bounded_vector<TrackSeed<7>>&, bounded_vector<int>&, bounded_vector<int>&);

Detectors/ITSMFT/ITS/tracking/src/TrackingInterface.cxx

@@ -467,6 +467,9 @@ void ITSTrackingInterface::finaliseCCDB(ConcreteDataMatcher& matcher, void* obj)
 
 void ITSTrackingInterface::printSummary() const
 {
+#ifdef ITS_ALLOW_DEBUG_TREES
+  mTimeFrame->closeDebugTree();
+#endif
   mVertexer->printSummary();
   mTracker->printSummary();
 }