From 7436f407c06b00f785d4fd652d8662c7efdea029 Mon Sep 17 00:00:00 2001 From: Andreas Stefl Date: Wed, 16 Oct 2024 14:51:44 +0200 Subject: [PATCH] fix: Implement `DirectNavigator` direction handling (#3702) Our current implementation of the `DirectNavigator` does not handle the propagation direction which is essential if we for example reverse direction to smooth track in the Kalman filter. Apart from that I added a search mechanism in the initialization for the starting surface. --- This pull request to `Core/include/Acts/Propagator/DirectNavigator.hpp` includes changes to improve the navigation logic in the `DirectNavigator` class. The most important changes include adding new utility methods for surface navigation, updating the surface index handling, and enhancing logging for better debugging. Enhancements to navigation logic: * Added utility methods `navSurface`, `nextSurface`, `endOfSurfaces`, and `remainingSurfaces` to the `DirectNavigator` class to streamline surface navigation. (`Core/include/Acts/Propagator/DirectNavigator.hpp`) [[1]](diffhunk://#diff-56afaecc4b563cfe416f8615f35dd5ac21f706869dd2c164dd326087b3edd655R77-R100) [[2]](diffhunk://#diff-56afaecc4b563cfe416f8615f35dd5ac21f706869dd2c164dd326087b3edd655L138-R192) [[3]](diffhunk://#diff-56afaecc4b563cfe416f8615f35dd5ac21f706869dd2c164dd326087b3edd655L166-R215) [[4]](diffhunk://#diff-56afaecc4b563cfe416f8615f35dd5ac21f706869dd2c164dd326087b3edd655L186-R229) [[5]](diffhunk://#diff-56afaecc4b563cfe416f8615f35dd5ac21f706869dd2c164dd326087b3edd655L227-R279) [[6]](diffhunk://#diff-56afaecc4b563cfe416f8615f35dd5ac21f706869dd2c164dd326087b3edd655L255-R299) Updates to surface index handling: * Changed the type of `surfaceIndex` from `std::size_t` to `int` to support bidirectional navigation. (`Core/include/Acts/Propagator/DirectNavigator.hpp`) Logging improvements: * Enhanced logging in the `initialize` method to provide more detailed information about the start surface and surface sequence. (`Core/include/Acts/Propagator/DirectNavigator.hpp`) --- .../Acts/Propagator/DirectNavigator.hpp | 99 +++++++++++++++---- .../Acts/Propagator/MultiStepperAborters.hpp | 5 +- Core/include/Acts/Propagator/Navigator.hpp | 2 - .../Acts/Propagator/TryAllNavigator.hpp | 4 - .../Acts/TrackFitting/GaussianSumFitter.hpp | 6 +- 5 files changed, 82 insertions(+), 34 deletions(-) diff --git a/Core/include/Acts/Propagator/DirectNavigator.hpp b/Core/include/Acts/Propagator/DirectNavigator.hpp index 8df2afdfe8e..eb9db91e71f 100644 --- a/Core/include/Acts/Propagator/DirectNavigator.hpp +++ b/Core/include/Acts/Propagator/DirectNavigator.hpp @@ -8,6 +8,7 @@ #pragma once +#include "Acts/Definitions/Direction.hpp" #include "Acts/Geometry/BoundarySurfaceT.hpp" #include "Acts/Geometry/Layer.hpp" #include "Acts/Geometry/TrackingGeometry.hpp" @@ -17,6 +18,7 @@ #include "Acts/Surfaces/BoundaryTolerance.hpp" #include "Acts/Surfaces/Surface.hpp" #include "Acts/Utilities/Intersection.hpp" +#include "Acts/Utilities/Logger.hpp" #include #include @@ -63,7 +65,7 @@ class DirectNavigator { Options options; /// Index of the next surface to try - std::size_t surfaceIndex = 0; + int surfaceIndex = 0; /// Navigation state - external interface: the current surface const Surface* currentSurface = nullptr; @@ -72,6 +74,35 @@ class DirectNavigator { bool targetReached = false; /// Navigation state - external interface: a break has been detected bool navigationBreak = false; + + const Surface* navSurface() const { + return options.surfaces.at(surfaceIndex); + } + + void nextSurface(Direction direction) { + if (direction == Direction::Forward) { + ++surfaceIndex; + } else { + --surfaceIndex; + } + } + + bool endOfSurfaces() const { + return surfaceIndex < 0 || + surfaceIndex >= static_cast(options.surfaces.size()); + } + + int remainingSurfaces(Direction direction) const { + if (direction == Direction::Forward) { + return options.surfaces.size() - surfaceIndex; + } + return surfaceIndex + 1; + } + + void resetSurfaceIndex(Direction direction) { + surfaceIndex = + direction == Direction::Forward ? 0 : options.surfaces.size() - 1; + } }; explicit DirectNavigator(std::unique_ptr _logger = @@ -135,17 +166,41 @@ class DirectNavigator { void initialize(propagator_state_t& state, const stepper_t& /*stepper*/) const { ACTS_VERBOSE("Initialize. Surface sequence for navigation:"); - for (auto surface : state.navigation.options.surfaces) { + for (const Surface* surface : state.navigation.options.surfaces) { ACTS_VERBOSE(surface->geometryId() << " - " << surface->center(state.geoContext).transpose()); } // We set the current surface to the start surface state.navigation.currentSurface = state.navigation.options.startSurface; - if (state.navigation.currentSurface) { + if (state.navigation.currentSurface != nullptr) { ACTS_VERBOSE("Current surface set to start surface " << state.navigation.currentSurface->geometryId()); + } else { + ACTS_VERBOSE("Current surface set to nullptr"); + } + + // Reset the surface index + state.navigation.resetSurfaceIndex(state.options.direction); + for (const Surface* surface : state.navigation.options.surfaces) { + // make sure we skip over the start surface + state.navigation.nextSurface(state.options.direction); + if (surface == state.navigation.currentSurface) { + break; + } + } + ACTS_VERBOSE("Start surface index set to " + << state.navigation.surfaceIndex); + if (state.navigation.endOfSurfaces()) { + ACTS_DEBUG( + "Did not find the start surface in the sequence. Assuming it is not " + "part of the sequence. Trusting the correctness of the input " + "sequence. Resetting the surface index."); + state.navigation.resetSurfaceIndex(state.options.direction); } + + state.navigation.navigationBreak = false; + state.navigation.targetReached = false; } /// @brief Navigator pre step call @@ -157,21 +212,25 @@ class DirectNavigator { /// @param [in] stepper Stepper in use template void preStep(propagator_state_t& state, const stepper_t& stepper) const { + if (state.navigation.navigationBreak) { + return; + } + ACTS_VERBOSE("pre step"); // Navigator target always resets the current surface state.navigation.currentSurface = nullptr; // Output the position in the sequence - ACTS_VERBOSE((state.navigation.options.surfaces.size() - - state.navigation.surfaceIndex) + ACTS_VERBOSE(state.navigation.remainingSurfaces(state.options.direction) << " out of " << state.navigation.options.surfaces.size() << " surfaces remain to try."); - if (state.navigation.surfaceIndex >= - state.navigation.options.surfaces.size()) { + if (state.navigation.endOfSurfaces()) { // Set the navigation break + ACTS_VERBOSE("End of surfaces reached, navigation break."); state.navigation.navigationBreak = true; + stepper.releaseStepSize(state.stepping, ConstrainedStep::actor); // If no externally provided target is given, the target is reached if (state.navigation.options.targetSurface == nullptr) { state.navigation.targetReached = true; @@ -183,8 +242,7 @@ class DirectNavigator { // Establish & update the surface status // TODO we do not know the intersection index - passing the closer one - const auto& surface = - *state.navigation.options.surfaces.at(state.navigation.surfaceIndex); + const auto& surface = *state.navigation.navSurface(); const double farLimit = std::numeric_limits::max(); const auto index = chooseIntersection( @@ -202,7 +260,7 @@ class DirectNavigator { "Surface not reachable anymore, switching to next one in " "sequence"); // Move the sequence to the next surface - ++state.navigation.surfaceIndex; + state.navigation.nextSurface(state.options.direction); } else { ACTS_VERBOSE("Navigation stepSize set to " << stepper.outputStepSize(state.stepping)); @@ -218,26 +276,27 @@ class DirectNavigator { /// @param [in] stepper Stepper in use template void postStep(propagator_state_t& state, const stepper_t& stepper) const { + if (state.navigation.navigationBreak) { + return; + } + ACTS_VERBOSE("post step"); // Navigator post step always resets the current surface state.navigation.currentSurface = nullptr; // Output the position in the sequence - ACTS_VERBOSE((state.navigation.options.surfaces.size() - - state.navigation.surfaceIndex) + ACTS_VERBOSE(state.navigation.remainingSurfaces(state.options.direction) << " out of " << state.navigation.options.surfaces.size() << " surfaces remain to try."); - if (state.navigation.surfaceIndex >= - state.navigation.options.surfaces.size()) { + if (state.navigation.endOfSurfaces()) { return; } // Establish the surface status // TODO we do not know the intersection index - passing the closer one - const auto& surface = - *state.navigation.options.surfaces.at(state.navigation.surfaceIndex); + const auto& surface = *state.navigation.navSurface(); const double farLimit = std::numeric_limits::max(); const auto index = chooseIntersection( @@ -252,14 +311,12 @@ class DirectNavigator { *m_logger); if (surfaceStatus == Intersection3D::Status::onSurface) { // Set the current surface - state.navigation.currentSurface = - state.navigation.options.surfaces.at(state.navigation.surfaceIndex); + state.navigation.currentSurface = state.navigation.navSurface(); ACTS_VERBOSE("Current surface set to " << state.navigation.currentSurface->geometryId()); // Move the sequence to the next surface - ++state.navigation.surfaceIndex; - if (state.navigation.surfaceIndex < - state.navigation.options.surfaces.size()) { + state.navigation.nextSurface(state.options.direction); + if (!state.navigation.endOfSurfaces()) { ACTS_VERBOSE("Next surface candidate is " << state.navigation.options.surfaces .at(state.navigation.surfaceIndex) diff --git a/Core/include/Acts/Propagator/MultiStepperAborters.hpp b/Core/include/Acts/Propagator/MultiStepperAborters.hpp index f675200975e..0e28375df94 100644 --- a/Core/include/Acts/Propagator/MultiStepperAborters.hpp +++ b/Core/include/Acts/Propagator/MultiStepperAborters.hpp @@ -71,8 +71,8 @@ struct MultiStepperSurfaceReached : public SurfaceReached { } ACTS_VERBOSE( - "MultiStepperSurfaceReached aborter | " - "Target intersection not found. Maybe next time?"); + "MultiStepperSurfaceReached aborter | Average distance to target: " + << sIntersection.pathLength()); } bool reached = true; @@ -84,6 +84,7 @@ struct MultiStepperSurfaceReached : public SurfaceReached { if (!SurfaceReached::checkAbort(singleState, singleStepper, navigator, logger)) { + cmp.status() = Acts::Intersection3D::Status::reachable; reached = false; } else { cmp.status() = Acts::Intersection3D::Status::onSurface; diff --git a/Core/include/Acts/Propagator/Navigator.hpp b/Core/include/Acts/Propagator/Navigator.hpp index b2bea757749..c2a9427257c 100644 --- a/Core/include/Acts/Propagator/Navigator.hpp +++ b/Core/include/Acts/Propagator/Navigator.hpp @@ -205,8 +205,6 @@ class Navigator { : m_cfg{std::move(cfg)}, m_logger{std::move(_logger)} {} State makeState(const Options& options) const { - assert(options.startSurface != nullptr && "Start surface must be set"); - State state; state.options = options; state.startSurface = options.startSurface; diff --git a/Core/include/Acts/Propagator/TryAllNavigator.hpp b/Core/include/Acts/Propagator/TryAllNavigator.hpp index 9a3a2d78475..83d7cbc64a4 100644 --- a/Core/include/Acts/Propagator/TryAllNavigator.hpp +++ b/Core/include/Acts/Propagator/TryAllNavigator.hpp @@ -96,8 +96,6 @@ class TryAllNavigatorBase { : m_cfg(std::move(cfg)), m_logger{std::move(_logger)} {} State makeState(const Options& options) const { - assert(options.startSurface != nullptr && "Start surface must be set"); - State state; state.options = options; state.startSurface = options.startSurface; @@ -587,8 +585,6 @@ class TryAllOverstepNavigator : public TryAllNavigatorBase { : TryAllNavigatorBase(std::move(cfg), std::move(logger)) {} State makeState(const Options& options) const { - assert(options.startSurface != nullptr && "Start surface must be set"); - State state; state.options = options; state.startSurface = options.startSurface; diff --git a/Core/include/Acts/TrackFitting/GaussianSumFitter.hpp b/Core/include/Acts/TrackFitting/GaussianSumFitter.hpp index efe21b58a60..f6c990a696d 100644 --- a/Core/include/Acts/TrackFitting/GaussianSumFitter.hpp +++ b/Core/include/Acts/TrackFitting/GaussianSumFitter.hpp @@ -123,15 +123,11 @@ struct GaussianSumFitter { using Actors = ActorList; using PropagatorOptions = typename propagator_t::template Options; - std::vector backwardSequence( - std::next(sSequence.rbegin()), sSequence.rend()); - backwardSequence.push_back(opts.referenceSurface); - PropagatorOptions propOptions(opts.geoContext, opts.magFieldContext); propOptions.setPlainOptions(opts.propagatorPlainOptions); - propOptions.navigation.surfaces = backwardSequence; + propOptions.navigation.surfaces = sSequence; propOptions.actorList.template get() .m_cfg.bethe_heitler_approx = &m_betheHeitlerApproximation;