Merge pull request #950 from niermann999/test-detray-main2

stephenswat · web-flow · commit 7911e66ca08e · 2025-04-30T17:46:40.000+02:00
Update to detray version 0.95.0
diff --git a/core/include/traccc/finding/actors/ckf_aborter.hpp b/core/include/traccc/finding/actors/ckf_aborter.hpp
@@ -43,10 +43,10 @@ struct ckf_aborter : detray::actor {
         abrt_state.count++;
         abrt_state.path_from_surface += stepping.step_size();
 
-        // Abort at the next sensitive surface
+        // Stop at the next sensitive surface
         if (navigation.is_on_sensitive() &&
             abrt_state.path_from_surface > abrt_state.min_step_length) {
-            prop_state._heartbeat &= navigation.abort();
+            prop_state._heartbeat &= navigation.pause();
             abrt_state.success = true;
         }
 
@@ -57,7 +57,9 @@ struct ckf_aborter : detray::actor {
         }
 
         if (abrt_state.count > abrt_state.max_count) {
-            prop_state._heartbeat &= navigation.abort();
+            prop_state._heartbeat &= navigation.abort(
+                "CKF: Maximum number of steps to reach next sensitive surface "
+                "exceeded");
         }
     }
 };
diff --git a/core/include/traccc/finding/details/find_tracks.hpp b/core/include/traccc/finding/details/find_tracks.hpp
@@ -73,7 +73,8 @@ track_candidate_container_types::host find_tracks(
 
     using actor_type = detray::actor_chain<
         detray::pathlimit_aborter<scalar_type>, transporter_type,
-        interaction_register<interactor_type>, interactor_type, ckf_aborter>;
+        interaction_register<interactor_type>, interactor_type,
+        detray::momentum_aborter<scalar_type>, ckf_aborter>;
 
     using propagator_type =
         detray::propagator<stepper_t, navigator_t, actor_type>;
@@ -314,12 +315,18 @@ track_candidate_container_types::host find_tracks(
                     config.propagation.stepping.step_constraint);
 
             typename detray::pathlimit_aborter<scalar_type>::state s0;
-            typename detray::parameter_transporter<algebra_type>::state s1;
-            typename interactor_type::state s3;
-            typename interaction_register<interactor_type>::state s2{s3};
+            typename interactor_type::state s2;
+            typename interaction_register<interactor_type>::state s1{s2};
+            typename detray::momentum_aborter<scalar_type>::state s3{};
             typename ckf_aborter::state s4;
+            // Update the actor config
             s4.min_step_length = config.min_step_length_for_next_surface;
             s4.max_count = config.max_step_counts_for_next_surface;
+            if (config.is_min_pT) {
+                s3.min_pT(static_cast<scalar_type>(config.min_p_mag));
+            } else {
+                s3.min_p(static_cast<scalar_type>(config.min_p_mag));
+            }
 
             // Propagate to the next surface
             propagator.propagate_sync(propagation,
@@ -328,6 +335,7 @@ track_candidate_container_types::host find_tracks(
             // If a surface found, add the parameter for the next
             // step
             if (s4.success) {
+                assert(propagation._navigation.is_on_sensitive());
                 out_params.push_back(propagation._stepping.bound_params());
                 param_to_link[step].push_back(link_id);
             }
diff --git a/core/include/traccc/finding/finding_config.hpp b/core/include/traccc/finding/finding_config.hpp
@@ -45,6 +45,10 @@ struct finding_config {
     /// Propagation configuration
     detray::propagation::config propagation{};
 
+    /// Minimum momentum for reconstructed tracks
+    bool is_min_pT = true;
+    float min_p_mag = 100.f * traccc::unit<float>::MeV;
+
     /// Particle hypothesis
     traccc::pdg_particle<traccc::scalar> ptc_hypothesis =
         traccc::muon<traccc::scalar>();
@@ -64,6 +68,20 @@ struct finding_config {
     /// @note This parameter affects GPU-based track finding only.
     unsigned int initial_links_per_seed = 100;
     /// @}
+
+    /// Set the momentum limit to @param p
+    TRACCC_HOST_DEVICE
+    inline void min_p(const float p) {
+        is_min_pT = false;
+        min_p_mag = p;
+    }
+
+    /// Set the transverse momentum limit to @param p
+    TRACCC_HOST_DEVICE
+    inline void min_pT(const float p) {
+        is_min_pT = true;
+        min_p_mag = p;
+    }
 };
 
 }  // namespace traccc
diff --git a/core/include/traccc/fitting/kalman_filter/kalman_actor.hpp b/core/include/traccc/fitting/kalman_filter/kalman_actor.hpp
@@ -109,7 +109,7 @@ struct kalman_actor : detray::actor {
 
         // If the iterator reaches the end, terminate the propagation
         if (actor_state.is_complete()) {
-            propagation._heartbeat &= navigation.abort();
+            propagation._heartbeat &= navigation.exit();
             return;
         }
 
@@ -153,7 +153,8 @@ struct kalman_actor : detray::actor {
 
             // Abort if the Kalman update fails
             if (res != kalman_fitter_status::SUCCESS) {
-                propagation._heartbeat &= navigation.abort();
+                propagation._heartbeat &=
+                    navigation.abort(fitter_debug_msg{res});
                 return;
             }
 
diff --git a/core/include/traccc/fitting/kalman_filter/kalman_step_aborter.hpp b/core/include/traccc/fitting/kalman_filter/kalman_step_aborter.hpp
@@ -56,7 +56,9 @@ struct kalman_step_aborter : public detray::actor {
         // Abort if the step count exceeds the maximum allowed
         if (++(abrt_state.step) > abrt_state.max_steps) {
             VECMEM_DEBUG_MSG(1, "Kalman fitter step aborter triggered");
-            prop_state._heartbeat &= navigation.abort();
+            prop_state._heartbeat &= navigation.abort(
+                "Kalman Fitter: Maximum number of steps to reach next "
+                "sensitive surface exceeded");
         }
     }
 
diff --git a/core/include/traccc/fitting/status_codes.hpp b/core/include/traccc/fitting/status_codes.hpp
@@ -8,6 +8,8 @@
 #pragma once
 
 #include <cstdint>
+#include <string>
+
 namespace traccc {
 enum class kalman_fitter_status : uint32_t {
     SUCCESS,
@@ -22,4 +24,47 @@ enum class kalman_fitter_status : uint32_t {
     ERROR_OTHER,
     MAX_STATUS
 };
+
+/// Convert a status code into a human readable string
+struct fitter_debug_msg {
+
+    TRACCC_HOST std::string operator()() const {
+        const std::string msg{"Kalman Fitter: "};
+        switch (m_error_code) {
+            using enum kalman_fitter_status;
+            case ERROR_QOP_ZERO: {
+                return msg + "Track qop is zero";
+            }
+            case ERROR_THETA_ZERO: {
+                return msg + "Track theta is zero";
+            }
+            case ERROR_INVERSION: {
+                return msg + "Failed matrix inversion";
+            }
+            case ERROR_SMOOTHER_CHI2_NEGATIVE: {
+                return msg + "Negative chi2 in smoother";
+            }
+            case ERROR_SMOOTHER_CHI2_NOT_FINITE: {
+                return msg + "Invalid chi2 in smoother";
+            }
+            case ERROR_UPDATER_CHI2_NEGATIVE: {
+                return msg + "Negative chi2 in gain matrix update";
+            }
+            case ERROR_UPDATER_CHI2_NOT_FINITE: {
+                return msg + "Invalid chi2 in gain matrix update";
+            }
+            case ERROR_BARCODE_SEQUENCE_OVERFLOW: {
+                return msg + "Barcode sequence overflow in direct navigator";
+            }
+            case ERROR_OTHER: {
+                return msg + "Unspecified error";
+            }
+            default: {
+                return "";
+            }
+        }
+    }
+
+    kalman_fitter_status m_error_code{kalman_fitter_status::MAX_STATUS};
+};
 }  // namespace traccc
diff --git a/core/include/traccc/seeding/impl/spacepoint_formation.ipp b/core/include/traccc/seeding/impl/spacepoint_formation.ipp
@@ -30,7 +30,7 @@ TRACCC_HOST_DEVICE inline void fill_pixel_spacepoint(edm::spacepoint<soa_t>& sp,
 
     // Get the global position of this silicon pixel measurement.
     const detray::tracking_surface sf{det, meas.surface_link};
-    const auto global = sf.bound_to_global({}, meas.local, {});
+    const auto global = sf.local_to_global({}, meas.local, {});
 
     // Fill the spacepoint with the global position and the measurement.
     sp.x() = global[0];
diff --git a/device/alpaka/include/traccc/alpaka/finding/finding_algorithm.hpp b/device/alpaka/include/traccc/alpaka/finding/finding_algorithm.hpp
@@ -54,7 +54,7 @@ class finding_algorithm
         detray::actor_chain<detray::pathlimit_aborter<scalar_type>,
                             detray::parameter_transporter<algebra_type>,
                             interaction_register<interactor>, interactor,
-                            ckf_aborter>;
+                            detray::momentum_aborter<scalar_type>, ckf_aborter>;
 
     using propagator_type =
         detray::propagator<stepper_t, navigator_t, actor_type>;
diff --git a/device/common/include/traccc/finding/device/impl/propagate_to_next_surface.ipp b/device/common/include/traccc/finding/device/impl/propagate_to_next_surface.ipp
@@ -21,6 +21,8 @@ TRACCC_HOST_DEVICE inline void propagate_to_next_surface(
     const global_index_t globalIndex, const finding_config& cfg,
     const propagate_to_next_surface_payload<propagator_t, bfield_t>& payload) {
 
+    using scalar_t = propagator_t::detector_type::scalar_type;
+
     if (globalIndex >= payload.n_in_params) {
         return;
     }
@@ -100,14 +102,23 @@ TRACCC_HOST_DEVICE inline void propagate_to_next_surface(
     typename detray::detail::tuple_element<2, actor_tuple_type>::type::state s2{
         s3};
     typename detray::detail::tuple_element<4, actor_tuple_type>::type::state s4;
-    s4.min_step_length = cfg.min_step_length_for_next_surface;
-    s4.max_count = cfg.max_step_counts_for_next_surface;
+
+    typename detray::detail::tuple_element<5, actor_tuple_type>::type::state s5;
+    s5.min_step_length = cfg.min_step_length_for_next_surface;
+    s5.max_count = cfg.max_step_counts_for_next_surface;
+    if (cfg.is_min_pT) {
+        s4.min_pT(static_cast<scalar_t>(cfg.min_p_mag));
+    } else {
+        s4.min_p(static_cast<scalar_t>(cfg.min_p_mag));
+    }
 
     // Propagate to the next surface
-    propagator.propagate_sync(propagation, detray::tie(s0, s2, s3, s4));
+    propagator.propagate_sync(propagation, detray::tie(s0, s2, s3, s4, s5));
 
     // If a surface found, add the parameter for the next step
-    if (s4.success) {
+    if (s5.success) {
+        assert(propagation._navigation.is_on_sensitive());
+
         params[param_id] = propagation._stepping.bound_params();
 
         if (payload.step == cfg.max_track_candidates_per_track - 1) {
diff --git a/device/cuda/include/traccc/cuda/finding/finding_algorithm.hpp b/device/cuda/include/traccc/cuda/finding/finding_algorithm.hpp
@@ -56,7 +56,7 @@ class finding_algorithm
         detray::actor_chain<detray::pathlimit_aborter<scalar_type>,
                             detray::parameter_transporter<algebra_type>,
                             interaction_register<interactor>, interactor,
-                            ckf_aborter>;
+                            detray::momentum_aborter<scalar_type>, ckf_aborter>;
 
     using propagator_type =
         detray::propagator<stepper_t, navigator_t, actor_type>;
diff --git a/device/sycl/src/finding/find_tracks.hpp b/device/sycl/src/finding/find_tracks.hpp
@@ -369,11 +369,11 @@ track_candidate_container_types::buffer find_tracks(
                 typename navigator_t::detector_type::scalar_type;
             using interactor_type =
                 detray::pointwise_material_interactor<algebra_type>;
-            using actor_type =
-                detray::actor_chain<detray::pathlimit_aborter<scalar_type>,
-                                    detray::parameter_transporter<algebra_type>,
-                                    interaction_register<interactor_type>,
-                                    interactor_type, ckf_aborter>;
+            using actor_type = detray::actor_chain<
+                detray::pathlimit_aborter<scalar_type>,
+                detray::parameter_transporter<algebra_type>,
+                interaction_register<interactor_type>, interactor_type,
+                detray::momentum_aborter<scalar_type>, ckf_aborter>;
             using propagator_type =
                 detray::propagator<stepper_t, navigator_t, actor_type>;
 
diff --git a/examples/options/include/traccc/options/track_finding.hpp b/examples/options/include/traccc/options/track_finding.hpp
@@ -26,6 +26,12 @@ class track_finding : public interface, public config_provider<finding_config> {
     /// Constructor
     track_finding();
 
+    /// Read/process the command line options
+    ///
+    /// @param vm The command line options to interpret/read
+    ///
+    void read(const boost::program_options::variables_map& vm) override;
+
     /// Configuration conversion operators
     operator finding_config() const override;
 
diff --git a/examples/options/src/track_finding.cpp b/examples/options/src/track_finding.cpp
@@ -63,6 +63,22 @@ track_finding::track_finding() : interface("Track Finding Options") {
     m_desc.add_options()("particle-hypothesis",
                          po::value(&m_pdg_number)->default_value(m_pdg_number),
                          "PDG number for the particle hypothesis");
+    m_desc.add_options()(
+        "min-total-momentum [GeV]",
+        po::value(&m_config.min_p_mag)->default_value(m_config.min_p_mag),
+        "Minimum total track momentum");
+    m_desc.add_options()(
+        "min-transverse-momentum [GeV]",
+        po::value(&m_config.min_p_mag)->default_value(m_config.min_p_mag),
+        "Minimum transverse track momentum");
+}
+
+void track_finding::read(const po::variables_map &vm) {
+
+    m_config.min_p_mag *= traccc::unit<float>::GeV;
+
+    // If not set as total momentum, interpret as transverse momentum
+    m_config.is_min_pT = vm["min-total-momentum"].defaulted();
 }
 
 track_finding::operator finding_config() const {
@@ -102,6 +118,16 @@ std::unique_ptr<configuration_printable> track_finding::as_printable() const {
         std::to_string(m_config.max_num_skipping_per_cand)));
     cat->add_child(std::make_unique<configuration_kv_pair>(
         "PDG number", std::to_string(m_pdg_number)));
+    // How to interpret the minimum track momentum value
+    if (m_config.is_min_pT) {
+        cat->add_child(std::make_unique<configuration_kv_pair>(
+            "Minimum transverse track momentum",
+            std::to_string(m_config.min_p_mag)));
+    } else {
+        cat->add_child(std::make_unique<configuration_kv_pair>(
+            "Minimum total track momentum",
+            std::to_string(m_config.min_p_mag)));
+    }
 
     return cat;
 }
diff --git a/extern/algebra-plugins/CMakeLists.txt b/extern/algebra-plugins/CMakeLists.txt
@@ -13,7 +13,7 @@ message( STATUS "Building Algebra Plugins as part of the TRACCC project" )
 
 # Declare where to get Algebra Plugins from.
 set( TRACCC_ALGEBRA_PLUGINS_SOURCE
-   "URL;https://github.com/acts-project/algebra-plugins/archive/refs/tags/v0.26.2.tar.gz;URL_MD5;607bd0b8beffcf7981b10c03f237e550"
+   "URL;https://github.com/acts-project/algebra-plugins/archive/refs/tags/v0.27.0.tar.gz;URL_MD5;d59270836aba4cebc27369acea6a995f"
    CACHE STRING "Source for Algebra Plugins, when built as part of this project" )
 mark_as_advanced( TRACCC_ALGEBRA_PLUGINS_SOURCE )
 FetchContent_Declare( AlgebraPlugins SYSTEM ${TRACCC_ALGEBRA_PLUGINS_SOURCE} )
diff --git a/extern/detray/CMakeLists.txt b/extern/detray/CMakeLists.txt
@@ -13,7 +13,7 @@ message( STATUS "Building Detray as part of the TRACCC project" )
 
 # Declare where to get Detray from.
 set( TRACCC_DETRAY_SOURCE
-   "URL;https://github.com/acts-project/detray/archive/refs/tags/v0.94.0.tar.gz;URL_MD5;b9f6c8f0fa5c6df22df1ebe1aff2ab4a"
+   "URL;https://github.com/acts-project/detray/archive/refs/tags/v0.95.0.tar.gz;URL_MD5;e5c89d7fc17f868eaf46e4eb3b8e630c"
    CACHE STRING "Source for Detray, when built as part of this project" )
 mark_as_advanced( TRACCC_DETRAY_SOURCE )
 FetchContent_Declare( Detray SYSTEM ${TRACCC_DETRAY_SOURCE} )
diff --git a/io/src/obj/write_track_candidates.cpp b/io/src/obj/write_track_candidates.cpp
@@ -8,6 +8,9 @@
 // Local include(s).
 #include "write_track_candidates.hpp"
 
+// Detray include(s)
+#include <detray/geometry/tracking_surface.hpp>
+
 // System include(s).
 #include <cassert>
 #include <fstream>
@@ -49,7 +52,7 @@ void write_track_candidates(
             const detray::tracking_surface surface{detector, m.surface_link};
 
             // Calculate a position for this measurement in global 3D space.
-            const auto global = surface.bound_to_global({}, m.local, {});
+            const auto global = surface.local_to_global({}, m.local, {});
 
             // Write the 3D coordinates of the measurement / spacepoint.
             assert(global.size() == 3);
diff --git a/io/src/read_detector_description.cpp b/io/src/read_detector_description.cpp
@@ -13,6 +13,9 @@
 #include "traccc/io/read_digitization_config.hpp"
 #include "traccc/io/utils.hpp"
 
+// Detray include(s)
+#include <detray/geometry/tracking_surface.hpp>
+
 // VecMem include(s).
 #include <vecmem/memory/host_memory_resource.hpp>
 
@@ -117,8 +120,7 @@ void read_json_dd(traccc::silicon_detector_description::host& dd,
         dd.resize(dd.size() + 1);
 
         // Construct a Detray surface object.
-        const detray::tracking_surface<traccc::default_detector::host> surface{
-            detector, surface_desc};
+        const detray::tracking_surface surface{detector, surface_desc};
 
         // Fill the new element with the geometry ID and the transformation of
         // the surface in question.
diff --git a/simulation/include/traccc/simulation/simulator.hpp b/simulation/include/traccc/simulation/simulator.hpp
diff --git a/tests/cuda/test_ckf_toy_detector.cpp b/tests/cuda/test_ckf_toy_detector.cpp

Original file line number	Diff line number	Diff line change
`@@ -43,10 +43,10 @@ struct ckf_aborter : detray::actor {`
`43`	`43`	`abrt_state.count++;`
`44`	`44`	`abrt_state.path_from_surface += stepping.step_size();`
`45`	`45`
`46`		`- // Abort at the next sensitive surface`
	`46`	`+ // Stop at the next sensitive surface`
`47`	`47`	`if (navigation.is_on_sensitive() &&`
`48`	`48`	`abrt_state.path_from_surface > abrt_state.min_step_length) {`
`49`		`- prop_state._heartbeat &= navigation.abort();`
	`49`	`+ prop_state._heartbeat &= navigation.pause();`
`50`	`50`	`abrt_state.success = true;`
`51`	`51`	`}`
`52`	`52`
`@@ -57,7 +57,9 @@ struct ckf_aborter : detray::actor {`
`57`	`57`	`}`
`58`	`58`
`59`	`59`	`if (abrt_state.count > abrt_state.max_count) {`
`60`		`- prop_state._heartbeat &= navigation.abort();`
	`60`	`+ prop_state._heartbeat &= navigation.abort(`
	`61`	`+ "CKF: Maximum number of steps to reach next sensitive surface "`
	`62`	`+ "exceeded");`
`61`	`63`	`}`
`62`	`64`	`}`
`63`	`65`	`};`
Original file line number	Diff line number	Diff line change
`@@ -109,7 +109,7 @@ struct kalman_actor : detray::actor {`
`109`	`109`
`110`	`110`	`// If the iterator reaches the end, terminate the propagation`
`111`	`111`	`if (actor_state.is_complete()) {`
`112`		`- propagation._heartbeat &= navigation.abort();`
	`112`	`+ propagation._heartbeat &= navigation.exit();`
`113`	`113`	`return;`
`114`	`114`	`}`
`115`	`115`
`@@ -153,7 +153,8 @@ struct kalman_actor : detray::actor {`
`153`	`153`
`154`	`154`	`// Abort if the Kalman update fails`
`155`	`155`	`if (res != kalman_fitter_status::SUCCESS) {`
`156`		`- propagation._heartbeat &= navigation.abort();`
	`156`	`+ propagation._heartbeat &=`
	`157`	`+ navigation.abort(fitter_debug_msg{res});`
`157`	`158`	`return;`
`158`	`159`	`}`
`159`	`160`
Original file line number	Diff line number	Diff line change
`@@ -56,7 +56,9 @@ struct kalman_step_aborter : public detray::actor {`
`56`	`56`	`// Abort if the step count exceeds the maximum allowed`
`57`	`57`	`if (++(abrt_state.step) > abrt_state.max_steps) {`
`58`	`58`	`VECMEM_DEBUG_MSG(1, "Kalman fitter step aborter triggered");`
`59`		`- prop_state._heartbeat &= navigation.abort();`
	`59`	`+ prop_state._heartbeat &= navigation.abort(`
	`60`	`+ "Kalman Fitter: Maximum number of steps to reach next "`
	`61`	`+ "sensitive surface exceeded");`
`60`	`62`	`}`
`61`	`63`	`}`
`62`	`64`