[WIP] Allow ACTS tracker in O2

Detectors/Upgrades/ALICE3/GlobalReconstruction/reconstruction/CMakeLists.txt

@@ -10,12 +10,18 @@
 # or submit itself to any jurisdiction.
 
 if(Acts_FOUND)
-  set(actsTarget Acts::Core)
+  set(actsSources src/TrackerACTS.cxx)
+  set(actsTarget
+        Acts::Core
+        O2::ACTSInterface
+        O2::ALICE3ACTS
+        ROOT::Tree
+        ROOT::RIO)
 endif()
 
 o2_add_library(ALICE3GlobalReconstruction
                TARGETVARNAME targetName
-               SOURCES src/TimeFrame.cxx
+               SOURCES src/TimeFrame.cxx ${actsSources}
                PUBLIC_LINK_LIBRARIES
                        O2::ITStracking
                        O2::GPUCommon

Detectors/Upgrades/ALICE3/GlobalReconstruction/reconstruction/include/ALICE3GlobalReconstruction/TimeFrameMixin.h

@@ -71,11 +71,12 @@ class TimeFrameMixin : public Base
   void deriveAndInitTiming(const std::array<gsl::span<const o2::trkft3::ROFRecord>, nLayers>& layerROFs);
 
   const o2::InteractionRecord& getTFAnchorIR() const noexcept { return mTFAnchorIR; }
+  const auto& getClusterSensorIDs() const noexcept { return mClusterSensorIDs; }
 
  protected:
   void initTimingTables(const std::array<o2::its::LayerTiming, nLayers>& timings);
   void updateHostROFVertexLookupTable();
-
+  std::array<std::vector<int>, nLayers> mClusterSensorIDs;
   bool mTimingTablesInitialised{false};
   o2::InteractionRecord mTFAnchorIR{0, 0};
 };
@@ -205,6 +206,7 @@ int TimeFrameMixin<nLayers, Base>::loadROFsFromHitTree(TTree* hitsTree, Geometry
     this->mUnsortedClusters[iLayer].clear();
     this->mTrackingFrameInfo[iLayer].clear();
     this->mClusterExternalIndices[iLayer].clear();
+    mClusterSensorIDs[iLayer].clear();
     this->mClusterSize[iLayer].clear();
   }
 
@@ -228,6 +230,7 @@ int TimeFrameMixin<nLayers, Base>::loadROFsFromHitTree(TTree* hitsTree, Geometry
     this->mUnsortedClusters[iLayer].reserve(clusterCountPerLayer[iLayer]);
     this->mTrackingFrameInfo[iLayer].reserve(clusterCountPerLayer[iLayer]);
     this->mClusterExternalIndices[iLayer].reserve(clusterCountPerLayer[iLayer]);
+    mClusterSensorIDs[iLayer].reserve(clusterCountPerLayer[iLayer]);
     this->mClusterSize[iLayer].reserve(clusterCountPerLayer[iLayer]);
   }
 
@@ -294,6 +297,7 @@ int TimeFrameMixin<nLayers, Base>::loadROFsFromHitTree(TTree* hitsTree, Geometry
                                         std::array<float, 2>{trkXYZ.y(), trkXYZ.z()},
                                         std::array<float, 3>{resolution[layer] * resolution[layer], 0., resolution[layer] * resolution[layer]});
       this->addClusterToLayer(layer, gloXYZ.x(), gloXYZ.y(), gloXYZ.z(), this->mUnsortedClusters[layer].size());
+      mClusterSensorIDs[layer].push_back(hit.GetDetectorID());
       const int layerHitCounter = hitCounterPerLayer[layer]++;
       this->addClusterExternalIndexToLayer(layer, layerHitCounter);
       this->mClusterSize[layer].push_back(1);
@@ -344,10 +348,12 @@ int TimeFrameMixin<nLayers, Base>::loadROFrameData(const std::array<gsl::span<co
     this->mUnsortedClusters[iLayer].clear();
     this->mTrackingFrameInfo[iLayer].clear();
     this->mClusterExternalIndices[iLayer].clear();
+    mClusterSensorIDs[iLayer].clear();
     this->mClusterSize[iLayer].clear();
     this->mUnsortedClusters[iLayer].reserve(layerClusters[iLayer].size());
     this->mTrackingFrameInfo[iLayer].reserve(layerClusters[iLayer].size());
     this->mClusterExternalIndices[iLayer].reserve(layerClusters[iLayer].size());
+    mClusterSensorIDs[iLayer].reserve(layerClusters[iLayer].size());
     this->mClusterSize[iLayer].reserve(layerClusters[iLayer].size());
   }
 
