diff --git a/CMakeLists.txt b/CMakeLists.txt index 57fcdfc781..2402d64e23 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -52,6 +52,10 @@ include( traccc-functions ) set( TRACCC_CUSTOM_SCALARTYPE "float" CACHE STRING "Scalar type to use in the TRACCC code" ) +# Temporary setting for the traccc device log level, until it can be removed. +set( TRACCC_DEVICE_LOG_LVL "NONE" CACHE STRING + "Log level for traccc and detray device code" ) + # Flags controlling which parts of traccc to build. option( TRACCC_BUILD_CUDA "Build the CUDA sources included in traccc" FALSE ) option( TRACCC_BUILD_HIP "Build the HIP sources included in traccc" FALSE) diff --git a/core/include/traccc/edm/track_parameters.hpp b/core/include/traccc/edm/track_parameters.hpp index 0a33dc3ffe..0f5ec78187 100644 --- a/core/include/traccc/edm/track_parameters.hpp +++ b/core/include/traccc/edm/track_parameters.hpp @@ -14,6 +14,7 @@ #include "traccc/definitions/qualifiers.hpp" #include "traccc/definitions/track_parametrization.hpp" #include "traccc/edm/container.hpp" +#include "traccc/utils/trigonometric_helpers.hpp" // detray include(s). #include @@ -48,19 +49,15 @@ using bound_track_parameters_collection_types = // Wrap the phi of track parameters to [-pi,pi] template -TRACCC_HOST_DEVICE inline void wrap_phi( +TRACCC_HOST_DEVICE constexpr void normalize_angles( bound_track_parameters& param) { + traccc::scalar phi; + traccc::scalar theta; + + std::tie(phi, theta) = detail::wrap_phi_theta(param.phi(), param.theta()); - traccc::scalar phi = param.phi(); - static constexpr traccc::scalar TWOPI = - 2.f * traccc::constant::pi; - phi = math::fmod(phi, TWOPI); - if (phi > traccc::constant::pi) { - phi -= TWOPI; - } else if (phi < -traccc::constant::pi) { - phi += TWOPI; - } param.set_phi(phi); + param.set_theta(theta); } /// Covariance inflation used for track fitting diff --git a/core/include/traccc/finding/actors/ckf_aborter.hpp b/core/include/traccc/finding/actors/ckf_aborter.hpp index 70d6fb8031..fa780667d4 100644 --- a/core/include/traccc/finding/actors/ckf_aborter.hpp +++ b/core/include/traccc/finding/actors/ckf_aborter.hpp @@ -10,8 +10,10 @@ // Project include(s) #include "traccc/definitions/primitives.hpp" #include "traccc/definitions/qualifiers.hpp" +#include "traccc/utils/logging.hpp" // detray include(s) +#include #include // System include(s) @@ -38,10 +40,10 @@ struct ckf_aborter : detray::actor { propagator_state_t &prop_state) const { auto &navigation = prop_state._navigation; - const auto &stepping = prop_state._stepping; abrt_state.count++; - abrt_state.path_from_surface += stepping.step_size(); + + TRACCC_VERBOSE_HOST_DEVICE("Checking CKF aborter"); // Stop at the next sensitive surface if (navigation.is_on_sensitive() && @@ -50,10 +52,8 @@ struct ckf_aborter : detray::actor { abrt_state.success = true; } - // Reset path from surface - if (navigation.is_on_sensitive()) { - abrt_state.path_from_surface = 0.f; - } + TRACCC_VERBOSE_HOST_DEVICE("-> Found sensitive surface: %d", + navigation.barcode().index()); if (abrt_state.count > abrt_state.max_count) { prop_state._heartbeat &= navigation.abort( diff --git a/core/include/traccc/finding/details/combinatorial_kalman_filter.hpp b/core/include/traccc/finding/details/combinatorial_kalman_filter.hpp index 27586ba155..30f7fcba36 100644 --- a/core/include/traccc/finding/details/combinatorial_kalman_filter.hpp +++ b/core/include/traccc/finding/details/combinatorial_kalman_filter.hpp @@ -15,6 +15,7 @@ #include "traccc/finding/candidate_link.hpp" #include "traccc/finding/details/combinatorial_kalman_filter_types.hpp" #include "traccc/finding/finding_config.hpp" +#include "traccc/finding/propagation_data.hpp" #include "traccc/fitting/kalman_filter/gain_matrix_updater.hpp" #include "traccc/fitting/kalman_filter/is_line_visitor.hpp" #include "traccc/fitting/status_codes.hpp" @@ -63,7 +64,7 @@ combinatorial_kalman_filter( const measurement_collection_types::const_view& measurements_view, const bound_track_parameters_collection_types::const_view& seeds_view, const finding_config& config, vecmem::memory_resource& mr, - const Logger& log) { + const Logger& /*log*/) { assert(config.min_step_length_for_next_surface > math::fabs(config.propagation.navigation.overstep_tolerance) && @@ -77,7 +78,8 @@ combinatorial_kalman_filter( using scalar_type = detray::dscalar; // Create a logger. - auto logger = [&log]() -> const Logger& { return log; }; + // @TODO: Turn back on, once detray can use the ACTS logger + // auto logger = [&log]() -> const Logger& { return log; }; /***************************************************************** * Measurement Operations @@ -126,8 +128,10 @@ combinatorial_kalman_filter( std::vector> tips; // Create propagator + auto prop_cfg{config.propagation}; + prop_cfg.navigation.estimate_scattering_noise = false; traccc::details::ckf_propagator_t propagator( - config.propagation); + prop_cfg); // Create the input seeds container. bound_track_parameters_collection_types::const_device seeds{seeds_view}; @@ -142,12 +146,13 @@ combinatorial_kalman_filter( for (unsigned int step = 0u; step < config.max_track_candidates_per_track; step++) { - TRACCC_VERBOSE("Starting step " - << step + 1 << " / " - << config.max_track_candidates_per_track); + TRACCC_VERBOSE_HOST("Starting step " + << step + 1 << " / " + << config.max_track_candidates_per_track); // Iterate over input parameters const std::size_t n_in_params = in_params.size(); + TRACCC_VERBOSE_HOST("No. Params: " << n_in_params); // Terminate if there is no parameter to proceed if (n_in_params == 0) { @@ -197,10 +202,13 @@ combinatorial_kalman_filter( : links[step - 1][param_to_link[step - 1][in_param_id]] .ndf_sum); - TRACCC_VERBOSE("Processing input parameter " - << in_param_id + 1 << " / " << n_in_params << ": " - << in_param << " (orig_param_id=" << orig_param_id - << ", skip_counter=" << skip_counter << ")"); + TRACCC_VERBOSE_HOST("Processing input parameter " + << in_param_id + 1 << " / " << n_in_params); + TRACCC_VERBOSE_HOST("-> orig_param_id=" + << orig_param_id << ", skip_counter=" + << skip_counter << "\nVec:\n" + << in_param.vector()); + TRACCC_DEBUG_HOST("Cov:\n" << in_param.covariance()); /************************* * Material interaction @@ -209,8 +217,8 @@ combinatorial_kalman_filter( // Get surface corresponding to bound params const detray::tracking_surface sf{det, in_param.surface_link()}; - TRACCC_VERBOSE( - " free params: " << sf.bound_to_free_vector({}, in_param)); + TRACCC_VERBOSE_HOST(" Free params:\n" + << sf.bound_to_free_vector({}, in_param)); // Apply interactor if (sf.has_material()) { @@ -239,7 +247,9 @@ combinatorial_kalman_filter( best_links; // Iterate over the measurements + TRACCC_VERBOSE_HOST("No. measurements: " << (up - lo)); for (unsigned int meas_id = lo; meas_id < up; meas_id++) { + TRACCC_VERBOSE_HOST("Testing measurement: " << meas_id); // The measurement on surface to handle. const measurement& meas = measurements.at(meas_id); @@ -248,7 +258,7 @@ combinatorial_kalman_filter( auto trk_state = edm::make_track_state(measurements, meas_id); - const bool is_line = sf.template visit_mask(); + const bool is_line = detail::is_line(sf); // Run the Kalman update on a copy of the track parameters const kalman_fitter_status res = @@ -257,9 +267,13 @@ combinatorial_kalman_filter( const traccc::scalar chi2 = trk_state.filtered_chi2(); + TRACCC_DEBUG_HOST("KF status: " << fitter_debug_msg{res}()); + // The chi2 from Kalman update should be less than chi2_max if (res == kalman_fitter_status::SUCCESS && - chi2 < config.chi2_max) { + (chi2 < config.chi2_max)) { + + TRACCC_VERBOSE_HOST("Found measurement: " << meas_id); best_links.push_back( {{.step = step, @@ -284,9 +298,9 @@ combinatorial_kalman_filter( const unsigned int n_branches = std::min(config.max_num_branches_per_surface, static_cast(best_links.size())); - TRACCC_VERBOSE("Found " << n_branches << " branches for step " - << step << " and input parameter " - << in_param_id); + TRACCC_VERBOSE_HOST("Found " << n_branches << " branches for step " + << step << " and input parameter " + << in_param_id + 1); for (unsigned int i = 0; i < n_branches; ++i) { const auto& [link, filtered_params] = best_links[i]; @@ -295,9 +309,9 @@ combinatorial_kalman_filter( // Add the updated parameter to the updated parameters updated_params.push_back(filtered_params); - TRACCC_VERBOSE("updated_params[" - << updated_params.size() - 1 - << "] = " << updated_params.back()); + TRACCC_DEBUG_HOST("updated_params[" + << updated_params.size() - 1 + << "] = " << updated_params.back()); } /***************************************************************** @@ -318,10 +332,12 @@ combinatorial_kalman_filter( .chi2_sum = prev_chi2_sum, .ndf_sum = prev_ndf_sum}); + TRACCC_VERBOSE_HOST("Hole state created"); + updated_params.push_back(in_param); - TRACCC_VERBOSE("updated_params[" - << updated_params.size() - 1 - << "] = " << updated_params.back()); + TRACCC_DEBUG_HOST("updated_params[" + << updated_params.size() - 1 + << "] = " << updated_params.back()); } } @@ -438,6 +454,8 @@ combinatorial_kalman_filter( } if (this_is_dominated) { + TRACCC_VERBOSE_HOST( + "Track is dead (deduplication)!"); param_liveness.at(tracks.at(i)) = 0u; break; } @@ -465,11 +483,15 @@ combinatorial_kalman_filter( // link to be a tip if (links.at(step).at(link_id).n_skipped > config.max_num_skipping_per_cand) { + TRACCC_WARNING_HOST( + "Create tip: Max no. of holes reached! Bound param:\n" + << updated_params[link_id].vector()); tips.push_back({step, link_id}); continue; } const auto& param = updated_params[link_id]; + // Create propagator state typename traccc::details::ckf_propagator_t< detector_t, bfield_t>::state propagation(param, field, det); @@ -484,10 +506,11 @@ combinatorial_kalman_filter( traccc::details::ckf_interactor_t::state s2; typename interaction_register< traccc::details::ckf_interactor_t>::state s1{s2}; - // typename detray::parameter_resetter< - // typename detector_t::algebra_type>::state s3{}; + typename detray::parameter_resetter< + typename detector_t::algebra_type>::state s3{prop_cfg}; typename detray::momentum_aborter::state s4{}; typename ckf_aborter::state s5; + // Update the actor config s4.min_pT(static_cast(config.min_pT)); s4.min_p(static_cast(config.min_p)); @@ -495,21 +518,53 @@ combinatorial_kalman_filter( s5.max_count = config.max_step_counts_for_next_surface; // Propagate to the next surface - propagator.propagate(propagation, detray::tie(s0, s1, s2, s4, s5)); + TRACCC_DEBUG_HOST("Propagating... "); + propagator.propagate(propagation, + detray::tie(s0, s1, s2, s3, s4, s5)); + TRACCC_DEBUG_HOST("Finished propagation: On surface " + << propagation._navigation.barcode()); // If a surface found, add the parameter for the next // step - if (s5.success) { + bool valid_track{s5.success}; + if (valid_track) { assert(propagation._navigation.is_on_sensitive()); assert(!propagation._stepping.bound_params().is_invalid()); - out_params.push_back(propagation._stepping.bound_params()); - param_to_link[step].push_back(link_id); + const auto& out_param = propagation._stepping.bound_params(); + + const scalar theta = out_param.theta(); + if (theta <= 0.f || + theta >= 2.f * constant::pi) { + TRACCC_ERROR_HOST("Theta is hit pole after propagation"); + valid_track = false; + } + + if (!std::isfinite(out_param.phi())) { + TRACCC_ERROR_HOST( + "Phi is infinite after propagation (Matrix inversion)"); + valid_track = false; + } + + if (math::fabs(out_param.qop()) == 0.f) { + TRACCC_ERROR_HOST("q over p is zero after propagation"); + valid_track = false; + } + + if (valid_track) { + out_params.push_back(out_param); + param_to_link[step].push_back(link_id); + } } // Unless the track found a surface, it is considered a // tip - else if (!s5.success && - (step >= (config.min_track_candidates_per_track - 1u))) { + if (!valid_track && + (step >= (config.min_track_candidates_per_track - 1u))) { + if (!s5.success) { + TRACCC_VERBOSE_HOST("Create tip: No next sensitive found"); + } else { + TRACCC_VERBOSE_HOST("Create tip: Encountered error"); + } tips.push_back({step, link_id}); } @@ -517,6 +572,7 @@ combinatorial_kalman_filter( // kept as a tip if (s5.success && (step == (config.max_track_candidates_per_track - 1u))) { + TRACCC_ERROR_HOST("Create tip: Max no. candidates"); tips.push_back({step, link_id}); } } diff --git a/core/include/traccc/finding/propagation_data.hpp b/core/include/traccc/finding/propagation_data.hpp new file mode 100644 index 0000000000..e22a1c4c01 --- /dev/null +++ b/core/include/traccc/finding/propagation_data.hpp @@ -0,0 +1,25 @@ +/** TRACCC library, part of the ACTS project (R&D line) + * + * (c) 2025 CERN for the benefit of the ACTS project + * + * Mozilla Public License Version 2.0 + */ + +#pragma once + +// Detray include(s). +#include + +namespace traccc { + +/// Data from the propagation loop that has to be kept between CKF steps +struct propagation_data { + /// The surface that was visited before the current one (overlaps + /// resolution) + detray::dindex prev_surface{detray::dindex_invalid}; + + /// Is the current surface hit in the extended tolerance band? + bool is_edge{false}; +}; + +} // namespace traccc diff --git a/core/include/traccc/fitting/details/kalman_fitting.hpp b/core/include/traccc/fitting/details/kalman_fitting.hpp index 2c0eab17bb..435f5b7ab1 100644 --- a/core/include/traccc/fitting/details/kalman_fitting.hpp +++ b/core/include/traccc/fitting/details/kalman_fitting.hpp @@ -88,7 +88,7 @@ typename edm::track_container::host kalman_fitting( result_tracks_device.at(result_tracks_device.size() - 1), typename edm::track_state_collection::device{ vecmem::get_data(result.states)}, - measurements, seqs_buffer); + measurements, seqs_buffer, fitter.config().propagation); // Run the fitter. The status that it returns is not used here. The main // failure modes are saved onto the fitted track itself. Not sure what diff --git a/core/include/traccc/fitting/fitting_config.hpp b/core/include/traccc/fitting/fitting_config.hpp index 9d956203a9..1e30ca0298 100644 --- a/core/include/traccc/fitting/fitting_config.hpp +++ b/core/include/traccc/fitting/fitting_config.hpp @@ -34,7 +34,7 @@ struct fitting_config { std::size_t barcode_sequence_size_factor = 5; std::size_t min_barcode_sequence_capacity = 100; traccc::scalar backward_filter_mask_tolerance = - 5.f * traccc::unit::mm; + 10.f * traccc::unit::m; }; } // namespace traccc diff --git a/core/include/traccc/fitting/kalman_filter/gain_matrix_updater.hpp b/core/include/traccc/fitting/kalman_filter/gain_matrix_updater.hpp index c96b2d4f76..03577636cc 100644 --- a/core/include/traccc/fitting/kalman_filter/gain_matrix_updater.hpp +++ b/core/include/traccc/fitting/kalman_filter/gain_matrix_updater.hpp @@ -14,6 +14,7 @@ #include "traccc/edm/measurement_helpers.hpp" #include "traccc/edm/track_state_collection.hpp" #include "traccc/fitting/status_codes.hpp" +#include "traccc/utils/logging.hpp" namespace traccc { @@ -62,6 +63,9 @@ struct gain_matrix_updater { static constexpr unsigned int D = 2; + TRACCC_VERBOSE_HOST_DEVICE("In gain-matrix-updater..."); + TRACCC_VERBOSE_HOST_DEVICE("Measurement dim: %d", dim); + assert(dim == 1u || dim == 2u); assert(!bound_params.is_invalid()); @@ -79,6 +83,8 @@ struct gain_matrix_updater { assert((dim > 1) || (getter::element(meas_local, 1u, 0u) == 0.f)); + TRACCC_DEBUG_HOST("Predicted param.: " << bound_params); + // Predicted vector of bound track parameters const bound_vector_type& predicted_vec = bound_params.vector(); @@ -109,6 +115,11 @@ struct gain_matrix_updater { getter::element(V, 1u, 1u) = 1.f; } + TRACCC_DEBUG_HOST("Measurement position: " << meas_local); + TRACCC_DEBUG_HOST("Measurement variance:\n" << V); + TRACCC_DEBUG_HOST("Predicted residual: " << meas_local - + H * predicted_vec); + const matrix_type projected_cov = algebra::matrix::transposed_product(predicted_cov, H); @@ -118,11 +129,49 @@ struct gain_matrix_updater { assert(matrix::determinant(M) != 0.f); const matrix_type<6, D> K = projected_cov * matrix::inverse(M); + TRACCC_DEBUG_HOST("H:\n" << H); + TRACCC_DEBUG_HOST("K:\n" << K); + // Calculate the filtered track parameters const matrix_type<6, 1> filtered_vec = predicted_vec + K * (meas_local - H * predicted_vec); const matrix_type<6, 6> filtered_cov = (I66 - K * H) * predicted_cov; + TRACCC_DEBUG_HOST("Filtered param:\n" << filtered_vec); + TRACCC_DEBUG_HOST("Filtered cov:\n" << filtered_cov); + + // Check the covariance for consistency + if (getter::element(filtered_cov, 0, 0) <= 0.f || + getter::element(filtered_cov, 1, 1) <= 0.f || + getter::element(filtered_cov, 2, 2) <= 0.f || + getter::element(filtered_cov, 3, 3) <= 0.f || + getter::element(filtered_cov, 4, 4) <= 0.f || + getter::element(filtered_cov, 5, 5) <= 0.f) { + TRACCC_ERROR_HOST_DEVICE( + "Invalid diagonal element in covariance matrix after " + "filtering"); + return kalman_fitter_status::ERROR_UPDATER_INVALID_COVARIANCE; + } + + // Return false if track is parallel to z-axis or phi is not finite + if (!std::isfinite(getter::element(filtered_vec, e_bound_theta, 0))) { + TRACCC_ERROR_HOST_DEVICE( + "Theta is infinite after filtering (Matrix inversion)"); + return kalman_fitter_status::ERROR_INVERSION; + } + + if (!std::isfinite(getter::element(filtered_vec, e_bound_phi, 0))) { + TRACCC_ERROR_HOST_DEVICE( + "Phi is infinite after filtering (Matrix inversion)"); + return kalman_fitter_status::ERROR_INVERSION; + } + + if (math::fabs(getter::element(filtered_vec, e_bound_qoverp, 0)) == + 0.f) { + TRACCC_ERROR_HOST_DEVICE("q/p is zero after filtering"); + return kalman_fitter_status::ERROR_QOP_ZERO; + } + // Residual between measurement and (projected) filtered vector const matrix_type residual = meas_local - H * filtered_vec; @@ -133,36 +182,38 @@ struct gain_matrix_updater { residual, matrix::inverse(R)) * residual; - // Return false if track is parallel to z-axis or phi is not finite - const scalar theta = bound_params.theta(); + const scalar chi2_val{getter::element(chi2, 0, 0)}; - if (theta <= 0.f || theta >= constant::pi) { - return kalman_fitter_status::ERROR_THETA_ZERO; - } + TRACCC_VERBOSE_HOST("Filtered residual: " << residual); + TRACCC_DEBUG_HOST("R:\n" << R); + TRACCC_DEBUG_HOST_DEVICE("det(R): %f", matrix::determinant(R)); + TRACCC_DEBUG_HOST("R_inv:\n" << matrix::inverse(R)); + TRACCC_VERBOSE_HOST_DEVICE("Chi2: %f", chi2_val); - if (!std::isfinite(bound_params.phi())) { - return kalman_fitter_status::ERROR_INVERSION; - } - - if (std::abs(bound_params.qop()) == 0.f) { - return kalman_fitter_status::ERROR_QOP_ZERO; - } - - if (getter::element(chi2, 0, 0) < 0.f) { + if (chi2_val < 0.f) { + TRACCC_ERROR_HOST_DEVICE("Chi2 negative"); return kalman_fitter_status::ERROR_UPDATER_CHI2_NEGATIVE; } - if (!std::isfinite(getter::element(chi2, 0, 0))) { + if (!std::isfinite(chi2_val)) { + TRACCC_ERROR_HOST_DEVICE("Chi2 infinite"); return kalman_fitter_status::ERROR_UPDATER_CHI2_NOT_FINITE; } - // Set the track state parameters + // Set the chi2 for this track and measurement trk_state.filtered_params().set_vector(filtered_vec); trk_state.filtered_params().set_covariance(filtered_cov); - trk_state.filtered_chi2() = getter::element(chi2, 0, 0); + trk_state.filtered_chi2() = chi2_val; + + // Wrap the phi and theta angles in their valid ranges + normalize_angles(trk_state.filtered_params()); - // Wrap the phi in the range of [-pi, pi] - wrap_phi(trk_state.filtered_params()); + const scalar theta = trk_state.filtered_params().theta(); + if (theta <= 0.f || theta >= 2.f * constant::pi) { + TRACCC_ERROR_HOST_DEVICE("Hit theta pole after filtering : %f", + theta); + return kalman_fitter_status::ERROR_THETA_POLE; + } assert(!trk_state.filtered_params().is_invalid()); diff --git a/core/include/traccc/fitting/kalman_filter/is_line_visitor.hpp b/core/include/traccc/fitting/kalman_filter/is_line_visitor.hpp index 34c3d5a14a..a4846c5774 100644 --- a/core/include/traccc/fitting/kalman_filter/is_line_visitor.hpp +++ b/core/include/traccc/fitting/kalman_filter/is_line_visitor.hpp @@ -7,20 +7,31 @@ #pragma once +#include #include +#include #include "traccc/definitions/qualifiers.hpp" -namespace traccc { +namespace traccc::detail { -struct is_line_visitor { - template - [[nodiscard]] TRACCC_HOST_DEVICE inline bool operator()( - const mask_group_t& /*mask_group*/, const index_t& /*index*/) const { - using shape_type = typename mask_group_t::value_type::shape; - return std::is_same_v> || - std::is_same_v>; +/// @returns true if the surface has "line" shape +template +[[nodiscard]] TRACCC_HOST_DEVICE bool constexpr is_line( + const detray::geometry::surface sf) { + using algebra_t = typename detector_t::algebra_type; + using straw_tube = detray::mask, algebra_t>; + using wire_cell = detray::mask, algebra_t>; + + if constexpr (detector_t::masks::template is_defined() || + detector_t::masks::template is_defined()) { + return (sf.shape_id() == + detector_t::masks::template get_id()) || + (sf.shape_id() == + detector_t::masks::template get_id()); + } else { + return false; } }; -} // namespace traccc +} // namespace traccc::detail diff --git a/core/include/traccc/fitting/kalman_filter/kalman_actor.hpp b/core/include/traccc/fitting/kalman_filter/kalman_actor.hpp index 94b70c1f91..61b5447363 100644 --- a/core/include/traccc/fitting/kalman_filter/kalman_actor.hpp +++ b/core/include/traccc/fitting/kalman_filter/kalman_actor.hpp @@ -16,6 +16,7 @@ #include "traccc/fitting/kalman_filter/is_line_visitor.hpp" #include "traccc/fitting/kalman_filter/two_filters_smoother.hpp" #include "traccc/fitting/status_codes.hpp" +#include "traccc/utils/logging.hpp" #include "traccc/utils/particle.hpp" // detray include(s). @@ -142,9 +143,17 @@ struct kalman_actor : detray::actor { // triggered only for sensitive surfaces if (navigation.is_on_sensitive()) { + TRACCC_VERBOSE_HOST_DEVICE("In Kalman actor..."); + TRACCC_DEBUG_HOST("-> on surface: " << navigation.get_surface()); + typename edm::track_state_collection::device::proxy_type trk_state = actor_state(); + TRACCC_DEBUG_HOST( + "-> expecting: " << actor_state.m_measurements + .at(trk_state.measurement_index()) + .surface_link); + // Increase the hole counts if the propagator fails to find the next // measurement if (navigation.barcode() != @@ -156,6 +165,8 @@ struct kalman_actor : detray::actor { return; } + TRACCC_VERBOSE_HOST_DEVICE("Found next track state to fit"); + // This track state is not a hole if (!actor_state.backward_mode) { trk_state.set_hole(false); @@ -164,7 +175,7 @@ struct kalman_actor : detray::actor { // Run Kalman Gain Updater const auto sf = navigation.get_surface(); - const bool is_line = sf.template visit_mask(); + const bool is_line = detail::is_line(sf); kalman_fitter_status res = kalman_fitter_status::SUCCESS; @@ -173,7 +184,11 @@ struct kalman_actor : detray::actor { kalman_actor_direction::FORWARD_ONLY || direction_e == kalman_actor_direction::BIDIRECTIONAL) { + // Wrap the phi and theta angles in their valid ranges + normalize_angles(propagation._stepping.bound_params()); + // Forward filter + TRACCC_DEBUG_HOST_DEVICE("Run filtering"); res = gain_matrix_updater{}( trk_state, actor_state.m_measurements, propagation._stepping.bound_params(), is_line); @@ -189,6 +204,7 @@ struct kalman_actor : detray::actor { direction_e == kalman_actor_direction::BIDIRECTIONAL) { // Backward filter for smoothing + TRACCC_DEBUG_HOST_DEVICE("Run smoothing"); res = two_filters_smoother{}( trk_state, actor_state.m_measurements, propagation._stepping.bound_params(), is_line); @@ -199,6 +215,16 @@ struct kalman_actor : detray::actor { // Abort if the Kalman update fails if (res != kalman_fitter_status::SUCCESS) { + if (actor_state.backward_mode) { + TRACCC_ERROR_DEVICE("Abort backward fit: KF status %d", + res); + TRACCC_ERROR_HOST( + "Abort backward fit: " << fitter_debug_msg{res}()); + } else { + TRACCC_ERROR_DEVICE("Abort forward fit: KF status %d", res); + TRACCC_ERROR_HOST( + "Abort forward fit: " << fitter_debug_msg{res}()); + } propagation._heartbeat &= navigation.abort(fitter_debug_msg{res}); return; @@ -214,8 +240,13 @@ struct kalman_actor : detray::actor { // Update iterator actor_state.next(); - // Flag renavigation of the current candidate - navigation.set_high_trust(); + // Flag renavigation of the current candidate (unless for overlap) + if (math::fabs(navigation()) > 1.f * unit::um) { + navigation.set_high_trust(); + } else { + TRACCC_DEBUG_HOST_DEVICE( + "Encountered overlap, jump to next surface"); + } } } }; diff --git a/core/include/traccc/fitting/kalman_filter/kalman_fitter.hpp b/core/include/traccc/fitting/kalman_filter/kalman_fitter.hpp index 06eddb4d85..ba5959084f 100644 --- a/core/include/traccc/fitting/kalman_filter/kalman_fitter.hpp +++ b/core/include/traccc/fitting/kalman_filter/kalman_fitter.hpp @@ -20,6 +20,7 @@ #include "traccc/fitting/kalman_filter/statistics_updater.hpp" #include "traccc/fitting/kalman_filter/two_filters_smoother.hpp" #include "traccc/fitting/status_codes.hpp" +#include "traccc/utils/logging.hpp" #include "traccc/utils/particle.hpp" #include "traccc/utils/prob.hpp" #include "traccc/utils/propagation.hpp" @@ -55,6 +56,7 @@ class kalman_fitter { // Actor types using aborter = detray::pathlimit_aborter; + using momentum_aborter = detray::momentum_aborter; using transporter = detray::parameter_transporter; using interactor = detray::pointwise_material_interactor; using forward_fit_actor = @@ -71,11 +73,13 @@ class kalman_fitter { using forward_actor_chain_type = detray::actor_chain; + resetter, momentum_aborter, barcode_sequencer, + kalman_step_aborter>; using backward_actor_chain_type = detray::actor_chain; + interactor, resetter, momentum_aborter, + kalman_step_aborter>; // Navigator type for backward propagator using direct_navigator_type = detray::direct_navigator; @@ -110,11 +114,13 @@ class kalman_fitter { track_states, const measurement_collection_types::const_device& measurements, vecmem::data::vector_view - sequence_buffer) + sequence_buffer, + const detray::propagation::config& prop_cfg) : m_fit_actor_state{track, track_states, measurements}, m_sequencer_state( vecmem::device_vector( sequence_buffer)), + m_parameter_resetter{prop_cfg}, m_fit_res{track}, m_sequence_buffer(sequence_buffer) {} @@ -122,7 +128,8 @@ class kalman_fitter { TRACCC_HOST_DEVICE typename forward_actor_chain_type::state_ref_tuple operator()() { return detray::tie(m_aborter_state, m_interactor_state, - m_fit_actor_state, m_sequencer_state, + m_fit_actor_state, m_parameter_resetter, + m_montum_aborter_state, m_sequencer_state, m_step_aborter_state); } @@ -131,7 +138,8 @@ class kalman_fitter { typename backward_actor_chain_type::state_ref_tuple backward_actor_state() { return detray::tie(m_aborter_state, m_fit_actor_state, - m_interactor_state, m_step_aborter_state); + m_interactor_state, m_parameter_resetter, + m_montum_aborter_state, m_step_aborter_state); } /// Individual actor states @@ -140,6 +148,8 @@ class kalman_fitter { typename forward_fit_actor::state m_fit_actor_state; typename barcode_sequencer::state m_sequencer_state; kalman_step_aborter::state m_step_aborter_state{}; + typename resetter::state m_parameter_resetter{}; + typename momentum_aborter::state m_montum_aborter_state{}; /// Fitting result per track typename edm::track_collection::device::proxy_type @@ -257,6 +267,7 @@ class kalman_fitter { smooth(state& fitter_state) const { if (fitter_state.m_sequencer_state.overflow) { + TRACCC_ERROR_HOST_DEVICE("Barcode sequence overlow"); return kalman_fitter_status::ERROR_BARCODE_SEQUENCE_OVERFLOW; } @@ -277,16 +288,26 @@ class kalman_fitter { } auto last = fitter_state.m_fit_actor_state(); - const scalar theta = last.filtered_params().theta(); - if (theta <= 0.f || theta >= constant::pi) { - return kalman_fitter_status::ERROR_THETA_ZERO; + + if (!std::isfinite(last.filtered_params().theta())) { + TRACCC_ERROR_HOST_DEVICE( + "Theta is infinite after filtering (Matrix inversion)"); + return kalman_fitter_status::ERROR_INVERSION; } if (!std::isfinite(last.filtered_params().phi())) { + TRACCC_ERROR_HOST_DEVICE( + "Phi is infinite after filtering (Matrix inversion)"); return kalman_fitter_status::ERROR_INVERSION; } + if (theta <= 0.f || theta >= 2.f * constant::pi) { + TRACCC_ERROR_HOST_DEVICE("Hit theta pole after filtering : %f", + theta); + return kalman_fitter_status::ERROR_THETA_POLE; + } + last.smoothed_params().set_parameter_vector(last.filtered_params()); last.smoothed_params().set_covariance( last.filtered_params().covariance()); @@ -324,6 +345,8 @@ class kalman_fitter { propagation._navigation.set_direction( detray::navigation::direction::e_backward); + propagation._navigation.safe_step_size = + 0.1f * traccc::unit::mm; // Synchronize the current barcode with the input track parameter while (propagation._navigation.get_target_barcode() != @@ -393,6 +416,7 @@ class kalman_fitter { auto trk_state = track_states.at(link.index); // Fitting fails if any of non-hole track states is not smoothed if (!trk_state.is_hole() && !trk_state.is_smoothed()) { + TRACCC_ERROR_HOST_DEVICE("Not all smoothed"); fit_res.fit_outcome() = track_fit_outcome::FAILURE_NOT_ALL_SMOOTHED; return; @@ -404,6 +428,7 @@ class kalman_fitter { return; } + TRACCC_ERROR_HOST_DEVICE("Negative NDF"); fit_res.fit_outcome() = track_fit_outcome::FAILURE_NON_POSITIVE_NDF; 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 index 6378f00d74..93aec62fe1 100644 --- a/core/include/traccc/fitting/kalman_filter/kalman_step_aborter.hpp +++ b/core/include/traccc/fitting/kalman_filter/kalman_step_aborter.hpp @@ -9,6 +9,7 @@ // Project include(s). #include "traccc/definitions/qualifiers.hpp" +#include "traccc/utils/logging.hpp" // Detray include(s). #include @@ -55,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"); + TRACCC_ERROR_HOST_DEVICE( + "Maximum number of steps to reach next " + "sensitive surface exceeded"); prop_state._heartbeat &= navigation.abort( "Kalman Fitter: Maximum number of steps to reach next " "sensitive surface exceeded"); diff --git a/core/include/traccc/fitting/kalman_filter/two_filters_smoother.hpp b/core/include/traccc/fitting/kalman_filter/two_filters_smoother.hpp index 81dcbe44d8..53f327ea7a 100644 --- a/core/include/traccc/fitting/kalman_filter/two_filters_smoother.hpp +++ b/core/include/traccc/fitting/kalman_filter/two_filters_smoother.hpp @@ -14,6 +14,7 @@ #include "traccc/edm/measurement_helpers.hpp" #include "traccc/edm/track_state_collection.hpp" #include "traccc/fitting/status_codes.hpp" +#include "traccc/utils/logging.hpp" namespace traccc { @@ -67,6 +68,8 @@ struct two_filters_smoother { // Do not smoothe if the forward pass produced an error if (trk_state.filtered_params().is_invalid()) { + TRACCC_ERROR_HOST_DEVICE("Filtered track state invalid"); + TRACCC_ERROR_HOST(trk_state.filtered_params()); return kalman_fitter_status::ERROR_INVALID_TRACK_STATE; } @@ -99,6 +102,19 @@ struct two_filters_smoother { const matrix_type smoothed_cov = matrix::inverse(smoothed_cov_inv); + // Check the covariance for consistency + if (getter::element(smoothed_cov, 0, 0) <= 0.f || + getter::element(smoothed_cov, 1, 1) <= 0.f || + getter::element(smoothed_cov, 2, 2) <= 0.f || + getter::element(smoothed_cov, 3, 3) <= 0.f || + getter::element(smoothed_cov, 4, 4) <= 0.f || + getter::element(smoothed_cov, 5, 5) <= 0.f) { + TRACCC_ERROR_HOST_DEVICE( + "Invalid diagonal element in covariance matrix after " + "filtering"); + return kalman_fitter_status::ERROR_SMOOTHER_INVALID_COVARIANCE; + } + // Eq (3.38) of "Pattern Recognition, Tracking and Vertex // Reconstruction in Particle Detectors" const matrix_type smoothed_vec = @@ -107,8 +123,30 @@ struct two_filters_smoother { predicted_cov_inv * predicted_vec); trk_state.smoothed_params().set_vector(smoothed_vec); + + // Return false if track is parallel to z-axis or phi is not finite + if (!std::isfinite(trk_state.smoothed_params().theta())) { + TRACCC_ERROR_HOST_DEVICE( + "Theta is infinite after smoothing (Matrix inversion)"); + return kalman_fitter_status::ERROR_INVERSION; + } + + if (!std::isfinite(trk_state.smoothed_params().phi())) { + TRACCC_ERROR_HOST_DEVICE( + "Phi is infinite after smoothing (Matrix inversion)"); + return kalman_fitter_status::ERROR_INVERSION; + } + + if (math::fabs(trk_state.smoothed_params().qop()) == 0.f) { + TRACCC_ERROR_HOST_DEVICE("q/p is zero after smoothing"); + return kalman_fitter_status::ERROR_QOP_ZERO; + } + trk_state.smoothed_params().set_covariance(smoothed_cov); + // Wrap the phi and theta angles in their valid ranges + normalize_angles(trk_state.smoothed_params()); + matrix_type H = measurements.at(trk_state.measurement_index()) .subs.template projector(); @@ -127,6 +165,11 @@ struct two_filters_smoother { getter::element(V, 1u, 1u) = 1.f; } + TRACCC_DEBUG_HOST("Measurement position: " << meas_local); + TRACCC_DEBUG_HOST("Measurement variance:\n" << V); + TRACCC_DEBUG_HOST("Predicted residual: " << meas_local - + H * predicted_vec); + // Eq (3.39) of "Pattern Recognition, Tracking and Vertex // Reconstruction in Particle Detectors" const matrix_type R_smt = @@ -141,11 +184,22 @@ struct two_filters_smoother { residual_smt, matrix::inverse(R_smt)) * residual_smt; - if (getter::element(chi2_smt, 0, 0) < 0.f) { + const scalar chi2_smt_value{getter::element(chi2_smt, 0, 0)}; + + TRACCC_VERBOSE_HOST("Smoothed residual: " << residual_smt); + TRACCC_DEBUG_HOST("R_smt:\n" << R_smt); + TRACCC_DEBUG_HOST_DEVICE("det(R_smt): %f", matrix::determinant(R_smt)); + TRACCC_DEBUG_HOST("R_smt_inv:\n" << matrix::inverse(R_smt)); + TRACCC_VERBOSE_HOST_DEVICE("Chi2: %f", chi2_smt_value); + + if (chi2_smt_value < 0.f) { + TRACCC_ERROR_HOST_DEVICE("Smoothed chi2 negative: %f", + chi2_smt_value); return kalman_fitter_status::ERROR_SMOOTHER_CHI2_NEGATIVE; } - if (!std::isfinite(getter::element(chi2_smt, 0, 0))) { + if (!std::isfinite(chi2_smt_value)) { + TRACCC_ERROR_HOST_DEVICE("Smoothed chi2 infinite"); return kalman_fitter_status::ERROR_SMOOTHER_CHI2_NOT_FINITE; } @@ -175,54 +229,94 @@ struct two_filters_smoother { assert(std::isfinite(matrix::determinant(M))); const matrix_type<6, D> K = projected_cov * matrix::inverse(M); + TRACCC_DEBUG_HOST("H:\n" << H); + TRACCC_DEBUG_HOST("K:\n" << K); + // Calculate the filtered track parameters const matrix_type<6, 1> filtered_vec = predicted_vec + K * (meas_local - H * predicted_vec); const matrix_type<6, 6> filtered_cov = (I66 - K * H) * predicted_cov; - // Residual between measurement and (projected) filtered vector - const matrix_type residual = meas_local - H * filtered_vec; - - // Calculate backward chi2 - const matrix_type R = (I_m - H * K) * V; - // assert(matrix::determinant(R) != 0.f); - assert(std::isfinite(matrix::determinant(R))); - const matrix_type<1, 1> chi2 = - algebra::matrix::transposed_product( - residual, matrix::inverse(R)) * - residual; + // Check the covariance for consistency + if (getter::element(filtered_cov, 0, 0) <= 0.f || + getter::element(filtered_cov, 1, 1) <= 0.f || + getter::element(filtered_cov, 2, 2) <= 0.f || + getter::element(filtered_cov, 3, 3) <= 0.f || + getter::element(filtered_cov, 4, 4) <= 0.f || + getter::element(filtered_cov, 5, 5) <= 0.f) { + TRACCC_ERROR_HOST_DEVICE( + "Invalid diagonal element in covariance matrix after " + "smoothing"); + return kalman_fitter_status::ERROR_SMOOTHER_INVALID_COVARIANCE; + } // Update the bound track parameters bound_params.set_vector(filtered_vec); - bound_params.set_covariance(filtered_cov); // Return false if track is parallel to z-axis or phi is not finite - const scalar theta = bound_params.theta(); - if (theta <= 0.f || theta >= constant::pi) { - return kalman_fitter_status::ERROR_THETA_ZERO; + if (!std::isfinite(bound_params.theta())) { + TRACCC_ERROR_HOST_DEVICE( + "Theta is infinite after filering in smoother (Matrix " + "inversion)"); + return kalman_fitter_status::ERROR_INVERSION; } if (!std::isfinite(bound_params.phi())) { + TRACCC_ERROR_HOST_DEVICE( + "Phi is infinite after filering in smoother (Matrix " + "inversion)"); return kalman_fitter_status::ERROR_INVERSION; } - if (std::abs(bound_params.qop()) == 0.f) { + if (math::fabs(bound_params.qop()) == 0.f) { + TRACCC_ERROR_HOST_DEVICE("q/p is zero after filering in smoother"); return kalman_fitter_status::ERROR_QOP_ZERO; } - if (getter::element(chi2, 0, 0) < 0.f) { - return kalman_fitter_status::ERROR_UPDATER_CHI2_NEGATIVE; + bound_params.set_covariance(filtered_cov); + + // Residual between measurement and (projected) filtered vector + const matrix_type residual = meas_local - H * filtered_vec; + + // Calculate backward chi2 + const matrix_type R = (I_m - H * K) * V; + // assert(matrix::determinant(R) != 0.f); // @TODO: This fails + assert(std::isfinite(matrix::determinant(R))); + const matrix_type<1, 1> chi2 = + algebra::matrix::transposed_product( + residual, matrix::inverse(R)) * + residual; + + const scalar chi2_val{getter::element(chi2, 0, 0)}; + + TRACCC_VERBOSE_HOST("Filtered residual: " << residual); + TRACCC_DEBUG_HOST("R:\n" << R); + TRACCC_DEBUG_HOST_DEVICE("det(R): %f", matrix::determinant(R)); + TRACCC_DEBUG_HOST("R_inv:\n" << matrix::inverse(R)); + TRACCC_VERBOSE_HOST_DEVICE("Chi2: %f", chi2_val); + + if (chi2_val < 0.f) { + TRACCC_ERROR_HOST_DEVICE("Chi2 negative: %f", chi2_val); + return kalman_fitter_status::ERROR_SMOOTHER_CHI2_NEGATIVE; } - if (!std::isfinite(getter::element(chi2, 0, 0))) { - return kalman_fitter_status::ERROR_UPDATER_CHI2_NOT_FINITE; + if (!std::isfinite(chi2_val)) { + TRACCC_ERROR_HOST_DEVICE("Chi2 infinite"); + return kalman_fitter_status::ERROR_SMOOTHER_CHI2_NOT_FINITE; } // Set backward chi2 - trk_state.backward_chi2() = getter::element(chi2, 0, 0); + trk_state.backward_chi2() = chi2_val; - // Wrap the phi in the range of [-pi, pi] - wrap_phi(bound_params); + // Wrap the phi and theta angles in their valid ranges + normalize_angles(bound_params); + + const scalar theta = bound_params.theta(); + if (theta <= 0.f || theta >= 2.f * constant::pi) { + TRACCC_ERROR_HOST_DEVICE("Hit theta pole after smoothing : %f", + theta); + return kalman_fitter_status::ERROR_THETA_POLE; + } trk_state.set_smoothed(); diff --git a/core/include/traccc/fitting/status_codes.hpp b/core/include/traccc/fitting/status_codes.hpp index ae7ce138af..05962e5add 100644 --- a/core/include/traccc/fitting/status_codes.hpp +++ b/core/include/traccc/fitting/status_codes.hpp @@ -14,12 +14,14 @@ namespace traccc { enum class kalman_fitter_status : uint32_t { SUCCESS, ERROR_QOP_ZERO, - ERROR_THETA_ZERO, + ERROR_THETA_POLE, ERROR_INVERSION, ERROR_SMOOTHER_CHI2_NEGATIVE, ERROR_SMOOTHER_CHI2_NOT_FINITE, + ERROR_SMOOTHER_INVALID_COVARIANCE, ERROR_UPDATER_CHI2_NEGATIVE, ERROR_UPDATER_CHI2_NOT_FINITE, + ERROR_UPDATER_INVALID_COVARIANCE, ERROR_BARCODE_SEQUENCE_OVERFLOW, ERROR_INVALID_TRACK_STATE, ERROR_OTHER, @@ -36,8 +38,8 @@ struct fitter_debug_msg { case ERROR_QOP_ZERO: { return msg + "Track qop is zero"; } - case ERROR_THETA_ZERO: { - return msg + "Track theta is zero"; + case ERROR_THETA_POLE: { + return msg + "Track theta hit pole"; } case ERROR_INVERSION: { return msg + "Failed matrix inversion"; @@ -49,10 +51,10 @@ struct fitter_debug_msg { return msg + "Invalid chi2 in smoother"; } case ERROR_UPDATER_CHI2_NEGATIVE: { - return msg + "Negative chi2 in gain matrix update"; + return msg + "Negative chi2 in gain matrix updater"; } case ERROR_UPDATER_CHI2_NOT_FINITE: { - return msg + "Invalid chi2 in gain matrix update"; + return msg + "Invalid chi2 in gain matrix updater"; } case ERROR_BARCODE_SEQUENCE_OVERFLOW: { return msg + "Barcode sequence overflow in direct navigator"; @@ -61,6 +63,12 @@ struct fitter_debug_msg { return msg + "Invalid track state in forward pass (skipped or error)"; } + case ERROR_UPDATER_INVALID_COVARIANCE: { + return msg + "Invalid track covariance after filtering"; + } + case ERROR_SMOOTHER_INVALID_COVARIANCE: { + return msg + "Invalid track covariance after smoothing"; + } case ERROR_OTHER: { return msg + "Unspecified error"; } diff --git a/core/include/traccc/seeding/doublet_finding_helper.hpp b/core/include/traccc/seeding/doublet_finding_helper.hpp index 2135cb2323..089955ed92 100644 --- a/core/include/traccc/seeding/doublet_finding_helper.hpp +++ b/core/include/traccc/seeding/doublet_finding_helper.hpp @@ -77,7 +77,7 @@ bool TRACCC_HOST_DEVICE doublet_finding_helper::isCompatible( (math::fabs(cotTheta) >= config.cotThetaMax * deltaR) || (zOrigin <= config.collisionRegionMin * deltaR) || (zOrigin >= config.collisionRegionMax * deltaR) || - std::abs(cotTheta) >= config.deltaZMax) { + math::fabs(cotTheta) >= config.deltaZMax) { return false; } diff --git a/core/include/traccc/utils/logging.hpp b/core/include/traccc/utils/logging.hpp index 225dbda6aa..4ce2e4cf8c 100644 --- a/core/include/traccc/utils/logging.hpp +++ b/core/include/traccc/utils/logging.hpp @@ -10,6 +10,9 @@ // Acts include(s). #include +// detray include(s) +#include + namespace traccc { /// Use the @c Acts::Logging namespace. @@ -46,3 +49,80 @@ const Logger& getDummyLogger(); #define TRACCC_WARNING(x) ACTS_WARNING(x) #define TRACCC_ERROR(x) ACTS_ERROR(x) #define TRACCC_FATAL(x) ACTS_FATAL(x) + +/// @TODO: Temporary usage of detray logger, until the ACTS logger can be +/// used/passed in host-device code + +// Exclude SYCL and HIP builds, which run into linker errors +#if defined(CL_SYCL_LANGUAGE_VERSION) || defined(SYCL_LANGUAGE_VERSION) || \ + defined(__HIP__) +#define __EXCLUDE_LOGS__ +#endif + +#if !defined(__EXCLUDE_LOGS__) +// Host log messages in host-device code (@TODO: This will currently claim to +// come from "DETRAY") +#define TRACCC_VERBOSE_HOST(x) DETRAY_VERBOSE_HOST(x) +#define TRACCC_DEBUG_HOST(x) DETRAY_DEBUG_HOST(x) +#define TRACCC_INFO_HOST(x) DETRAY_INFO_HOST(x) +#define TRACCC_WARNING_HOST(x) DETRAY_WARN_HOST(x) +#define TRACCC_ERROR_HOST(x) DETRAY_ERROR_HOST(x) +#define TRACCC_FATAL_HOST(x) DETRAY_FATAL_HOST(x) + +// Host-device log messages for host-device code +#define TRACCC_WARNING_HOST_DEVICE(x, ...) \ + DETRAY_LOG_PRINTF("TRACCC", "WARNING", x, __VA_ARGS__) +#define TRACCC_ERROR_HOST_DEVICE(x, ...) \ + DETRAY_LOG_PRINTF("TRACCC", "ERROR", x, __VA_ARGS__) +#define TRACCC_FATAL_HOST_DEVICE(x, ...) \ + DETRAY_LOG_PRINTF("TRACCC", "FATAL", x, __VA_ARGS__) +#else +#define TRACCC_VERBOSE_HOST(x) +#define TRACCC_DEBUG_HOST(x) +#define TRACCC_INFO_HOST(x) +#define TRACCC_WARNING_HOST(x) +#define TRACCC_ERROR_HOST(x) +#define TRACCC_FATAL_HOST(x) + +#define TRACCC_WARNING_HOST_DEVICE(x, ...) +#define TRACCC_ERROR_HOST_DEVICE(x, ...) +#define TRACCC_FATAL_HOST_DEVICE(x, ...) +#endif // defined(__EXCLUDE_LOGS__) + +#if DETRAY_LOG_LVL > 0 && !defined(__EXCLUDE_LOGS__) +#define TRACCC_INFO_HOST_DEVICE(x, ...) \ + DETRAY_LOG_PRINTF("TRACCC", "INFO", x, __VA_ARGS__) +#else +#define TRACCC_INFO_HOST_DEVICE(x, ...) +#endif + +#if DETRAY_LOG_LVL > 1 && !defined(__EXCLUDE_LOGS__) +#define TRACCC_VERBOSE_HOST_DEVICE(x, ...) \ + DETRAY_LOG_PRINTF("TRACCC", "VERBOSE", x, __VA_ARGS__) +#else +#define TRACCC_VERBOSE_HOST_DEVICE(x, ...) +#endif + +#if DETRAY_LOG_LVL > 2 && !defined(__EXCLUDE_LOGS__) +#define TRACCC_DEBUG_HOST_DEVICE(x, ...) \ + DETRAY_LOG_PRINTF("TRACCC", "DEBUG", x, __VA_ARGS__) +#else +#define TRACCC_DEBUG_HOST_DEVICE(x, ...) +#endif + +// Device log messages in device/host-device code +#if defined(__DEVICE_LOGGING__) && !defined(__EXCLUDE_LOGS__) +#define TRACCC_WARNING_DEVICE(x, ...) TRACCC_WARNING_HOST_DEVICE(x, __VA_ARGS__) +#define TRACCC_ERROR_DEVICE(x, ...) TRACCC_ERROR_HOST_DEVICE(x, __VA_ARGS__) +#define TRACCC_FATAL_DEVICE(x, ...) TRACCC_FATAL_HOST_DEVICE(x, __VA_ARGS__) +#define TRACCC_INFO_DEVICE(x, ...) TRACCC_INFO_HOST_DEVICE(x, __VA_ARGS__) +#define TRACCC_VERBOSE_DEVICE(x, ...) TRACCC_VERBOSE_HOST_DEVICE(x, __VA_ARGS__) +#define TRACCC_DEBUG_DEVICE(x, ...) TRACCC_DEBUG_HOST_DEVICE(x, __VA_ARGS__) +#else +#define TRACCC_WARNING_DEVICE(x, ...) +#define TRACCC_ERROR_DEVICE(x, ...) +#define TRACCC_FATAL_DEVICE(x, ...) +#define TRACCC_INFO_DEVICE(x, ...) +#define TRACCC_VERBOSE_DEVICE(x, ...) +#define TRACCC_DEBUG_DEVICE(x, ...) +#endif diff --git a/core/include/traccc/utils/trigonometric_helpers.hpp b/core/include/traccc/utils/trigonometric_helpers.hpp new file mode 100644 index 0000000000..e9a15b4a27 --- /dev/null +++ b/core/include/traccc/utils/trigonometric_helpers.hpp @@ -0,0 +1,83 @@ +/** TRACCC library, part of the ACTS project (R&D line) + * + * (c) 2025 CERN for the benefit of the ACTS project + * + * Mozilla Public License Version 2.0 + */ + +#pragma once + +// Project include(s). +#include "traccc/definitions/common.hpp" +#include "traccc/definitions/math.hpp" +#include "traccc/definitions/primitives.hpp" +#include "traccc/definitions/qualifiers.hpp" + +/// @see +/// https://github.com/acts-project/acts/blob/8098e6953ac35771c34a4e3b13dbfee50869a0c1/Core/include/Acts/Utilities/detail/periodic.hpp +namespace traccc::detail { + +/// Wrap a periodic value back into the nominal range. +TRACCC_HOST_DEVICE +constexpr traccc::scalar wrap_periodic(traccc::scalar value, + traccc::scalar start, + traccc::scalar range) { + // only wrap if really necessary + const traccc::scalar diff{value - start}; + return ((0 <= diff) && (diff < range)) + ? value + : (value - range * math::floor(diff / range)); +} + +/// Calculate the equivalent angle in the [-pi, pi) range. +TRACCC_HOST_DEVICE +constexpr traccc::scalar wrap_phi(traccc::scalar phi) { + constexpr traccc::scalar PI{traccc::constant::pi}; + constexpr traccc::scalar TWOPI{2.f * traccc::constant::pi}; + return wrap_periodic(phi, -PI, TWOPI); +} + +/// Calculate the equivalent angle in the [0, 2*pi) range. +TRACCC_HOST_DEVICE +constexpr traccc::scalar wrap_theta(traccc::scalar theta) { + constexpr traccc::scalar TWOPI{2.f * traccc::constant::pi}; + return wrap_periodic(theta, 0.f, TWOPI); +} + +/// Ensure both phi and theta direction angles are within the allowed range. +/// +/// @param[in] phi Transverse direction angle +/// @param[in] theta Longitudinal direction angle +/// @return pair containing the updated angles +/// +/// The phi angle is truly cyclic, i.e. all values outside the nominal range +/// [-pi,pi) have a corresponding value inside nominal range, independent from +/// the theta angle. The theta angle is more complicated. Imagine that the two +/// angles describe a position on the unit sphere. If theta moves outside its +/// nominal range [0,pi], we are moving over one of the two poles of the unit +/// sphere along the great circle defined by phi. The angles still describe a +/// valid position on the unit sphere, but to describe it with angles within +/// their nominal range, both phi and theta need to be updated; when moving over +/// the poles, phi needs to be flipped by 180degree to allow theta to remain +/// within its nominal range. +constexpr std::pair wrap_phi_theta( + traccc::scalar phi, traccc::scalar theta) { + constexpr traccc::scalar PI{traccc::constant::pi}; + constexpr traccc::scalar TWOPI{2.f * traccc::constant::pi}; + + // wrap to [0,2pi). while the nominal range of theta is [0,pi], it is + // periodic, i.e. describes identical positions, in the full [0,2pi) range. + // moving it first to the periodic range simplifies further steps as the + // possible range of theta becomes fixed. + theta = wrap_theta(theta); + if (PI < theta) { + // theta is in the second half of the great circle and outside its + // nominal range. need to change both phi and theta to be within range. + phi += PI; + theta = TWOPI - theta; + } + + return {wrap_phi(phi), theta}; +} + +} // namespace traccc::detail diff --git a/device/alpaka/src/finding/combinatorial_kalman_filter.hpp b/device/alpaka/src/finding/combinatorial_kalman_filter.hpp index b5176c0379..68bc6fe93d 100644 --- a/device/alpaka/src/finding/combinatorial_kalman_filter.hpp +++ b/device/alpaka/src/finding/combinatorial_kalman_filter.hpp @@ -28,7 +28,6 @@ #include "traccc/finding/device/propagate_to_next_surface.hpp" #include "traccc/finding/device/remove_duplicates.hpp" #include "traccc/finding/finding_config.hpp" -#include "traccc/utils/logging.hpp" #include "traccc/utils/memory_resource.hpp" #include "traccc/utils/projections.hpp" #include "traccc/utils/propagation.hpp" diff --git a/device/common/include/traccc/finding/device/impl/find_tracks.ipp b/device/common/include/traccc/finding/device/impl/find_tracks.ipp index 52ca4aec37..4e806e3b76 100644 --- a/device/common/include/traccc/finding/device/impl/find_tracks.ipp +++ b/device/common/include/traccc/finding/device/impl/find_tracks.ipp @@ -26,6 +26,7 @@ #include "traccc/fitting/kalman_filter/gain_matrix_updater.hpp" #include "traccc/fitting/kalman_filter/is_line_visitor.hpp" #include "traccc/fitting/status_codes.hpp" +#include "traccc/utils/logging.hpp" // Detray include(s) #include @@ -104,6 +105,8 @@ TRACCC_HOST_DEVICE inline void find_tracks( in_params.at(in_param_id).surface_link().index()}; init_meas = sf_idx == 0u ? 0u : meas_ranges[sf_idx - 1]; num_meas = meas_ranges[sf_idx] - init_meas; + TRACCC_VERBOSE_DEVICE("No. measurements: %d", num_meas); + TRACCC_VERBOSE_DEVICE("Testing measurement: %d", init_meas); /* * If we cannot find any corresponding measurements, set the number of @@ -209,13 +212,15 @@ TRACCC_HOST_DEVICE inline void find_tracks( const detray::tracking_surface sf{det, in_par.surface_link()}; - const bool is_line = sf.template visit_mask(); + const bool is_line = detail::is_line(sf); // Run the Kalman update const kalman_fitter_status res = gain_matrix_updater{}( trk_state, measurements, in_par, is_line); + TRACCC_DEBUG_DEVICE("KF status: %d", res); + /* * The $\chi^2$ value from the Kalman update should be less than * `chi2_max`, and the fit should have succeeded. If both @@ -229,11 +234,12 @@ TRACCC_HOST_DEVICE inline void find_tracks( */ if (const traccc::scalar chi2 = trk_state.filtered_chi2(); res != kalman_fitter_status::SUCCESS || - chi2 >= cfg.chi2_max) { + !(chi2 <= cfg.chi2_max)) { use_measurement = false; } if (use_measurement) { + TRACCC_VERBOSE_DEVICE("Found measurement: %d", meas_idx); result.emplace(trk_state, meas_idx, owner_local_thread_id); } } @@ -527,8 +533,8 @@ TRACCC_HOST_DEVICE inline void find_tracks( } /* - * Compute the offset at which this block will write, as well as the index - * at which this block will write. + * Compute the offset at which this block will write, as well as the + * index at which this block will write. */ if (in_param_is_live) { local_num_params = std::get<1>(decode_insertion_mutex( @@ -600,6 +606,8 @@ TRACCC_HOST_DEVICE inline void find_tracks( .chi2_sum = prev_chi2_sum, .ndf_sum = prev_ndf_sum}; + TRACCC_VERBOSE_DEVICE("Hole state created"); + unsigned int param_pos = out_offset - payload.curr_links_idx; out_params.at(param_pos) = in_params.at(in_param_id); @@ -609,6 +617,12 @@ TRACCC_HOST_DEVICE inline void find_tracks( const unsigned int n_cands = payload.step - n_skipped; if (n_cands >= cfg.min_track_candidates_per_track) { + if (!in_param_can_create_hole) { + TRACCC_WARNING_DEVICE("Create tip: Max no. holes"); + } else { + TRACCC_WARNING_DEVICE( + "Create tip: No next sensitive found"); + } auto tip_pos = tips.push_back(prev_link_idx); tip_lengths.at(tip_pos) = n_cands; } @@ -633,6 +647,7 @@ TRACCC_HOST_DEVICE inline void find_tracks( if (last_step && n_cands >= cfg.min_track_candidates_per_track) { + TRACCC_ERROR_DEVICE("Create tip: Max no. candidates"); auto tip_pos = tips.push_back(out_offset); tip_lengths.at(tip_pos) = n_cands; } 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 index eb4ddc91db..54d15271c1 100644 --- 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 @@ -8,6 +8,7 @@ #pragma once // Project include(s). +#include "traccc/utils/logging.hpp" #include "traccc/utils/particle.hpp" // Detray include(s). @@ -62,7 +63,9 @@ TRACCC_HOST_DEVICE inline void propagate_to_next_surface( const bound_track_parameters<> in_par = params.at(param_id); // Create propagator - propagator_t propagator(cfg.propagation); + auto prop_cfg{cfg.propagation}; + prop_cfg.navigation.estimate_scattering_noise = false; + propagator_t propagator(prop_cfg); // Create propagator state typename propagator_t::state propagation(in_par, payload.field_data, det); @@ -87,21 +90,20 @@ TRACCC_HOST_DEVICE inline void propagate_to_next_surface( typename detray::detail::tuple_element<2, actor_tuple_type>::type::state s2{ s3}; // Parameter resetter - // typename detray::detail::tuple_element<4, actor_tuple_type>::type::state - // s4{ - // cfg.propagation}; + typename detray::detail::tuple_element<4, actor_tuple_type>::type::state s4{ + prop_cfg}; // Momentum aborter typename detray::detail::tuple_element<5, actor_tuple_type>::type::state s5; // CKF aborter typename detray::detail::tuple_element<6, actor_tuple_type>::type::state s6; - s6.min_step_length = cfg.min_step_length_for_next_surface; - s6.max_count = cfg.max_step_counts_for_next_surface; s5.min_pT(static_cast(cfg.min_pT)); s5.min_p(static_cast(cfg.min_p)); + s6.min_step_length = cfg.min_step_length_for_next_surface; + s6.max_count = cfg.max_step_counts_for_next_surface; // Propagate to the next surface - propagator.propagate(propagation, detray::tie(s0, s2, s3, s5, s6)); + propagator.propagate(propagation, detray::tie(s0, s2, s3, s4, s5, s6)); // If a surface found, add the parameter for the next step if (s6.success) { @@ -110,13 +112,32 @@ TRACCC_HOST_DEVICE inline void propagate_to_next_surface( params[param_id] = propagation._stepping.bound_params(); params_liveness[param_id] = 1u; + + const scalar theta = params[param_id].theta(); + if (theta <= 0.f || theta >= 2.f * constant::pi) { + TRACCC_ERROR_DEVICE("Theta is zero after propagation"); + params_liveness[param_id] = 0u; + } + + if (!std::isfinite(params[param_id].phi())) { + TRACCC_ERROR_DEVICE( + "Phi is infinite after propagation (Matrix inversion)"); + params_liveness[param_id] = 0u; + } + + if (math::fabs(params[param_id].qop()) == 0.f) { + TRACCC_ERROR_DEVICE("q/p is zero after propagation"); + params_liveness[param_id] = 0u; + } } else { params_liveness[param_id] = 0u; + } - if (n_cands >= cfg.min_track_candidates_per_track) { - auto tip_pos = tips.push_back(link_idx); - tip_lengths.at(tip_pos) = n_cands; - } + if (params_liveness[param_id] == 0 && + n_cands >= cfg.min_track_candidates_per_track) { + TRACCC_VERBOSE_DEVICE("Create tip: No next sensitive found"); + auto tip_pos = tips.push_back(link_idx); + tip_lengths.at(tip_pos) = n_cands; } } diff --git a/device/common/include/traccc/fitting/device/fit_backward.hpp b/device/common/include/traccc/fitting/device/fit_backward.hpp index 70305decf9..83d2c90078 100644 --- a/device/common/include/traccc/fitting/device/fit_backward.hpp +++ b/device/common/include/traccc/fitting/device/fit_backward.hpp @@ -39,7 +39,8 @@ TRACCC_HOST_DEVICE inline void fit_backward( if (param_liveness.at(param_id) > 0u) { typename fitter_t::state fitter_state( track, tracks.states, tracks.measurements, - *(payload.barcodes_view.ptr() + param_id)); + *(payload.barcodes_view.ptr() + param_id), + fitter.config().propagation); kalman_fitter_status fit_status = fitter.smooth(fitter_state); diff --git a/device/common/include/traccc/fitting/device/fit_forward.hpp b/device/common/include/traccc/fitting/device/fit_forward.hpp index 5b3495e05c..c67c035fb8 100644 --- a/device/common/include/traccc/fitting/device/fit_forward.hpp +++ b/device/common/include/traccc/fitting/device/fit_forward.hpp @@ -42,7 +42,7 @@ TRACCC_HOST_DEVICE inline void fit_forward( typename fitter_t::state fitter_state( track, tracks.states, tracks.measurements, - *(payload.barcodes_view.ptr() + param_id)); + *(payload.barcodes_view.ptr() + param_id), fitter.config().propagation); kalman_fitter_status fit_status = fitter.filter(params, fitter_state); diff --git a/examples/options/src/track_propagation.cpp b/examples/options/src/track_propagation.cpp index c1221517f5..881d544754 100644 --- a/examples/options/src/track_propagation.cpp +++ b/examples/options/src/track_propagation.cpp @@ -52,6 +52,25 @@ track_propagation::track_propagation() "search-window", po::value(&m_search_window)->default_value(m_search_window), "Size of the grid surface search window"); + m_desc.add_options()( + "mask-tolerance-scaling", + po::value(&(m_config.navigation.mask_tolerance_scalor)) + ->default_value(m_config.navigation.mask_tolerance_scalor), + "Scale factor between min. and max. mask tolerance with surface " + "distance"); + + m_desc.add_options()( + "accumulated-noise-factor", + po::value(&(m_config.navigation.accumulated_error)) + ->default_value(m_config.navigation.accumulated_error), + "Scale factor on the total track path length to model accumualted " + "noise [%]"); + + m_desc.add_options()( + "scattering-stddevs", + po::value(&(m_config.navigation.n_scattering_stddev)) + ->default_value(m_config.navigation.n_scattering_stddev), + "Number of angle standard deviations to use for the noise modelling"); m_desc.add_options()("rk-tolerance-mm", po::value(&(m_config.stepping.rk_error_tol)) ->default_value(m_config.stepping.rk_error_tol / @@ -101,6 +120,7 @@ void track_propagation::read(const po::variables_map &) { m_config.navigation.min_mask_tolerance *= traccc::unit::mm; m_config.navigation.max_mask_tolerance *= traccc::unit::mm; m_config.navigation.search_window = m_search_window; + m_config.navigation.accumulated_error /= 100.f; m_config.stepping.min_stepsize *= traccc::unit::mm; m_config.stepping.path_limit *= traccc::unit::m; @@ -125,7 +145,7 @@ std::unique_ptr track_propagation::as_printable() traccc::unit::mm) + " mm")); cat_nav->add_child(std::make_unique( - "Mask tolerance scalar", + "Mask tolerance scaling", std::to_string(m_config.navigation.mask_tolerance_scalor))); cat_nav->add_child(std::make_unique( "Path tolerance", std::to_string(m_config.navigation.path_tolerance / @@ -140,6 +160,12 @@ std::unique_ptr track_propagation::as_printable() "Search window", std::to_string(m_config.navigation.search_window[0]) + " x " + std::to_string(m_config.navigation.search_window[1]))); + cat_nav->add_child(std::make_unique( + "Scale factor for accumulated noise", + std::to_string(m_config.navigation.accumulated_error * 100.f) + " %")); + cat_nav->add_child(std::make_unique( + "# scattering stddevs to assume", + std::to_string(m_config.navigation.n_scattering_stddev))); auto cat_tsp = std::make_unique("Transport"); diff --git a/examples/run/cuda/truth_finding_example_cuda.cpp b/examples/run/cuda/truth_finding_example_cuda.cpp index 136517a9bc..1e5248dfa7 100644 --- a/examples/run/cuda/truth_finding_example_cuda.cpp +++ b/examples/run/cuda/truth_finding_example_cuda.cpp @@ -190,6 +190,8 @@ int seq_run(const traccc::opts::track_finding& finding_opts, seeds.push_back(truth_track_candidates.tracks.at(i_trk).params()); } + std::cout << "Number of seeds: " << seeds.size() << std::endl; + traccc::bound_track_parameters_collection_types::buffer seeds_buffer{ static_cast(seeds.size()), mr.main}; async_copy.setup(seeds_buffer)->wait(); diff --git a/extern/algebra-plugins/CMakeLists.txt b/extern/algebra-plugins/CMakeLists.txt index 15869a9522..7707510a75 100644 --- 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.28.0.tar.gz;URL_MD5;24fa671f564a332858599df60fb2228f" + "URL;https://github.com/acts-project/algebra-plugins/archive/refs/tags/v0.30.0.tar.gz;URL_MD5;d1ec191c6fea5bd0e73f62d4b0319b8c" 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 index c3ce19a867..17263f5da3 100644 --- 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.103.0.tar.gz;URL_MD5;7df38072d676b63ee3f0f88bd84c0106" + "URL;https://github.com/acts-project/detray/archive/refs/tags/v0.104.1.tar.gz;URL_MD5;2233a3b9945122e1b3f16ac973c7d1be" 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} ) @@ -22,6 +22,10 @@ FetchContent_Declare( Detray SYSTEM ${TRACCC_DETRAY_SOURCE} ) set( DETRAY_CUSTOM_SCALARTYPE "${TRACCC_CUSTOM_SCALARTYPE}" CACHE STRING "Scalar type to use in the Detray code" ) +set( DETRAY_SET_LOGGING "${TRACCC_DEVICE_LOG_LVL}" CACHE STRING + "Log level to use in device code" FORCE) +message( STATUS "Setting DETRAY and TRACCC device log level to ${DETRAY_SET_LOGGING}" ) + set( DETRAY_BUILD_UNITTESTS FALSE CACHE BOOL "Turn off the build of the Detray unit tests" ) set( DETRAY_BUILD_INTEGRATIONTESTS FALSE CACHE BOOL diff --git a/simulation/include/traccc/simulation/simulator.hpp b/simulation/include/traccc/simulation/simulator.hpp index d6ba2f8183..bbc9043fe0 100644 --- a/simulation/include/traccc/simulation/simulator.hpp +++ b/simulation/include/traccc/simulation/simulator.hpp @@ -82,6 +82,15 @@ struct simulator { m_cfg.ptc_type = ptc_type; m_track_generator->config().charge(ptc_type.charge()); + + // Turn off tracking features + m_cfg.propagation.stepping.do_covariance_transport = false; + m_cfg.propagation.stepping.use_eloss_gradient = false; + m_cfg.propagation.stepping.use_field_gradient = false; + m_cfg.propagation.navigation.estimate_scattering_noise = false; + + m_resetter = typename detray::parameter_resetter::state{ + m_cfg.propagation}; } config& get_config() { return m_cfg; } @@ -108,8 +117,8 @@ struct simulator { m_scatterer.set_seed(event_id); writer_state.set_seed(event_id); - auto actor_states = - detray::tie(m_aborter_state, m_scatterer, writer_state); + auto actor_states = detray::tie(m_aborter_state, m_scatterer, + m_resetter, writer_state); for (auto track : *m_track_generator.get()) { @@ -149,6 +158,7 @@ struct simulator { /// Actor states typename detray::momentum_aborter::state m_aborter_state{}; typename detray::random_scatterer::state m_scatterer{}; + typename detray::parameter_resetter::state m_resetter{}; }; } // namespace traccc diff --git a/tests/cpu/test_kalman_fitter_momentum_resolution.cpp b/tests/cpu/test_kalman_fitter_momentum_resolution.cpp index 314653cb78..24c609bb97 100644 --- a/tests/cpu/test_kalman_fitter_momentum_resolution.cpp +++ b/tests/cpu/test_kalman_fitter_momentum_resolution.cpp @@ -187,9 +187,9 @@ TEST_P(KalmanFittingMomentumResolutionTests, Run) { // n_trakcs = 100 ASSERT_GE(static_cast(n_tracks), - 0.95 * static_cast(n_truth_tracks)); + 0.94 * static_cast(n_truth_tracks)); ASSERT_GE(static_cast(n_fitted_tracks), - 0.95 * static_cast(n_truth_tracks)); + 0.94 * static_cast(n_truth_tracks)); for (std::size_t i_trk = 0; i_trk < n_tracks; i_trk++) { @@ -242,7 +242,7 @@ TEST_P(KalmanFittingMomentumResolutionTests, Run) { float success_rate = static_cast(n_success) / static_cast(n_truth_tracks * n_events); - ASSERT_GE(success_rate, 0.98f); + ASSERT_GE(success_rate, 0.97f); } // Muon with 1, 10, 100 GeV/c, no materials diff --git a/tests/cpu/test_kalman_fitter_wire_chamber.cpp b/tests/cpu/test_kalman_fitter_wire_chamber.cpp index 4fa6e077d4..38ac359b18 100644 --- a/tests/cpu/test_kalman_fitter_wire_chamber.cpp +++ b/tests/cpu/test_kalman_fitter_wire_chamber.cpp @@ -132,6 +132,8 @@ TEST_P(KalmanFittingWireChamberTests, Run) { fit_cfg.propagation.navigation.min_mask_tolerance = static_cast(mask_tolerance); fit_cfg.propagation.navigation.search_window = search_window; + // TODO: Disable until overlaps are handled correctly + fit_cfg.propagation.navigation.estimate_scattering_noise = false; fit_cfg.ptc_hypothesis = ptc; traccc::host::kalman_fitting_algorithm fitting(fit_cfg, host_mr, copy);