Back to home page

EIC code displayed by LXR

 
 

    


File indexing completed on 2026-05-27 07:24:14

0001 // This file is part of the ACTS project.
0002 //
0003 // Copyright (C) 2016 CERN for the benefit of the ACTS project
0004 //
0005 // This Source Code Form is subject to the terms of the Mozilla Public
0006 // License, v. 2.0. If a copy of the MPL was not distributed with this
0007 // file, You can obtain one at https://mozilla.org/MPL/2.0/.
0008 
0009 #pragma once
0010 
0011 // TODO: Remove this when gcc fixes their false positives.
0012 #if defined(__GNUC__) && !defined(__clang__)
0013 #pragma GCC diagnostic warning "-Wmaybe-uninitialized"
0014 #endif
0015 
0016 // Project include(s)
0017 #include "detray/definitions/algebra.hpp"
0018 #include "detray/definitions/math.hpp"
0019 #include "detray/navigation/detail/print_state.hpp"
0020 #include "detray/navigation/intersection/intersection.hpp"
0021 #include "detray/navigation/navigation_config.hpp"
0022 #include "detray/propagator/base_actor.hpp"
0023 #include "detray/propagator/base_stepper.hpp"
0024 #include "detray/propagator/detail/print_stepper_state.hpp"
0025 #include "detray/propagator/stepping_config.hpp"
0026 #include "detray/tracks/ray.hpp"
0027 #include "detray/utils/tuple_helpers.hpp"
0028 
0029 // System include(s)
0030 #include <iomanip>
0031 #include <limits>
0032 #include <sstream>
0033 #include <string>
0034 #include <tuple>
0035 #include <type_traits>
0036 #include <utility>
0037 
0038 namespace detray {
0039 
0040 /// An inspector that aggregates a number of different inspectors.
0041 template <typename... Inspectors>
0042 struct aggregate_inspector {
0043   using view_type = dmulti_view<detail::get_view_t<Inspectors>...>;
0044   using const_view_type = dmulti_view<detail::get_view_t<const Inspectors>...>;
0045 
0046   using inspector_tuple_t = std::tuple<Inspectors...>;
0047   inspector_tuple_t _inspectors{};
0048 
0049   /// Default constructor
0050   aggregate_inspector() = default;
0051 
0052   /// Construct from the inspector @param view type. Mainly used device-side.
0053   template <concepts::device_view view_t>
0054   DETRAY_HOST_DEVICE explicit aggregate_inspector(view_t &view)
0055       : _inspectors(unroll_views(
0056             view, std::make_index_sequence<sizeof...(Inspectors)>{})) {}
0057 
0058   /// Inspector interface
0059   template <unsigned int current_id = 0, typename state_type,
0060             concepts::point3D point3_t, concepts::vector3D vector3_t,
0061             typename... Args>
0062   DETRAY_HOST_DEVICE auto operator()(state_type &state,
0063                                      const navigation::config &cfg,
0064                                      const point3_t &pos, const vector3_t &dir,
0065                                      const char *message, Args &&...args) {
0066     // Call inspector
0067     std::get<current_id>(_inspectors)(state, cfg, pos, dir, message,
0068                                       std::forward<Args>(args)...);
0069 
0070     // Next inspector
0071     if constexpr (current_id < std::tuple_size<inspector_tuple_t>::value - 1) {
0072       return operator()<current_id + 1>(state, cfg, pos, dir, message,
0073                                         std::forward<Args>(args)...);
0074     }
0075   }
0076 
0077   /// Inspector interface
0078   template <unsigned int current_id = 0, typename state_type, typename... Args>
0079   DETRAY_HOST_DEVICE auto operator()(state_type &state, const char *message,
0080                                      Args &&...args) {
0081     // Call inspector
0082     std::get<current_id>(_inspectors)(state, message);
0083 
0084     // Next inspector
0085     if constexpr (current_id < std::tuple_size<inspector_tuple_t>::value - 1) {
0086       return operator()<current_id + 1>(state, message,
0087                                         std::forward<Args>(args)...);
0088     }
0089   }
0090 
0091   /// @returns a specific inspector by type
0092   template <typename inspector_t>
0093   DETRAY_HOST_DEVICE constexpr decltype(auto) get() {
0094     return std::get<inspector_t>(_inspectors);
0095   }
0096 
0097   /// @returns a specific inspector by type - const
0098   template <typename inspector_t>
0099   DETRAY_HOST_DEVICE constexpr decltype(auto) get() const {
0100     return std::get<inspector_t>(_inspectors);
0101   }
0102 
0103   /// @returns a specific inspector by index
0104   template <std::size_t I>
0105   DETRAY_HOST_DEVICE constexpr decltype(auto) get() {
0106     return std::get<I>(_inspectors);
0107   }
0108 
0109   /// @returns a specific inspector by index - const
0110   template <std::size_t I>
0111   DETRAY_HOST_DEVICE constexpr decltype(auto) get() const {
0112     return std::get<I>(_inspectors);
0113   }
0114 
0115   /// @returns a tuple constructed from the inspector @param view s.
0116   template <concepts::device_view view_t, std::size_t... I>
0117   DETRAY_HOST_DEVICE constexpr auto unroll_views(
0118       view_t &view, std::index_sequence<I...> /*seq*/) {
0119     return detail::make_tuple<std::tuple>(
0120         Inspectors(detail::get<I>(view.m_view))...);
0121   }
0122 };
0123 
0124 namespace navigation {
0125 
0126 namespace detail {
0127 
0128 /// Record of a surface intersection along a track
0129 template <typename intersetion_t>
0130 struct candidate_record {
0131   using algebra_type = typename intersetion_t::algebra_type;
0132   using scalar_type = dscalar<algebra_type>;
0133   using point3_type = dpoint3D<algebra_type>;
0134   using vector3_type = dvector3D<algebra_type>;
0135   using intersection_type = intersetion_t;
0136 
0137   constexpr candidate_record() = default;
0138 
0139   /// The particle charge is not known in the navigation, but might be
0140   /// provided in a different context
0141   DETRAY_HOST_DEVICE
0142   constexpr candidate_record(
0143       const point3_type &position, const vector3_type &direction,
0144       const intersection_type &intr,
0145       const scalar_type q = detray::detail::invalid_value<scalar_type>(),
0146       const scalar_type p = detray::detail::invalid_value<scalar_type>())
0147       : pos{position},
0148         dir{direction},
0149         intersection{intr},
0150         charge{q},
0151         p_mag{p} {}
0152 
0153   /// Current global track position
0154   point3_type pos{0.f, 0.f, 0.f};
0155   /// Current global track direction
0156   vector3_type dir{0.f, 0.f, 1.f};
0157   /// The intersection result
0158   intersetion_t intersection{};
0159   /// Charge hypothesis of the particle (invalid value if not known)
0160   scalar_type charge{detray::detail::invalid_value<scalar_type>()};
0161   /// Current momentum magnitude of the particle
0162   scalar_type p_mag{1.f};
0163 };
0164 
0165 }  // namespace detail
0166 
0167 /// A navigation inspector that relays information about the encountered
0168 /// objects whenever the navigator reaches one or more status flags
0169 template <typename candidate_t, template <typename...> class vector_t = dvector,
0170           status... navigation_status>
0171 struct object_tracer {
0172   using candidate_record_t = detail::candidate_record<candidate_t>;
0173   using scalar_t = typename candidate_record_t::scalar_type;
0174 
0175   using view_type = dvector_view<candidate_record_t>;
0176   using const_view_type = dvector_view<const candidate_record_t>;
0177 
0178   /// Default constructor
0179   object_tracer() = default;
0180 
0181   /// Device-side construction from a vecmem based view type
0182   DETRAY_HOST_DEVICE explicit object_tracer(
0183       dvector_view<candidate_record_t> &view)
0184       : object_trace(view) {}
0185 
0186   // Record all objects the navigator encounters
0187   vector_t<candidate_record_t> object_trace;
0188   dindex current_vol{dindex_invalid};
0189   const scalar_t inv_pos{detray::detail::invalid_value<scalar_t>()};
0190   typename candidate_record_t::point3_type last_pos = {inv_pos, inv_pos,
0191                                                        inv_pos};
0192   typename candidate_record_t::vector3_type last_dir = {0.f, 0.f, 0.f};
0193 
0194   /// Inspector interface
0195   template <typename state_type, concepts::point3D point3_t,
0196             concepts::vector3D vector3_t, typename... Args>
0197   DETRAY_HOST_DEVICE auto operator()(const state_type &state,
0198                                      const navigation::config & /*unused*/,
0199                                      const point3_t &pos, const vector3_t &dir,
0200                                      const char * /*message*/,
0201                                      Args &&.../*unused*/) {
0202     // Record the candidate of an encountered object
0203     if ((is_status(state.status(), navigation_status) || ...)) {
0204       // Reached a new position: log it
0205       // Also log volume switches that happen without position update
0206       if ((vector::norm(last_pos - pos) >=
0207            10.f * std::numeric_limits<scalar_t>::epsilon()) ||
0208           (state.is_on_portal() && current_vol != state.volume()) ||
0209           state.current().surface().identifier() !=
0210               object_trace.back().intersection.surface().identifier()) {
0211         object_trace.push_back({pos, dir, state.current()});
0212         last_pos = pos;
0213         last_dir = dir;
0214         current_vol = state.volume();
0215       }
0216     }
0217   }
0218 
0219   /// Inspector interface
0220   template <typename state_type>
0221   DETRAY_HOST_DEVICE auto operator()(const state_type & /*state*/,
0222                                      const char * /*message*/) { /* Do nothing*/
0223   }
0224 
0225   /// @returns a specific candidate from the trace
0226   DETRAY_HOST_DEVICE
0227   constexpr const candidate_record_t &operator[](std::size_t i) const {
0228     return object_trace[i];
0229   }
0230 
0231   /// @returns the full object trace
0232   DETRAY_HOST_DEVICE
0233   constexpr const auto &trace() const { return object_trace; }
0234 
0235   /// Compares a navigation status with the tracers references
0236   DETRAY_HOST_DEVICE
0237   constexpr bool is_status(const status &nav_stat, const status &ref_stat) {
0238     return (nav_stat == ref_stat);
0239   }
0240 };
0241 
0242 /// A navigation inspector that prints information about the current navigation
0243 /// state. Meant for debugging.
0244 struct print_inspector {
0245   using view_type = dvector_view<char>;
0246   using const_view_type = dvector_view<const char>;
0247 
0248   struct void_generator {};
0249 
0250   /// Default constructor
0251   print_inspector() = default;
0252 
0253   /// Copy constructor ensures that the string stream is set up identically
0254   print_inspector(const print_inspector &other)
0255       : debug_stream(other.debug_stream.str()) {}
0256 
0257   /// Move constructor
0258   print_inspector(print_inspector &&other) = default;
0259 
0260   /// Default destructor
0261   ~print_inspector() = default;
0262 
0263   /// Copy assignemten operator ensures that the string stream is set up
0264   /// identically
0265   print_inspector &operator=(const print_inspector &other) {
0266     // Reset
0267     debug_stream.str(std::string());
0268     debug_stream.clear();
0269     // Move new debug string in
0270     debug_stream << other.debug_stream.str();
0271     return *this;
0272   }
0273 
0274   /// Move assignemten operator
0275   print_inspector &operator=(print_inspector &&other) = default;
0276 
0277   /// Inspector interface. Gathers detailed information during navigation
0278   template <typename state_type, concepts::point3D point3_t,
0279             concepts::vector3D vector3_t,
0280             typename message_generator_t = void_generator>
0281   auto operator()(const state_type &state, const navigation::config &cfg,
0282                   const point3_t &track_pos, const vector3_t &track_dir,
0283                   const char *message,
0284                   const message_generator_t &msg_gen = {}) {
0285     std::string msg(message);
0286     debug_stream << msg;
0287     if constexpr (!std::same_as<message_generator_t, void_generator>) {
0288       debug_stream << msg_gen();
0289 
0290       if (state.status() == navigation::status::e_abort) {
0291         fata_error_msg = msg_gen();
0292       }
0293     }
0294     debug_stream << std::endl;
0295     debug_stream << "----------------------------------------" << std::endl;
0296 
0297     debug_stream << navigation::print_state(state);
0298     debug_stream << navigation::print_candidates(state, cfg, track_pos,
0299                                                  track_dir);
0300 
0301     debug_stream << std::endl << std::endl;
0302   }
0303 
0304   /// Inspector interface. Print basic state information
0305   template <typename state_type, typename message_generator_t = void_generator>
0306   auto operator()(const state_type &state, const char *message,
0307                   const message_generator_t &msg_gen = {}) {
0308     std::string msg(message);
0309     debug_stream << msg;
0310     if constexpr (!std::same_as<message_generator_t, void_generator>) {
0311       debug_stream << msg_gen();
0312 
0313       if (state.status() == navigation::status::e_abort) {
0314         fata_error_msg = msg_gen();
0315       }
0316     }
0317     debug_stream << std::endl;
0318     debug_stream << "----------------------------------------" << std::endl;
0319 
0320     debug_stream << navigation::print_state(state);
0321 
0322     debug_stream << std::endl << std::endl;
0323   }
0324 
0325   /// @returns a string representation of the gathered information
0326   std::string to_string() const { return debug_stream.str(); }
0327 
0328   /// Gathers navigation information across navigator update calls
0329   std::stringstream debug_stream{};
0330   /// Special message that is collected if the navigator hits a fatal error
0331   std::string fata_error_msg{""};
0332 };
0333 
0334 }  // namespace navigation
0335 
0336 namespace stepping {
0337 
0338 /// A stepper inspector that prints information about the current stepper
0339 /// state. Meant for debugging.
0340 struct print_inspector {
0341   /// Default constructor
0342   print_inspector() = default;
0343 
0344   /// Copy constructor ensures that the string stream is set up identically
0345   print_inspector(const print_inspector &other)
0346       : debug_stream(other.debug_stream.str()) {}
0347 
0348   /// Copy assignemten operator ensures that the string stream is set up
0349   /// identically
0350   print_inspector &operator=(const print_inspector &other) {
0351     // Reset
0352     debug_stream.str(std::string());
0353     debug_stream.clear();
0354     // Move new debug string in
0355     debug_stream << other.debug_stream.str();
0356     return *this;
0357   }
0358 
0359   /// Move assignemten operator
0360   print_inspector &operator=(print_inspector &&other) = default;
0361 
0362   /// Gathers stepping information from inside the stepper methods
0363   std::stringstream debug_stream{};
0364 
0365   /// Inspector interface. Gathers detailed information during stepping
0366   template <typename state_type, concepts::scalar scalar_t>
0367   void operator()(const state_type &state, const stepping::config & /*unused*/,
0368                   const char *message, const scalar_t dist) {
0369     std::string msg(message);
0370 
0371     debug_stream << msg << std::endl;
0372     debug_stream << detray::stepping::print_state(state, dist);
0373   }
0374 
0375   /// Inspector interface. Gathers detailed information during stepping
0376   template <typename state_type, concepts::scalar scalar_t>
0377   void operator()(const state_type &state, const stepping::config & /*unused*/,
0378                   const char *message, const scalar_t /*dist*/,
0379                   const std::size_t n_trials, const scalar_t step_scalor) {
0380     std::string msg(message);
0381 
0382     debug_stream << msg << std::endl;
0383     debug_stream << detray::stepping::print_state(state, n_trials, step_scalor);
0384   }
0385 
0386   /// @returns a string representation of the gathered information
0387   std::string to_string() const { return debug_stream.str(); }
0388 };
0389 
0390 }  // namespace stepping
0391 
0392 }  // namespace detray