@@ -436,6 +442,7 @@ int TimeFrameMixin<nLayers, Base>::loadROFrameData(const std::array<gsl::span<co
                                           std::array<float, 2>{trkXYZ.y(), trkXYZ.z()},
                                           std::array<float, 3>{sigmaY2, 0.f, sigmaZ2});
         this->addClusterToLayer(iLayer, gloXYZ.x(), gloXYZ.y(), gloXYZ.z(), this->mUnsortedClusters[iLayer].size());
+        mClusterSensorIDs[iLayer].push_back(c.chipID);
         this->addClusterExternalIndexToLayer(iLayer, clusterId);
         this->mClusterSize[iLayer].push_back(std::clamp(static_cast<unsigned int>(c.size), 0u, 255u));
       }

Detectors/Upgrades/ALICE3/GlobalReconstruction/reconstruction/include/ALICE3GlobalReconstruction/TrackerACTS.h

@@ -19,119 +19,201 @@
 #ifndef ALICE3_GLOBALRECONSTRUCTION_INCLUDE_TRACKERACTS_H_
 #define ALICE3_GLOBALRECONSTRUCTION_INCLUDE_TRACKERACTS_H_
 
-#include "Acts/Definitions/Units.hpp"
-#include "Framework/Logger.h"
-
-#include "ITStracking/TimeFrame.h"
-#include "TH2F.h"
-
+#include <Acts/Definitions/Units.hpp>
+#include <chrono>
+#include <cmath>
+#include <iomanip>
+#include <memory>
+#include <string>
+#include <string_view>
+#include <Framework/Logger.h>
+#include <array>
+#include <vector>
+#include <unordered_set>
+#include "DataFormatsITS/TrackITS.h"
+#include <ITStracking/TimeFrame.h>
+#include <TH2F.h>
+
+class TFile;
+class TTree;
+namespace Acts
+{
+class Surface;
+class BoundTrackParameters;
+} // namespace Acts
 namespace o2::trk
 {
 
-/// Configuration for the ACTS-based tracker
+/// Configuration for ACTS seeding, track finding, and optional diagnostic output.
+///
+/// Length and momentum-valued cuts use ACTS units (mm and GeV, respectively).
 struct TrackerACTSConfig {
-  // Seeding parameters
-  float minPt = 0.4f * Acts::UnitConstants::GeV;             ///< Minimum pT for seeds
-  float maxImpactParameter = 10.f * Acts::UnitConstants::mm; ///< Maximum impact parameter
-  float cotThetaMax = std::sinh(4.0f);                       ///< Maximum cot(theta), corresponds to eta ~4
-
-  // Delta R cuts for doublet/triplet formation
+  /// Minimum transverse momentum accepted by the seeding algorithms.
+  float minPt = 0.4f * Acts::UnitConstants::GeV; ///< Minimum pT for seeds
+  /// Maximum transverse distance of closest approach used in seeding.
+  float maxImpactParameter = 10.f * Acts::UnitConstants::mm;
+  /// Maximum absolute cot(theta); the default corresponds approximately to |eta| = 4.
+  float cotThetaMax = std::sinh(4.0f);
+
+  /// Minimum and maximum radial separation for bottom-to-middle space-point pairs.
   float deltaRMinBottom = 5.f * Acts::UnitConstants::mm;   ///< Min deltaR for bottom-middle
   float deltaRMaxBottom = 200.f * Acts::UnitConstants::mm; ///< Max deltaR for bottom-middle
-  float deltaRMinTop = 5.f * Acts::UnitConstants::mm;      ///< Min deltaR for middle-top
-  float deltaRMaxTop = 200.f * Acts::UnitConstants::mm;    ///< Max deltaR for middle-top
+  /// Minimum and maximum radial separation for middle-to-top space-point pairs.
+  float deltaRMinTop = 5.f * Acts::UnitConstants::mm;   ///< Min deltaR for middle-top
+  float deltaRMaxTop = 200.f * Acts::UnitConstants::mm; ///< Max deltaR for middle-top
 
-  // Z cuts
+  /// Longitudinal range used to select space points for seeding.
   float zMin = -3000.f * Acts::UnitConstants::mm;
   float zMax = 3000.f * Acts::UnitConstants::mm;
 
-  // Collision region
+  /// Longitudinal range of the collision region used when constructing seeds.
   float collisionRegionMin = -150.f * Acts::UnitConstants::mm;
   float collisionRegionMax = 150.f * Acts::UnitConstants::mm;
 
-  // Quality cuts
+  /// Maximum number of seeds retained for each middle space point.
   float maxSeedsPerMiddleSP = 2;
-  float deltaPhiMax = 0.1f; ///< Maximum phi difference for doublets
+  /// Reserved maximum azimuthal angle difference for doublets (radians); currently unused.
+  float deltaPhiMax = 0.1f;
+
+  /// Enable diagnostic ROOT trees and histograms; disabled by default.
+  bool writeDebugTrees = false;
+  /// Write a JSON snapshot of the tracking geometry on the first use.
+  bool dumpTrackingGeometry = false;
+  /// Number of time frames with debug output; 0 disables it and -1 means all.
+  int debugMaxTimeFrames = 1;
+  /// Destination ROOT file for diagnostic output.
+  std::string debugFileName = "TrackerACTS_debug.root";
+  /// Destination JSON file for the tracking-geometry snapshot.
+  std::string geometryDumpFileName = "TrackerACTS_geometry.json";
 };
 
-/// Space point representation for tracking
+/// Cluster-derived space point passed to the ACTS seeding algorithms.
 struct SpacePoint {
+  /// Cartesian position in the global tracking frame, in mm.
   float x{0.f};
   float y{0.f};
   float z{0.f};
+  /// Detector layer and cluster indices identifying the source measurement.
   int layer{-1};
   int clusterId{-1};
+  /// Readout-frame and sensor indices identifying the source measurement.
   int rof{-1};
+  int sensorId{-1};
+  /// Geometry surface associated with this measurement; not owned by the space point.
+  const Acts::Surface* surface{nullptr};
 
-  // Derived quantities
+  /// Transverse radius in mm.
   float r() const { return std::hypot(x, y); }
-  float radius() const { return r(); } // required by Acts::CylindricalGridElement concept
+  /// Radius accessor required by the ACTS cylindrical-grid element concept.
+  float radius() const { return r(); }
+  /// Azimuthal angle in radians.
   float phi() const { return std::atan2(y, x); }
 
-  // Variance estimates (can be refined based on cluster properties)
-  float varianceR{0.01f}; // ~100 um resolution squared
+  /// Estimated radial and longitudinal position variances, in mm^2.
+  float varianceR{0.01f};
   float varianceZ{0.01f};
 };
 
-/// Seed (triplet of space points)
+/// Three-space-point seed ordered from the innermost to outermost measurement.
 struct SeedACTS {
   const SpacePoint* bottom{nullptr};
   const SpacePoint* middle{nullptr};
   const SpacePoint* top{nullptr};
+  /// ACTS seed quality value.
   float quality{0.f};
 };
 
-/// TRK Tracker using ACTS algorithms for seeding and track finding
+/// TRK tracker using ACTS algorithms for seeding and track finding.
+///
+/// The caller supplies a TimeFrame and, optionally, sensor-ID mappings before
+/// calling clustersToTracks(). The adopted data is referenced rather than
+/// copied, so it must remain alive while the tracker processes it.
 template <int nLayers>
 class TrackerACTS
 {
  public:
   TrackerACTS();
-  ~TrackerACTS()
-  {
-    if (mHistSpacePoints) {
-      mHistSpacePoints->SaveAs("/tmp/mHistSpacePoints.C");
-      delete mHistSpacePoints;
-      mHistSpacePoints = nullptr;
-    }
-  }
+  ~TrackerACTS();
 
-  /// Adopt a TimeFrame for processing
+  /// Set the TimeFrame whose clusters will be reconstructed.
+  /// The TimeFrame is borrowed and must outlive its use by this tracker.
   void adoptTimeFrame(o2::its::TimeFrame<nLayers>& tf);
+  /// Set per-layer sensor-ID mappings used to associate clusters with geometry.
+  /// The mapping array is borrowed and must remain alive while it is used.
+  void adoptSensorIDs(const std::array<std::vector<int>, nLayers>& sensorIDs)
+  {
+    mSensorIDs = &sensorIDs;
+  }
 
-  /// Main tracking entry point: convert clusters to tracks
+  /// Build space points, find tracks, and publish accepted tracks in the TimeFrame.
   void clustersToTracks();
 
-  /// Configuration
+  /// Replace the current tracker configuration.
   void setConfig(const TrackerACTSConfig& cfg) { mConfig = cfg; }
+  /// Access the tracker configuration for modification.
   TrackerACTSConfig& getConfig() { return mConfig; }
+  /// Read-only access to the tracker configuration.
   const TrackerACTSConfig& getConfig() const { return mConfig; }
 
-  /// Set the magnetic field strength
+  /// Set the magnetic field along z, in Tesla.
   void setBz(float bz) { mBz = bz; }
 
-  /// Get the magnetic field strength
+  /// Get the magnetic field along z, in Tesla.
   float getBz() const { return mBz; }
 
-  /// Print tracking summary
+  /// Print the accumulated tracking timing and processing summary.
   void printSummary() const;
 
  private:
-  TH2F* mHistSpacePoints = nullptr;
+  // Debug resources are allocated only if requested and closed at destruction.
+  std::unique_ptr<TFile> mDebugFile;
+  TTree* mClusterTree = nullptr; // Owned by mDebugFile
+  TTree* mSeedTree = nullptr;    // Owned by mDebugFile
+  TTree* mROFTree = nullptr;     // Owned by mDebugFile
+  std::unique_ptr<TH2F> mHistSpacePoints;
+  bool mGeometryDumped = false;
+
+  struct DebugClusterRecord {
+    int tf = 0, rof = 0, layer = 0, clusterId = -1;
+    int externalId = -1, sensorId = -1, mappingStatus = 0;
+    unsigned long long geometryId = 0;
+    float x = 0.f, y = 0.f, z = 0.f;
+    float r = 0.f, phi = 0.f;
+    float varianceR = 0.f, varianceZ = 0.f;
+  } mDebugCluster;
+
+  struct DebugSeedRecord {
+    int tf = 0, rof = 0, seedId = 0;
+    float quality = 0.f;
+    int layer[3] = {}, clusterId[3] = {}, sensorId[3] = {};
+    unsigned long long geometryId[3] = {};
+    float x[3] = {}, y[3] = {}, z[3] = {}, r[3] = {}, phi[3] = {};
+  } mDebugSeed;
+
+  struct DebugROFRecord {
+    int tf = 0, rof = 0, nSpacePoints = 0, nSeeds = 0;
+    int nAcceptedTracks = 0;
+  } mDebugROF;
+
+  void initDebugOutput();
+  void dumpGeometryOnce();
+  bool isDebugTimeFrame() const;
+  void fillDebugSeeds(int rof);
 
   /// Build space points from clusters in the TimeFrame
   void buildSpacePoints(int rof);
 
   /// Create seeds (triplets) from space points using ACTS SeedFinder
   void createSeeds();
 
-  /// Estimate track parameters from a seed using ACTS
-  bool estimateTrackParams(const SeedACTS& seed, o2::its::TrackITSExt& track) const;
-
   /// Run track finding from seeds
   void findTracks();
 
-  /// Assign MC labels to tracks
+  /// Convert fitted ACTS states to a TrackITSExt and publish it in the TimeFrame.
+  bool acceptTrack(const Acts::BoundTrackParameters& inner,
+                   const Acts::BoundTrackParameters& outer,
+                   const SeedACTS& seed, float chi2);
+  /// Assign one MC label per accepted track, after all ROFs have been processed.
   void computeTracksMClabels();
 
   /// Helper: time a task
@@ -150,12 +232,17 @@ class TrackerACTS
 
   // Seeds
   std::vector<SeedACTS> mSeeds;
+  // Avoid accepting seeds which reuse a cluster already assigned to a track.
+  std::unordered_set<uint64_t> mAcceptedClusters;
+  // Layer-local indices in mUnsortedClusters, in the same order as getTracks().
+  // Keep them separate from the external indices stored in TrackITSExt.
+  std::vector<std::array<int, nLayers>> mAcceptedTrackClusterIndices;
 
   // Tracking state
   float mBz{0.5f}; ///< Magnetic field in Tesla
   unsigned int mTimeFrameCounter{0};
   double mTotalTime{0.};
-
+  const std::array<std::vector<int>, nLayers>* mSensorIDs = nullptr;
   // Tracking states for logging
   enum State {
     SpacePointBuilding = 0,