diff --git a/Core/include/Acts/Fitter/KalmanFitter.hpp b/Core/include/Acts/Fitter/KalmanFitter.hpp
index ddb52be9e2a13d3466dfb3bdbb808d7bfe7de0c5..56c8d1da901fe93070552ce8e6886808c06c376f 100644
--- a/Core/include/Acts/Fitter/KalmanFitter.hpp
+++ b/Core/include/Acts/Fitter/KalmanFitter.hpp
@@ -20,8 +20,8 @@
#include "Acts/MagneticField/MagneticFieldContext.hpp"
#include "Acts/Propagator/AbortList.hpp"
#include "Acts/Propagator/ActionList.hpp"
+#include "Acts/Propagator/ConstrainedStep.hpp"
#include "Acts/Propagator/Propagator.hpp"
-#include "Acts/Propagator/detail/ConstrainedStep.hpp"
#include "Acts/Propagator/detail/StandardAborters.hpp"
#include "Acts/Utilities/CalibrationContext.hpp"
#include "Acts/Utilities/Definitions.hpp"
@@ -396,7 +396,7 @@ class KalmanFitter {
stepper.update(state.stepping, smoothedPars);
// Reverse the propagation direction
state.stepping.stepSize =
- detail::ConstrainedStep(-1. * state.options.maxStepSize);
+ ConstrainedStep(-1. * state.options.maxStepSize);
state.options.direction = backward;
return Result<void>::success();
diff --git a/Core/include/Acts/Propagator/AtlasStepper.hpp b/Core/include/Acts/Propagator/AtlasStepper.hpp
index c1f923dbac2aa879a64ac4cbf4ea642ae2441472..36c5b9b3886c0074ddefd53e99d1d799e684c2a6 100644
--- a/Core/include/Acts/Propagator/AtlasStepper.hpp
+++ b/Core/include/Acts/Propagator/AtlasStepper.hpp
@@ -12,7 +12,8 @@
#include <functional>
#include "Acts/EventData/TrackParameters.hpp"
#include "Acts/Geometry/GeometryContext.hpp"
-#include "Acts/Propagator/detail/ConstrainedStep.hpp"
+#include "Acts/Propagator/ConstrainedStep.hpp"
+#include "Acts/Propagator/detail/SteppingHelper.hpp"
#include "Acts/Surfaces/Surface.hpp"
#include "Acts/Utilities/Definitions.hpp"
#include "Acts/Utilities/Intersection.hpp"
@@ -28,8 +29,6 @@ using namespace Acts::UnitLiterals;
template <typename bfield_t>
class AtlasStepper {
public:
- using cstep = detail::ConstrainedStep;
-
using Jacobian = BoundMatrix;
using Covariance = BoundSymMatrix;
using BoundState = std::tuple<BoundParameters, Jacobian, double>;
@@ -51,11 +50,13 @@ class AtlasStepper {
/// @param[in] pars Input parameters
/// @param[in] ndir The navigation direction w.r.t. parameters
/// @param[in] ssize the steps size limitation
+ /// @param [in] stolerance is the stepping tolerance
template <typename Parameters>
State(std::reference_wrapper<const GeometryContext> gctx,
std::reference_wrapper<const MagneticFieldContext> mctx,
const Parameters& pars, NavigationDirection ndir = forward,
- double ssize = std::numeric_limits<double>::max())
+ double ssize = std::numeric_limits<double>::max(),
+ double stolerance = s_onSurfaceTolerance)
: navDir(ndir),
useJacobian(false),
step(0.),
@@ -67,6 +68,7 @@ class AtlasStepper {
covariance(nullptr),
t0(pars.time()),
stepSize(ndir * std::abs(ssize)),
+ tolerance(stolerance),
fieldCache(mctx),
geoContext(gctx) {
// The rest of this constructor is copy&paste of AtlasStepper::update() -
@@ -271,8 +273,14 @@ class AtlasStepper {
// Starting time
const double t0;
- // adaptive step size of the runge-kutta integration
- cstep stepSize = std::numeric_limits<double>::max();
+ // Adaptive step size of the runge-kutta integration
+ ConstrainedStep stepSize = std::numeric_limits<double>::max();
+
+ // Previous step size for overstep estimation
+ double previousStepSize = 0.;
+
+ /// The tolerance for the stepping
+ double tolerance = s_onSurfaceTolerance;
/// It caches the current magnetic field cell and stays (and interpolates)
/// within as long as this is valid. See step() code for details.
@@ -330,15 +338,59 @@ class AtlasStepper {
/// Time access
double time(const State& state) const { return state.t0 + state.pVector[3]; }
- /// Tests if the state reached a surface
+ /// Update surface status
+ ///
+ /// This method intersect the provided surface and update the navigation
+ /// step estimation accordingly (hence it changes the state). It also
+ /// returns the status of the intersection to trigger onSurface in case
+ /// the surface is reached.
+ ///
+ /// @param state [in,out] The stepping state (thread-local cache)
+ /// @param surface [in] The surface provided
+ /// @param bcheck [in] The boundary check for this status update
+ Intersection::Status updateSurfaceStatus(State& state, const Surface& surface,
+ const BoundaryCheck& bcheck) const {
+ return detail::updateSingleSurfaceStatus<AtlasStepper>(*this, state,
+ surface, bcheck);
+ }
+
+ /// Update step size
+ ///
+ /// It checks the status to the reference surface & updates
+ /// the step size accordingly
+ ///
+ /// @param state [in,out] The stepping state (thread-local cache)
+ /// @param oIntersection [in] The ObjectIntersection to layer, boundary, etc
+ /// @param release [in] boolean to trigger step size release
+ template <typename object_intersection_t>
+ void updateStepSize(State& state, const object_intersection_t& oIntersection,
+ bool release = true) const {
+ detail::updateSingleStepSize<AtlasStepper>(state, oIntersection, release);
+ }
+
+ /// Set Step size - explicitely with a double
+ ///
+ /// @param state [in,out] The stepping state (thread-local cache)
+ /// @param stepSize [in] The step size value
+ /// @param stype [in] The step size type to be set
+ void setStepSize(State& state, double stepSize,
+ ConstrainedStep::Type stype = ConstrainedStep::actor) const {
+ state.previousStepSize = state.stepSize;
+ state.stepSize.update(stepSize, stype, true);
+ }
+
+ /// Release the Step size
///
- /// @param [in] state State that is tests
- /// @param [in] surface Surface that is tested
+ /// @param state [in,out] The stepping state (thread-local cache)
+ void releaseStepSize(State& state) const {
+ state.stepSize.release(ConstrainedStep::actor);
+ }
+
+ /// Output the Step Size - single component
///
- /// @return Boolean statement if surface is reached by state
- bool surfaceReached(const State& state, const Surface* surface) const {
- return surface->isOnSurface(state.geoContext, position(state),
- direction(state), true);
+ /// @param state [in,out] The stepping state (thread-local cache)
+ std::string outputStepSize(const State& state) const {
+ return state.stepSize.toString();
}
/// Create and return the bound state at the current position
diff --git a/Core/include/Acts/Propagator/detail/ConstrainedStep.hpp b/Core/include/Acts/Propagator/ConstrainedStep.hpp
similarity index 89%
rename from Core/include/Acts/Propagator/detail/ConstrainedStep.hpp
rename to Core/include/Acts/Propagator/ConstrainedStep.hpp
index 5826341bbd7d317bd0e2423b6208ce3e8e6538a5..e92fbc8c8b2fd0cf7f34b6611408d26666763a8a 100644
--- a/Core/include/Acts/Propagator/detail/ConstrainedStep.hpp
+++ b/Core/include/Acts/Propagator/ConstrainedStep.hpp
@@ -17,7 +17,6 @@
#include "Acts/Utilities/Definitions.hpp"
namespace Acts {
-namespace detail {
/// A constrained step class for the steppers
struct ConstrainedStep {
@@ -129,12 +128,28 @@ struct ConstrainedStep {
inline std::string ConstrainedStep::toString() const {
std::stringstream dstream;
- dstream << "(" << std::setw(5) << values[accuracy];
- dstream << ", " << std::setw(5) << values[actor];
- dstream << ", " << std::setw(5) << values[aborter];
- dstream << ", " << std::setw(5) << values[user] << " )";
+
+ // Helper method to avoid unreadable screen output
+ auto streamValue = [&](ConstrainedStep cstep) -> void {
+ double val = values[cstep];
+ dstream << std::setw(5);
+ if (std::abs(val) == std::numeric_limits<double>::max()) {
+ dstream << (val > 0 ? "+∞" : "-∞");
+ } else {
+ dstream << val;
+ }
+ };
+
+ dstream << "(";
+ streamValue(accuracy);
+ dstream << ", ";
+ streamValue(actor);
+ dstream << ", ";
+ streamValue(aborter);
+ dstream << ", ";
+ streamValue(user);
+ dstream << " )";
return dstream.str();
}
-} // namespace detail
} // namespace Acts
diff --git a/Core/include/Acts/Propagator/DirectNavigator.hpp b/Core/include/Acts/Propagator/DirectNavigator.hpp
index a71a1a41002f89dab8ba2a40782d6904507d0425..bcf26f9a82660ce748c97f05d8c22efe808320d7 100644
--- a/Core/include/Acts/Propagator/DirectNavigator.hpp
+++ b/Core/include/Acts/Propagator/DirectNavigator.hpp
@@ -17,14 +17,12 @@
#include "Acts/Geometry/Layer.hpp"
#include "Acts/Geometry/TrackingGeometry.hpp"
#include "Acts/Geometry/TrackingVolume.hpp"
+#include "Acts/Propagator/ConstrainedStep.hpp"
#include "Acts/Propagator/Propagator.hpp"
-#include "Acts/Propagator/detail/ConstrainedStep.hpp"
#include "Acts/Surfaces/Surface.hpp"
namespace Acts {
-using Cstep = detail::ConstrainedStep;
-
/// DirectNavigator class
///
/// This is a fully guided navigator that progresses through
@@ -36,7 +34,7 @@ class DirectNavigator {
public:
/// The sequentially crossed surfaces
using SurfaceSequence = std::vector<const Surface*>;
- using SurfaceIter = std::vector<const Surface*>::const_iterator;
+ using SurfaceIter = std::vector<const Surface*>::iterator;
/// Defaulted Constructed
DirectNavigator() = default;
@@ -76,8 +74,6 @@ class DirectNavigator {
state.navigation.surfaceSequence = surfaceSequence;
state.navigation.nextSurfaceIter =
state.navigation.surfaceSequence.begin();
- state.navigation.endSurfaceIter =
- state.navigation.surfaceSequence.end();
r.initialized = true;
}
}
@@ -100,8 +96,6 @@ class DirectNavigator {
/// Iterator the the next surface
SurfaceIter nextSurfaceIter = surfaceSequence.begin();
- /// Iterator to the end for the end sequence trigger
- SurfaceIter endSurfaceIter = surfaceSequence.end();
/// Navigation state - external interface: the start surface
const Surface* startSurface = nullptr;
@@ -130,20 +124,50 @@ class DirectNavigator {
// Navigator status always resets the current surface
state.navigation.currentSurface = nullptr;
+ // Output the position in the sequence
+ debugLog(state, [&] {
+ std::stringstream dstream;
+ dstream << std::distance(state.navigation.nextSurfaceIter,
+ state.navigation.surfaceSequence.end());
+ dstream << " out of " << state.navigation.surfaceSequence.size();
+ dstream << " surfaces remain to try.";
+ return dstream.str();
+ });
// Check if we are on surface
- if (state.navigation.nextSurfaceIter != state.navigation.endSurfaceIter &&
- stepper.surfaceReached(state.stepping,
- *state.navigation.nextSurfaceIter)) {
- // Set the current surface
- state.navigation.currentSurface = *state.navigation.nextSurfaceIter;
- debugLog(state, [&] {
- std::stringstream dstream;
- dstream << "Current surface set to "
- << state.navigation.currentSurface->geoID();
- return dstream.str();
- });
- // Move the sequence to the next surface
- ++state.navigation.nextSurfaceIter;
+ if (state.navigation.nextSurfaceIter !=
+ state.navigation.surfaceSequence.end()) {
+ // Establish the surface status
+ auto surfaceStatus = stepper.updateSurfaceStatus(
+ state.stepping, **state.navigation.nextSurfaceIter, false);
+ if (surfaceStatus == Intersection::Status::onSurface) {
+ // Set the current surface
+ state.navigation.currentSurface = *state.navigation.nextSurfaceIter;
+ debugLog(state, [&] {
+ std::stringstream dstream;
+ dstream << "Current surface set to "
+ << state.navigation.currentSurface->geoID();
+ return dstream.str();
+ });
+ // Move the sequence to the next surface
+ ++state.navigation.nextSurfaceIter;
+ if (state.navigation.nextSurfaceIter !=
+ state.navigation.surfaceSequence.end()) {
+ debugLog(state, [&] {
+ std::stringstream dstream;
+ dstream << "Next surface candidate is "
+ << (*state.navigation.nextSurfaceIter)->geoID();
+ return dstream.str();
+ });
+ stepper.releaseStepSize(state.stepping);
+ }
+ } else if (surfaceStatus == Intersection::Status::reachable) {
+ debugLog(state, [&] {
+ std::stringstream dstream;
+ dstream << "Next surface reachable at distance "
+ << stepper.outputStepSize(state.stepping);
+ return dstream.str();
+ });
+ }
}
}
@@ -161,29 +185,37 @@ class DirectNavigator {
// Navigator target always resets the current surface
state.navigation.currentSurface = nullptr;
-
- if (state.navigation.nextSurfaceIter != state.navigation.endSurfaceIter) {
- // take the target intersection
- auto nextIntersection =
- (*state.navigation.nextSurfaceIter)
- ->intersect(state.geoContext, stepper.position(state.stepping),
- stepper.direction(state.stepping), false);
-
- // Intersect the next surface and go
- double navStep = nextIntersection.intersection.pathLength;
- double overstepLimit = stepper.overstepLimit(state.stepping);
- if (navStep < overstepLimit and nextIntersection.alternative) {
- navStep = nextIntersection.alternative.pathLength;
+ // Output the position in the sequence
+ debugLog(state, [&] {
+ std::stringstream dstream;
+ dstream << std::distance(state.navigation.nextSurfaceIter,
+ state.navigation.surfaceSequence.end());
+ dstream << " out of " << state.navigation.surfaceSequence.size();
+ dstream << " surfaces remain to try.";
+ return dstream.str();
+ });
+ if (state.navigation.nextSurfaceIter !=
+ state.navigation.surfaceSequence.end()) {
+ // Establish & update the surface status
+ auto surfaceStatus = stepper.updateSurfaceStatus(
+ state.stepping, **state.navigation.nextSurfaceIter, false);
+ if (surfaceStatus == Intersection::Status::unreachable) {
+ debugLog(state, [&] {
+ std::stringstream dstream;
+ dstream << "Surface not reachable anymore, switching to next one in "
+ "sequence";
+ return dstream.str();
+ });
+ // Move the sequence to the next surface
+ ++state.navigation.nextSurfaceIter;
+ } else {
+ debugLog(state, [&] {
+ std::stringstream dstream;
+ dstream << "Navigation stepSize set to ";
+ dstream << stepper.outputStepSize(state.stepping);
+ return dstream.str();
+ });
}
-
- // Set the step size - this has to be outsourced to the Stepper
- state.stepping.stepSize.update(navStep, Cstep::actor, true);
- debugLog(state, [&] {
- std::stringstream dstream;
- dstream << "Navigation stepSize released and updated to ";
- dstream << state.stepping.stepSize.toString();
- return dstream.str();
- });
} else {
// Set the navigation break
state.navigation.navigationBreak = true;
diff --git a/Core/include/Acts/Propagator/EigenStepper.hpp b/Core/include/Acts/Propagator/EigenStepper.hpp
index 5be22e9f7f2d367e3683610e755f0d63f23f1ed3..2281518c1d7ca3c8f8ffda7bd4a43d8091108441 100644
--- a/Core/include/Acts/Propagator/EigenStepper.hpp
+++ b/Core/include/Acts/Propagator/EigenStepper.hpp
@@ -18,6 +18,7 @@
#include "Acts/Propagator/EigenStepperError.hpp"
#include "Acts/Propagator/StepperExtensionList.hpp"
#include "Acts/Propagator/detail/Auctioneer.hpp"
+#include "Acts/Propagator/detail/SteppingHelper.hpp"
#include "Acts/Utilities/Intersection.hpp"
#include "Acts/Utilities/Result.hpp"
#include "Acts/Utilities/Units.hpp"
@@ -43,8 +44,6 @@ template <typename bfield_t,
typename auctioneer_t = detail::VoidAuctioneer>
class EigenStepper {
public:
- using cstep = detail::ConstrainedStep;
-
/// Jacobian, Covariance and State defintions
using Jacobian = BoundMatrix;
using Covariance = BoundSymMatrix;
@@ -67,13 +66,15 @@ class EigenStepper {
/// @param [in] par The track parameters at start
/// @param [in] ndir The navigation direciton w.r.t momentum
/// @param [in] ssize is the maximum step size
+ /// @param [in] stolerance is the stepping tolerance
///
/// @note the covariance matrix is copied when needed
template <typename parameters_t>
explicit State(std::reference_wrapper<const GeometryContext> gctx,
std::reference_wrapper<const MagneticFieldContext> mctx,
const parameters_t& par, NavigationDirection ndir = forward,
- double ssize = std::numeric_limits<double>::max())
+ double ssize = std::numeric_limits<double>::max(),
+ double stolerance = s_onSurfaceTolerance)
: pos(par.position()),
dir(par.momentum().normalized()),
p(par.momentum().norm()),
@@ -81,11 +82,13 @@ class EigenStepper {
t0(par.time()),
navDir(ndir),
stepSize(ndir * std::abs(ssize)),
+ tolerance(stolerance),
fieldCache(mctx),
geoContext(gctx) {
// remember the start parameters
startPos = pos;
startDir = dir;
+
// Init the jacobian matrix if needed
if (par.covariance()) {
// Get the reference surface for navigation
@@ -142,11 +145,17 @@ class EigenStepper {
bool covTransport = false;
Covariance cov = Covariance::Zero();
- /// accummulated path length state
+ /// Accummulated path length state
double pathAccumulated = 0.;
- /// adaptive step size of the runge-kutta integration
- cstep stepSize{std::numeric_limits<double>::max()};
+ /// Adaptive step size of the runge-kutta integration
+ ConstrainedStep stepSize{std::numeric_limits<double>::max()};
+
+ /// Last performed step (for overstep limit calculation)
+ double previousStepSize = 0.;
+
+ /// The tolerance for the stepping
+ double tolerance = s_onSurfaceTolerance;
/// This caches the current magnetic field cell and stays
/// (and interpolates) within it as long as this is valid.
@@ -210,20 +219,67 @@ class EigenStepper {
/// @param state [in] The stepping state (thread-local cache)
double time(const State& state) const { return state.t0 + state.dt; }
- /// Overstep limit
+ /// Update surface status
///
- /// @param state [in] The stepping state (thread-local cache)
- double overstepLimit(const State& /*state*/) const { return m_overstepLimit; }
+ /// It checks the status to the reference surface & updates
+ /// the step size accordingly
+ ///
+ /// @param state [in,out] The stepping state (thread-local cache)
+ /// @param surface [in] The surface provided
+ /// @param bcheck [in] The boundary check for this status update
+ Intersection::Status updateSurfaceStatus(State& state, const Surface& surface,
+ const BoundaryCheck& bcheck) const {
+ return detail::updateSingleSurfaceStatus<EigenStepper>(*this, state,
+ surface, bcheck);
+ }
+
+ /// Update step size
+ ///
+ /// This method intersect the provided surface and update the navigation
+ /// step estimation accordingly (hence it changes the state). It also
+ /// returns the status of the intersection to trigger onSurface in case
+ /// the surface is reached.
+ ///
+ /// @param state [in,out] The stepping state (thread-local cache)
+ /// @param oIntersection [in] The ObjectIntersection to layer, boundary, etc
+ /// @param release [in] boolean to trigger step size release
+ template <typename object_intersection_t>
+ void updateStepSize(State& state, const object_intersection_t& oIntersection,
+ bool release = true) const {
+ detail::updateSingleStepSize<EigenStepper>(state, oIntersection, release);
+ }
+
+ /// Set Step size - explicitely with a double
+ ///
+ /// @param state [in,out] The stepping state (thread-local cache)
+ /// @param stepSize [in] The step size value
+ /// @param stype [in] The step size type to be set
+ void setStepSize(State& state, double stepSize,
+ ConstrainedStep::Type stype = ConstrainedStep::actor) const {
+ state.previousStepSize = state.stepSize;
+ state.stepSize.update(stepSize, stype, true);
+ }
- /// Tests if the state reached a surface
+ /// Release the Step size
///
- /// @param [in] state State that is tests
- /// @param [in] surface Surface that is tested
+ /// @param state [in,out] The stepping state (thread-local cache)
+ void releaseStepSize(State& state) const {
+ state.stepSize.release(ConstrainedStep::actor);
+ }
+
+ /// Output the Step Size - single component
///
- /// @return Boolean statement if surface is reached by state
- bool surfaceReached(const State& state, const Surface* surface) const {
- return surface->isOnSurface(state.geoContext, position(state),
- direction(state), true);
+ /// @param state [in,out] The stepping state (thread-local cache)
+ std::string outputStepSize(const State& state) const {
+ return state.stepSize.toString();
+ }
+
+ /// Overstep limit
+ ///
+ /// @param state [in] The stepping state (thread-local cache)
+ double overstepLimit(const State& /*state*/) const {
+ // A dynamic overstep limit could sit here
+ return -m_overstepLimit;
}
/// Create and return the bound state at the current position
@@ -319,7 +375,7 @@ class EigenStepper {
BField m_bField;
/// Overstep limit: could/should be dynamic
- double m_overstepLimit = -50_um;
+ double m_overstepLimit = 100_um;
};
} // namespace Acts
diff --git a/Core/include/Acts/Propagator/EigenStepper.ipp b/Core/include/Acts/Propagator/EigenStepper.ipp
index e30f2f8e25a768a8427ab6e2429c0612082730f6..109160aa92a0c7e654d7f4893fd3f7ccc64f5575 100644
--- a/Core/include/Acts/Propagator/EigenStepper.ipp
+++ b/Core/include/Acts/Propagator/EigenStepper.ipp
@@ -260,9 +260,7 @@ Acts::Result<double> Acts::EigenStepper<B, E, A>::step(
return (error_estimate <= state.options.tolerance);
};
- double stepSizeScaling;
-
- unsigned int stepAttempts = 0;
+ double stepSizeScaling = 1.;
// Select and adjust the appropriate Runge-Kutta step size as given
// ATL-SOFT-PUB-2009-001
while (!tryRungeKuttaStep(state.stepping.stepSize)) {
@@ -276,14 +274,6 @@ Acts::Result<double> Acts::EigenStepper<B, E, A>::step(
}
state.stepping.stepSize = state.stepping.stepSize * stepSizeScaling;
- // @HOTFIX to break the RK step trying, @TODO will be replaced
- // by proper overstepping mechanism
- if (++stepAttempts == 100) {
- // step in mm steps, costly but should do
- state.stepping.stepSize = state.stepping.navDir * 1_mm;
- break;
- }
-
// If step size becomes too small the particle remains at the initial
// place
if (state.stepping.stepSize * state.stepping.stepSize <
diff --git a/Core/include/Acts/Propagator/Navigator.hpp b/Core/include/Acts/Propagator/Navigator.hpp
index 5e9352c5d9f3b632a4c30a5eda4179f3656c46cd..de498b9e43084eeab69d38db54c7fa485d586f99 100644
--- a/Core/include/Acts/Propagator/Navigator.hpp
+++ b/Core/include/Acts/Propagator/Navigator.hpp
@@ -17,14 +17,13 @@
#include "Acts/Geometry/Layer.hpp"
#include "Acts/Geometry/TrackingGeometry.hpp"
#include "Acts/Geometry/TrackingVolume.hpp"
+#include "Acts/Propagator/ConstrainedStep.hpp"
#include "Acts/Propagator/Propagator.hpp"
-#include "Acts/Propagator/detail/ConstrainedStep.hpp"
#include "Acts/Surfaces/Surface.hpp"
#include "Acts/Utilities/Units.hpp"
namespace Acts {
-using Cstep = detail::ConstrainedStep;
using namespace Acts::UnitLiterals;
/// @brief struct for the Navigation options that are forwarded to
@@ -97,8 +96,8 @@ struct NavigationOptions {
/// This navigation actor thus always needs to run first!
/// It does two things: it figures out the order of volumes, layers and
/// surfaces. For each propagation step, the operator() runs, which checks if
-/// the current surface (or layer/volume boundary) is reached via
-/// surfaceReached.
+/// the current surface (or layer/volume boundary) is reached.
+///
/// The current target surface is the surface pointed to by of the iterators
/// for the surfaces, layers or volume boundaries.
/// If a surface is found, the state.navigation.currentSurface
@@ -106,7 +105,7 @@ struct NavigationOptions {
/// actor uses the ordered iterators to figure out which surface, layer or
/// volume boundary is _supposed_ to be hit next. It then sets the maximum
/// step size to the path length found out by straight line intersection. If
-/// the surfaceReached call fails, it also re-computes the step size, to make
+/// the state is not on surface, it also re-computes the step size, to make
/// sure we end up at the desired surface.
///
class Navigator {
@@ -325,7 +324,7 @@ class Navigator {
});
// Navigation break & release navigation stepping
state.navigation.navigationBreak = true;
- state.stepping.stepSize.release(Cstep::actor);
+ stepper.releaseStepSize(state.stepping);
return;
} else {
debugLog(state, [&] { return std::string("Volume updated."); });
@@ -345,7 +344,7 @@ class Navigator {
});
// Set navigation break and release the navigation step size
state.navigation.navigationBreak = true;
- state.stepping.stepSize.release(Cstep::actor);
+ stepper.releaseStepSize(state.stepping);
} else {
debugLog(state, [&] {
return std::string("Status could not be determined - good luck.");
@@ -399,7 +398,7 @@ class Navigator {
});
// Set navigation break and release the navigation step size
state.navigation.navigationBreak = true;
- state.stepping.stepSize.release(Cstep::actor);
+ stepper.releaseStepSize(state.stepping);
}
// Navigator target always resets the current surface
@@ -518,7 +517,9 @@ class Navigator {
// Check if we are at a surface
// If we are on the surface pointed at by the iterator, we can make
// it the current one to pass it to the other actors
- if (stepper.surfaceReached(state.stepping, surface)) {
+ auto surfaceStatus =
+ stepper.updateSurfaceStatus(state.stepping, *surface, true);
+ if (surfaceStatus == Intersection::Status::onSurface) {
debugLog(state, [&] {
return std::string("Status Surface successfully hit, storing it.");
});
@@ -572,7 +573,7 @@ class Navigator {
});
// set the navigation break
state.navigation.navigationBreak = true;
- state.stepping.stepSize.release(Cstep::actor);
+ stepper.releaseStepSize(state.stepping);
}
return startResolved;
}
@@ -586,13 +587,7 @@ class Navigator {
});
return false;
}
-
- // Create the navigaton options
- NavigationOptions<Surface> navOpts(state.stepping.navDir, true);
- navOpts.pathLimit = state.stepping.stepSize.value(Cstep::aborter);
- navOpts.overstepLimit = stepper.overstepLimit(state.stepping);
-
- // Loop over the navigation surfaces
+ // Loop over the remaining navigation surfaces
while (state.navigation.navSurfaceIter !=
state.navigation.navSurfaces.end()) {
// Screen output how much is left to try
@@ -613,26 +608,20 @@ class Navigator {
dstream << surface->geoID();
return dstream.str();
});
- // Now intersect (should exclude punch-through)
- auto surfaceIntersect = surface->intersect(
- state.geoContext, stepper.position(state.stepping),
- stepper.direction(state.stepping), navOpts.boundaryCheck);
- double surfaceDistance = surfaceIntersect.intersection.pathLength;
- if (!surfaceIntersect) {
+ // Estimate the surface status
+ auto surfaceStatus =
+ stepper.updateSurfaceStatus(state.stepping, *surface, true);
+ if (surfaceStatus == Intersection::Status::reachable) {
debugLog(state, [&] {
std::stringstream dstream;
- dstream << "Surface intersection at path length ";
- dstream << surfaceDistance;
- dstream << " is not valid.";
+ dstream << "Surface reachable, step size updated to ";
+ dstream << stepper.outputStepSize(state.stepping);
return dstream.str();
});
- ++state.navigation.navSurfaceIter;
- continue;
+ return true;
}
- // Update the step for the stepper
- updateStep(state, surfaceDistance, true);
- // Return to the propagator
- return true;
+ ++state.navigation.navSurfaceIter;
+ continue;
}
// Reached the end of the surface iteration
@@ -730,7 +719,7 @@ class Navigator {
if (!protoNavSurfaces.empty()) {
// did we find any surfaces?
- // are we on the first surface?
+ // Check: are we on the first surface?
if (state.navigation.currentSurface == nullptr ||
state.navigation.currentSurface !=
protoNavSurfaces.front().object) {
@@ -741,9 +730,15 @@ class Navigator {
state.navigation.navSurfaces.begin();
state.navigation.navLayers = {};
state.navigation.navLayerIter = state.navigation.navLayers.end();
- updateStep(state,
- state.navigation.navSurfaceIter->intersection.pathLength,
- true);
+ // The stepper updates the step size ( single / multi component)
+ stepper.updateStepSize(state.stepping,
+ *state.navigation.navSurfaceIter, true);
+ debugLog(state, [&] {
+ std::stringstream dstream;
+ dstream << "Navigation stepSize updated to ";
+ dstream << stepper.outputStepSize(state.stepping);
+ return dstream.str();
+ });
return true;
}
}
@@ -772,24 +767,24 @@ class Navigator {
continue;
}
}
- // Otherwise try to step towards it
- NavigationOptions<Surface> navOpts(state.stepping.navDir, true);
- navOpts.overstepLimit = stepper.overstepLimit(state.stepping);
- auto layerIntersect = layerSurface->intersect(
- state.geoContext, stepper.position(state.stepping),
- stepper.direction(state.stepping), navOpts.boundaryCheck);
- // check if the intersect is invalid
- if (!layerIntersect) {
+ // Try to step towards it
+ auto layerStatus =
+ stepper.updateSurfaceStatus(state.stepping, *layerSurface, true);
+ if (layerStatus == Intersection::Status::reachable) {
debugLog(state, [&] {
- return std::string("Layer intersection not valid, skipping it.");
+ std::stringstream dstream;
+ dstream << "Layer reachable, step size updated to ";
+ dstream << stepper.outputStepSize(state.stepping);
+ return dstream.str();
});
- ++state.navigation.navLayerIter;
- } else {
- // update the navigation step size, release the former first
- updateStep(state, layerIntersect.intersection.pathLength, true);
return true;
}
+ debugLog(state, [&] {
+ return std::string("Layer intersection not valid, skipping it.");
+ });
+ ++state.navigation.navLayerIter;
}
+
// Re-initialize target at last layer, only in case it is the target volume
// This avoids a wrong target volume estimation
if (state.navigation.currentVolume == state.navigation.targetVolume) {
@@ -851,7 +846,7 @@ class Navigator {
"No sufficient information to resolve boundary, "
"stopping navigation.");
});
- state.stepping.stepSize.release(Cstep::actor);
+ stepper.releaseStepSize(state.stepping);
return false;
} else if (state.navigation.currentVolume ==
state.navigation.targetVolume) {
@@ -860,7 +855,7 @@ class Navigator {
"In target volume: no need to resolve boundary, "
"stopping navigation.");
state.navigation.navigationBreak = true;
- state.stepping.stepSize.release(Cstep::actor);
+ stepper.releaseStepSize(state.stepping);
});
return true;
}
@@ -870,7 +865,7 @@ class Navigator {
// The navigation options
NavigationOptions<Surface> navOpts(state.stepping.navDir, true);
- navOpts.pathLimit = state.stepping.stepSize.value(Cstep::aborter);
+ navOpts.pathLimit = state.stepping.stepSize.value(ConstrainedStep::aborter);
navOpts.overstepLimit = stepper.overstepLimit(state.stepping);
// If you came until here, and you might not have boundaries
@@ -902,6 +897,18 @@ class Navigator {
});
// Set the begin iterator
state.navigation.navBoundaryIter = state.navigation.navBoundaries.begin();
+ if (not state.navigation.navBoundaries.empty()) {
+ // Set to the first and return to the stepper
+ stepper.updateStepSize(state.stepping,
+ *state.navigation.navBoundaryIter, true);
+ debugLog(state, [&] {
+ std::stringstream dstream;
+ dstream << "Navigation stepSize updated to ";
+ dstream << stepper.outputStepSize(state.stepping);
+ return dstream.str();
+ });
+ return true;
+ }
}
// Loop over the boundary surface
@@ -909,38 +916,20 @@ class Navigator {
state.navigation.navBoundaries.end()) {
// That is the current boundary surface
auto boundarySurface = state.navigation.navBoundaryIter->representation;
- // Step towards the boundary surface
- auto boundaryIntersect = boundarySurface->intersect(
- state.geoContext, stepper.position(state.stepping),
- stepper.direction(state.stepping), navOpts.boundaryCheck);
- // Distance
- auto distance = boundaryIntersect.intersection.pathLength;
- // Check the boundary is properly intersected: we are in target mode
- if (!boundaryIntersect or
- distance * distance < s_onSurfaceTolerance * s_onSurfaceTolerance) {
+ // Step towards the boundary surfrace
+ auto boundaryStatus = stepper.updateSurfaceStatus(
+ state.stepping, *boundarySurface, navOpts.boundaryCheck);
+ if (boundaryStatus == Intersection::Status::reachable) {
debugLog(state, [&] {
- std::stringstream ss;
- ss << "Boundary intersection with:\n";
- boundaryIntersect.object->toStream(state.geoContext, ss);
- ss << "\n";
- ss << "Boundary intersection not valid, skipping it.\n";
- ss << "valid: " << bool(boundaryIntersect) << "\n";
- ss << "pathLength: " << distance << "\n";
- if (distance < 0 && std::abs(distance) < 0.01) {
- ss << "Very likely overstepped over boundary surface! \n";
- }
- return ss.str();
+ std::stringstream dstream;
+ dstream << "Boundary reachable, step size updated to ";
+ dstream << stepper.outputStepSize(state.stepping);
+ return dstream.str();
});
- // Increase the iterator to the next one
- ++state.navigation.navBoundaryIter;
-
- } else {
- debugLog(state,
- [&] { return std::string("Boundary intersection valid."); });
- // This is a new navigation stream, release the former first
- updateStep(state, distance, true);
return true;
}
+ // Increase the iterator to the next one
+ ++state.navigation.navBoundaryIter;
}
// Could not do anything
return false;
@@ -969,7 +958,8 @@ class Navigator {
template <typename propagator_state_t, typename stepper_t>
void initializeTarget(propagator_state_t& state,
const stepper_t& stepper) const {
- if (state.navigation.targetVolume && state.stepping.pathAccumulated == 0.) {
+ if (state.navigation.targetVolume and
+ state.stepping.pathAccumulated == 0.) {
debugLog(state, [&] {
return std::string("Re-initialzing cancelled as it is the first step.");
});
@@ -1002,39 +992,35 @@ class Navigator {
}
// Slow navigation initialization for target:
// target volume and layer search through global search
- NavigationOptions<Surface> navOpts(state.stepping.navDir, false,
- resolveSensitive, resolveMaterial,
- resolvePassive);
- navOpts.overstepLimit = stepper.overstepLimit(state.stepping);
-
- // take the target intersection
auto targetIntersection = state.navigation.targetSurface->intersect(
state.geoContext, stepper.position(state.stepping),
- stepper.direction(state.stepping), navOpts.boundaryCheck);
- debugLog(state, [&] {
- std::stringstream dstream;
- dstream << "Target estimate position (";
- dstream << targetIntersection.intersection.position.x() << ", ";
- dstream << targetIntersection.intersection.position.y() << ", ";
- dstream << targetIntersection.intersection.position.z() << ")";
- return dstream.str();
- });
- /// get the target volume from the intersection
- auto tPosition = targetIntersection.intersection.position;
- state.navigation.targetVolume =
- trackingGeometry->lowestTrackingVolume(state.geoContext, tPosition);
- state.navigation.targetLayer =
- state.navigation.targetVolume
- ? state.navigation.targetVolume->associatedLayer(state.geoContext,
- tPosition)
- : nullptr;
- if (state.navigation.targetVolume) {
+ state.stepping.navDir * stepper.direction(state.stepping), false);
+ if (targetIntersection) {
debugLog(state, [&] {
std::stringstream dstream;
- dstream << "Target volume estimated : ";
- dstream << state.navigation.targetVolume->volumeName();
+ dstream << "Target estimate position (";
+ dstream << targetIntersection.intersection.position.x() << ", ";
+ dstream << targetIntersection.intersection.position.y() << ", ";
+ dstream << targetIntersection.intersection.position.z() << ")";
return dstream.str();
});
+ /// get the target volume from the intersection
+ auto tPosition = targetIntersection.intersection.position;
+ state.navigation.targetVolume =
+ trackingGeometry->lowestTrackingVolume(state.geoContext, tPosition);
+ state.navigation.targetLayer =
+ state.navigation.targetVolume
+ ? state.navigation.targetVolume->associatedLayer(
+ state.geoContext, tPosition)
+ : nullptr;
+ if (state.navigation.targetVolume) {
+ debugLog(state, [&] {
+ std::stringstream dstream;
+ dstream << "Target volume estimated : ";
+ dstream << state.navigation.targetVolume->volumeName();
+ return dstream.str();
+ });
+ }
}
}
}
@@ -1064,7 +1050,7 @@ class Navigator {
state.stepping.navDir, true, resolveSensitive, resolveMaterial,
resolvePassive, startSurface, state.navigation.targetSurface);
// Check the limit
- navOpts.pathLimit = state.stepping.stepSize.value(Cstep::aborter);
+ navOpts.pathLimit = state.stepping.stepSize.value(ConstrainedStep::aborter);
// No overstepping on start layer, otherwise ask the stepper
navOpts.overstepLimit = (cLayer != nullptr)
? s_onSurfaceTolerance
@@ -1087,11 +1073,15 @@ class Navigator {
});
// set the iterator
state.navigation.navSurfaceIter = state.navigation.navSurfaces.begin();
- // Update the navigation step size before you return to the stepper
- // This is a new navigation stream, release the former step size first
- updateStep(state,
- state.navigation.navSurfaceIter->intersection.pathLength,
- true);
+ // The stepper updates the step size ( single / multi component)
+ stepper.updateStepSize(state.stepping, *state.navigation.navSurfaceIter,
+ true);
+ debugLog(state, [&] {
+ std::stringstream dstream;
+ dstream << "Navigation stepSize updated to ";
+ dstream << stepper.outputStepSize(state.stepping);
+ return dstream.str();
+ });
return true;
}
state.navigation.navSurfaceIter = state.navigation.navSurfaces.end();
@@ -1133,7 +1123,7 @@ class Navigator {
resolvePassive, startLayer, nullptr);
// Set also the target surface
navOpts.targetSurface = state.navigation.targetSurface;
- navOpts.pathLimit = state.stepping.stepSize.value(Cstep::aborter);
+ navOpts.pathLimit = state.stepping.stepSize.value(ConstrainedStep::aborter);
navOpts.overstepLimit = stepper.overstepLimit(state.stepping);
// Request the compatible layers
state.navigation.navLayers =
@@ -1160,11 +1150,15 @@ class Navigator {
state.navigation.navLayerIter->object !=
state.navigation.startLayer) {
debugLog(state, [&] { return std::string("Target at layer."); });
- // update the navigation step size before you return
- updateStep(state,
- state.navigation.navLayerIter->intersection.pathLength,
- true);
- // Trigger the return to the propagator
+ // The stepper updates the step size ( single / multi component)
+ stepper.updateStepSize(state.stepping, *state.navigation.navLayerIter,
+ true);
+ debugLog(state, [&] {
+ std::stringstream dstream;
+ dstream << "Navigation stepSize updated to ";
+ dstream << stepper.outputStepSize(state.stepping);
+ return dstream.str();
+ });
return true;
}
}
@@ -1176,33 +1170,11 @@ class Navigator {
debugLog(state, [&] {
return std::string("No compatible layer candidates found.");
});
- // Update the navigation step to the target step to trigger
- // step modification when requested
- updateStep(state, state.stepping.stepSize.value(Cstep::aborter));
-
+ // Release the step size
+ stepper.releaseStepSize(state.stepping);
return false;
}
- /// This method updates the constrained step size
- ///
- /// @tparam propagator_state_t is the state type
- ///
- /// @param[in,out] state The state object for the step length
- /// @param[in] step the step size
- /// @param[in] release flag steers if the step is released first
- template <typename propagator_state_t>
- void updateStep(propagator_state_t& state, double navigationStep,
- bool release = false) const { // update the step
- state.stepping.stepSize.update(navigationStep, Cstep::actor, release);
- debugLog(state, [&] {
- std::stringstream dstream;
- std::string releaseFlag = release ? "released and " : "";
- dstream << "Navigation stepSize " << releaseFlag << "updated to ";
- dstream << state.stepping.stepSize.toString();
- return dstream.str();
- });
- }
-
/// --------------------------------------------------------------------
/// Inactive
///
@@ -1238,12 +1210,12 @@ class Navigator {
if (state.navigation.targetReached || !state.navigation.targetSurface) {
return true;
}
+ auto targetStatus = stepper.updateSurfaceStatus(
+ state.stepping, *state.navigation.targetSurface, true);
// the only advance could have been to the target
- if (stepper.surfaceReached(state.stepping,
- state.navigation.targetSurface)) {
+ if (targetStatus == Intersection::Status::onSurface) {
// set the target surface
state.navigation.currentSurface = state.navigation.targetSurface;
-
debugLog(state, [&] {
std::stringstream dstream;
dstream << "Current surface set to target surface ";
diff --git a/Core/include/Acts/Propagator/Propagator.hpp b/Core/include/Acts/Propagator/Propagator.hpp
index 57822b9a9101f0b94b5af5530d0f057b418a42fd..eb55c63bf94b3170321791c5c66f01c1edf863da 100644
--- a/Core/include/Acts/Propagator/Propagator.hpp
+++ b/Core/include/Acts/Propagator/Propagator.hpp
@@ -94,9 +94,9 @@ struct PropagatorOptions {
eoptions.absPdgCode = absPdgCode;
eoptions.mass = mass;
eoptions.maxSteps = maxSteps;
- eoptions.maxStepSize = maxStepSize;
+ eoptions.maxStepSize = direction * std::abs(maxStepSize);
eoptions.targetTolerance = targetTolerance;
- eoptions.pathLimit = pathLimit;
+ eoptions.pathLimit = direction * std::abs(pathLimit);
eoptions.loopProtection = loopProtection;
eoptions.loopFraction = loopFraction;
// Output option
@@ -241,7 +241,7 @@ class Propagator final {
State(const parameters_t& start, const propagator_options_t& topts)
: options(topts),
stepping(topts.geoContext, topts.magFieldContext, start,
- topts.direction, topts.maxStepSize),
+ topts.direction, topts.maxStepSize, topts.tolerance),
geoContext(topts.geoContext) {
// Setting the start surface
navigation.startSurface = &start.referenceSurface();
diff --git a/Core/include/Acts/Propagator/Propagator.ipp b/Core/include/Acts/Propagator/Propagator.ipp
index b4c684c4d686ef9ca58763108a3471a7c0d31c32..f5f188081e57b4f83fbc3ced3fa1808adf3144f0 100644
--- a/Core/include/Acts/Propagator/Propagator.ipp
+++ b/Core/include/Acts/Propagator/Propagator.ipp
@@ -107,6 +107,8 @@ auto Acts::Propagator<S, N>::propagate(
// Expand the abort list with a path aborter
path_aborter_t pathAborter;
+ pathAborter.internalLimit = options.pathLimit;
+
auto abortList = options.abortList.append(pathAborter);
// The expanded options (including path limit)
@@ -166,6 +168,7 @@ auto Acts::Propagator<S, N>::propagate(
// Type of provided options
target_aborter_t targetAborter;
path_aborter_t pathAborter;
+ pathAborter.internalLimit = options.pathLimit;
auto abortList = options.abortList.append(targetAborter, pathAborter);
// Create the extended options and declare their type
diff --git a/Core/include/Acts/Propagator/StepperConcept.hpp b/Core/include/Acts/Propagator/StepperConcept.hpp
index 585d33281f4c857d0ed531a95743f4a6a9cd8ae9..7dc850663c11be5dd8e809fd18cf2a58984303d7 100644
--- a/Core/include/Acts/Propagator/StepperConcept.hpp
+++ b/Core/include/Acts/Propagator/StepperConcept.hpp
@@ -9,8 +9,11 @@
#pragma once
#include "Acts/EventData/TrackParameters.hpp"
-#include "Acts/Propagator/detail/ConstrainedStep.hpp"
+#include "Acts/Propagator/ConstrainedStep.hpp"
+#include "Acts/Surfaces/BoundaryCheck.hpp"
+#include "Acts/Surfaces/Surface.hpp"
#include "Acts/Utilities/Definitions.hpp"
+#include "Acts/Utilities/Intersection.hpp"
#include "Acts/Utilities/TypeTraits.hpp"
namespace Acts {
@@ -40,12 +43,15 @@ namespace concept {
METHOD_TRAIT(charge_t, charge);
METHOD_TRAIT(time_t, time);
METHOD_TRAIT(overstep_t, overstepLimit);
- METHOD_TRAIT(surface_reached_t, surfaceReached);
METHOD_TRAIT(bound_state_method_t, boundState);
METHOD_TRAIT(curvilinear_state_method_t, curvilinearState);
METHOD_TRAIT(update_t, update);
METHOD_TRAIT(covariance_transport_t, covarianceTransport);
METHOD_TRAIT(step_t, step);
+ METHOD_TRAIT(update_surface_status_t, updateSurfaceStatus);
+ METHOD_TRAIT(set_step_size_t, setStepSize);
+ METHOD_TRAIT(release_step_size_t, releaseStepSize);
+ METHOD_TRAIT(output_step_size_t, outputStepSize);
template <typename T>
using cov_transport_t = decltype(std::declval<T>().covTransport);
@@ -65,7 +71,7 @@ namespace concept {
has_member<S, cov_t, BoundSymMatrix>,
has_member<S, nav_dir_t, NavigationDirection>,
has_member<S, path_accumulated_t, double>,
- has_member<S, step_size_t, detail::ConstrainedStep>
+ has_member<S, step_size_t, ConstrainedStep>
>;
// clang-format on
@@ -98,8 +104,6 @@ namespace concept {
static_assert(time_exists, "time method not found");
constexpr static bool overstep_exists = has_method<const S, double, overstep_t, const state&>;
static_assert(overstep_exists, "overstepLimit method not found");
- constexpr static bool surface_reached_exists = has_method<const S, bool, surface_reached_t, const state&, const Surface*>;
- static_assert(surface_reached_exists, "surfaceReached method not found");
constexpr static bool bound_state_method_exists= has_method<const S, typename S::BoundState, bound_state_method_t, state&, const Surface&, bool>;
static_assert(bound_state_method_exists, "boundState method not found");
constexpr static bool curvilinear_state_method_exists = has_method<const S, typename S::CurvilinearState, curvilinear_state_method_t, state&, bool>;
@@ -110,6 +114,14 @@ namespace concept {
constexpr static bool covariance_transport_exists = require<has_method<const S, void, covariance_transport_t, state&, bool>,
has_method<const S, void, covariance_transport_t, state&, const Surface&, bool>>;
static_assert(covariance_transport_exists, "covarianceTransport method not found");
+ constexpr static bool update_surface_exists = has_method<const S, Intersection::Status, update_surface_status_t, state&, const Surface&, const BoundaryCheck&>;
+ static_assert(update_surface_exists, "updateSurfaceStatus method not found");
+ constexpr static bool set_step_size_exists = has_method<const S, void, set_step_size_t, state&, double, ConstrainedStep::Type>;
+ static_assert(set_step_size_exists, "setStepSize method not found");
+ constexpr static bool release_step_size_exists = has_method<const S, void, release_step_size_t, state&>;
+ static_assert(release_step_size_exists, "releaseStepSize method not found");
+ constexpr static bool output_step_size_exists = has_method<const S, std::string, output_step_size_t, const state&>;
+ static_assert(output_step_size_exists, "outputStepSize method not found");
constexpr static bool value = require<state_exists,
jacobian_exists,
@@ -123,11 +135,14 @@ namespace concept {
momentum_exists,
charge_exists,
time_exists,
- surface_reached_exists,
bound_state_method_exists,
curvilinear_state_method_exists,
update_method_exists,
- covariance_transport_exists>;
+ covariance_transport_exists,
+ update_surface_exists,
+ set_step_size_exists,
+ release_step_size_exists,
+ output_step_size_exists>;
};
// clang-format on
} // namespace Stepper
diff --git a/Core/include/Acts/Propagator/StraightLineStepper.hpp b/Core/include/Acts/Propagator/StraightLineStepper.hpp
index df6cfa061d6a247b8a17749e5223f2bbf29c8062..df07a82071cee025051778f589f78da867c21f0b 100644
--- a/Core/include/Acts/Propagator/StraightLineStepper.hpp
+++ b/Core/include/Acts/Propagator/StraightLineStepper.hpp
@@ -15,7 +15,8 @@
#include "Acts/Geometry/GeometryContext.hpp"
#include "Acts/MagneticField/MagneticFieldContext.hpp"
#include "Acts/MagneticField/NullBField.hpp"
-#include "Acts/Propagator/detail/ConstrainedStep.hpp"
+#include "Acts/Propagator/ConstrainedStep.hpp"
+#include "Acts/Propagator/detail/SteppingHelper.hpp"
#include "Acts/Surfaces/Surface.hpp"
#include "Acts/Utilities/Definitions.hpp"
#include "Acts/Utilities/Intersection.hpp"
@@ -30,8 +31,6 @@ namespace Acts {
/// used for simple material mapping, navigation validation
class StraightLineStepper {
public:
- using cstep = detail::ConstrainedStep;
-
using Jacobian = BoundMatrix;
using Covariance = BoundSymMatrix;
using BoundState = std::tuple<BoundParameters, Jacobian, double>;
@@ -53,11 +52,13 @@ class StraightLineStepper {
/// @param [in] par The track parameters at start
/// @param [in] ndir is the navigation direction
/// @param [in] ssize is the (absolute) maximum step size
+ /// @param [in] stolerance is the stepping tolerance
template <typename parameters_t>
explicit State(std::reference_wrapper<const GeometryContext> gctx,
std::reference_wrapper<const MagneticFieldContext> /*mctx*/,
const parameters_t& par, NavigationDirection ndir = forward,
- double ssize = std::numeric_limits<double>::max())
+ double ssize = std::numeric_limits<double>::max(),
+ double stolerance = s_onSurfaceTolerance)
: pos(par.position()),
dir(par.momentum().normalized()),
p(par.momentum().norm()),
@@ -65,6 +66,7 @@ class StraightLineStepper {
t0(par.time()),
navDir(ndir),
stepSize(ndir * std::abs(ssize)),
+ tolerance(stolerance),
geoContext(gctx) {
if (par.covariance()) {
// Get the reference surface for navigation
@@ -119,7 +121,13 @@ class StraightLineStepper {
double pathAccumulated = 0.;
/// adaptive step size of the runge-kutta integration
- cstep stepSize = std::numeric_limits<double>::max();
+ ConstrainedStep stepSize = std::numeric_limits<double>::max();
+
+ // Previous step size for overstep estimation (ignored for SL stepper)
+ double previousStepSize = 0.;
+
+ /// The tolerance for the stepping
+ double tolerance = s_onSurfaceTolerance;
// Cache the geometry context of this propagation
std::reference_wrapper<const GeometryContext> geoContext;
@@ -174,15 +182,60 @@ class StraightLineStepper {
return s_onSurfaceTolerance;
}
- /// Tests if the state reached a surface
+ /// Update surface status
+ ///
+ /// This method intersect the provided surface and update the navigation
+ /// step estimation accordingly (hence it changes the state). It also
+ /// returns the status of the intersection to trigger onSurface in case
+ /// the surface is reached.
+ ///
+ /// @param state [in,out] The stepping state (thread-local cache)
+ /// @param surface [in] The surface provided
+ /// @param bcheck [in] The boundary check for this status update
+ Intersection::Status updateSurfaceStatus(State& state, const Surface& surface,
+ const BoundaryCheck& bcheck) const {
+ return detail::updateSingleSurfaceStatus<StraightLineStepper>(
+ *this, state, surface, bcheck);
+ }
+
+ /// Update step size
+ ///
+ /// It checks the status to the reference surface & updates
+ /// the step size accordingly
+ ///
+ /// @param state [in,out] The stepping state (thread-local cache)
+ /// @param oIntersection [in] The ObjectIntersection to layer, boundary, etc
+ /// @param release [in] boolean to trigger step size release
+ template <typename object_intersection_t>
+ void updateStepSize(State& state, const object_intersection_t& oIntersection,
+ bool release = true) const {
+ detail::updateSingleStepSize<StraightLineStepper>(state, oIntersection,
+ release);
+ }
+
+ /// Set Step size - explicitely with a double
+ ///
+ /// @param state [in,out] The stepping state (thread-local cache)
+ /// @param stepSize [in] The step size value
+ /// @param stype [in] The step size type to be set
+ void setStepSize(State& state, double stepSize,
+ ConstrainedStep::Type stype = ConstrainedStep::actor) const {
+ state.previousStepSize = state.stepSize;
+ state.stepSize.update(stepSize, stype, true);
+ }
+
+ /// Release the Step size
///
- /// @param [in] state State that is tests
- /// @param [in] surface Surface that is tested
+ /// @param state [in,out] The stepping state (thread-local cache)
+ void releaseStepSize(State& state) const {
+ state.stepSize.release(ConstrainedStep::actor);
+ }
+
+ /// Output the Step Size - single component
///
- /// @return Boolean statement if surface is reached by state
- bool surfaceReached(const State& state, const Surface* surface) const {
- return surface->isOnSurface(state.geoContext, position(state),
- direction(state), true);
+ /// @param state [in,out] The stepping state (thread-local cache)
+ std::string outputStepSize(const State& state) const {
+ return state.stepSize.toString();
}
/// Create and return the bound state at the current position
diff --git a/Core/include/Acts/Propagator/detail/LoopProtection.hpp b/Core/include/Acts/Propagator/detail/LoopProtection.hpp
index 58d92d1483c56589a11cc731f73377b36388a049..53bade8e34f62836a8902b0891d1e041549c95bd 100644
--- a/Core/include/Acts/Propagator/detail/LoopProtection.hpp
+++ b/Core/include/Acts/Propagator/detail/LoopProtection.hpp
@@ -9,6 +9,7 @@
#pragma once
#include "Acts/Utilities/Definitions.hpp"
+#include "Acts/Utilities/Helpers.hpp"
namespace Acts {
namespace detail {
@@ -22,20 +23,65 @@ struct LoopProtection {
///
/// @param [in,out] state State object provided for the call
/// @param [in] stepper Stepper used
- template <typename state_t, typename stepper_t>
- void operator()(state_t& state, const stepper_t& stepper) const {
+ template <typename propagator_state_t, typename stepper_t>
+ void operator()(propagator_state_t& state, const stepper_t& stepper) const {
// Estimate the loop protection limit
if (state.options.loopProtection) {
// Get the field at the start position
Vector3D field =
stepper.getField(state.stepping, stepper.position(state.stepping));
- const double p = stepper.momentum(state.stepping);
const double B = field.norm();
- const double helixPath = 2 * M_PI * p / B;
- // now set the new loop limit
- auto& pathAborter =
- state.options.abortList.template get<path_arborter_t>();
- pathAborter.internalLimit = state.options.loopFraction * helixPath;
+ if (B != 0.) {
+ // Transverse component at start is taken for the loop protection
+ const double p = stepper.momentum(state.stepping);
+ // Calculate the full helix path
+ const double helixPath = state.stepping.navDir * 2 * M_PI * p / B;
+ // And set it as the loop limit if it overwrites the internal limit
+ auto& pathAborter =
+ state.options.abortList.template get<path_arborter_t>();
+ double loopLimit = state.options.loopFraction * helixPath;
+ double pathLimit = pathAborter.internalLimit;
+ if (loopLimit * loopLimit < pathLimit * pathLimit) {
+ pathAborter.internalLimit = loopLimit;
+ debugLog(state, [&] {
+ std::stringstream dstream;
+ dstream << "Path aborter limit set to ";
+ dstream << loopLimit << " (full helix = " << helixPath << ")";
+ return dstream.str();
+ });
+ }
+ }
+ }
+ }
+
+ /// The private loop protection debug logging
+ ///
+ /// It needs to be fed by a lambda function that returns a string,
+ /// that guarantees that the lambda is only called in the
+ /// state.options.debug == true
+ /// case in order not to spend time when not needed.
+ ///
+ /// @tparam propagator_state_t Type of the propagator state
+ ///
+ /// @param[in,out] state the propagator state for the debug flag,
+ /// prefix and length
+ /// @param logAction is a callable function that returns a streamable object
+ template <typename propagator_state_t>
+ void debugLog(propagator_state_t& state,
+ const std::function<std::string()>& logAction) const {
+ if (state.options.debug) {
+ std::vector<std::string> lines;
+ std::string input = logAction();
+ boost::split(lines, input, boost::is_any_of("\n"));
+ for (const auto& line : lines) {
+ std::stringstream dstream;
+ dstream << '\n';
+ dstream << " ∞ " << std::setw(state.options.debugPfxWidth)
+ << " loop protection "
+ << " | ";
+ dstream << std::setw(state.options.debugMsgWidth) << line << '\n';
+ state.options.debugString += dstream.str();
+ }
}
}
};
diff --git a/Core/include/Acts/Propagator/detail/StandardAborters.hpp b/Core/include/Acts/Propagator/detail/StandardAborters.hpp
index 93ec9c9d85378b373abe8bc18ca441612d3ba460..0bec804b675cb925d60063caeddb21ff3ec1d704 100644
--- a/Core/include/Acts/Propagator/detail/StandardAborters.hpp
+++ b/Core/include/Acts/Propagator/detail/StandardAborters.hpp
@@ -11,7 +11,7 @@
#include <limits>
#include <sstream>
#include <string>
-#include "Acts/Propagator/detail/ConstrainedStep.hpp"
+#include "Acts/Propagator/ConstrainedStep.hpp"
#include "Acts/Surfaces/BoundaryCheck.hpp"
#include "Acts/Surfaces/Surface.hpp"
#include "Acts/Utilities/Definitions.hpp"
@@ -80,10 +80,8 @@ struct PathLimitReached {
if (state.navigation.targetReached) {
return true;
}
-
// Check if the maximum allowed step size has to be updated
- double limit = std::min(internalLimit, state.options.pathLimit);
- double distance = limit - state.stepping.pathAccumulated;
+ double distance = internalLimit - state.stepping.pathAccumulated;
double tolerance = state.options.targetTolerance;
state.stepping.stepSize.update(distance, ConstrainedStep::aborter);
bool limitReached = (distance * distance < tolerance * tolerance);
diff --git a/Core/include/Acts/Propagator/detail/SteppingHelper.hpp b/Core/include/Acts/Propagator/detail/SteppingHelper.hpp
new file mode 100644
index 0000000000000000000000000000000000000000..0356245298b3eeba07ecf32c9ed1e882f015241d
--- /dev/null
+++ b/Core/include/Acts/Propagator/detail/SteppingHelper.hpp
@@ -0,0 +1,83 @@
+// This file is part of the Acts project.
+//
+// Copyright (C) 2019 CERN for the benefit of the Acts project
+//
+// This Source Code Form is subject to the terms of the Mozilla Public
+// License, v. 2.0. If a copy of the MPL was not distributed with this
+// file, You can obtain one at http://mozilla.org/MPL/2.0/.
+
+#pragma once
+
+#include "Acts/Propagator/ConstrainedStep.hpp"
+#include "Acts/Surfaces/BoundaryCheck.hpp"
+#include "Acts/Surfaces/Surface.hpp"
+#include "Acts/Utilities/Definitions.hpp"
+#include "Acts/Utilities/Intersection.hpp"
+
+namespace Acts {
+
+namespace detail {
+
+/// Update surface status - Single component
+///
+/// This method intersect the provided surface and update the navigation
+/// step estimation accordingly (hence it changes the state). It also
+/// returns the status of the intersection to trigger onSurface in case
+/// the surface is reached.
+///
+/// @param state [in,out] The stepping state (thread-local cache)
+/// @param surface [in] The surface provided
+/// @param bcheck [in] The boundary check for this status update
+template <typename stepper_t>
+Acts::Intersection::Status updateSingleSurfaceStatus(
+ const stepper_t& stepper, typename stepper_t::State& state,
+ const Surface& surface, const BoundaryCheck& bcheck) {
+ auto sIntersection =
+ surface.intersect(state.geoContext, stepper.position(state),
+ state.navDir * stepper.direction(state), bcheck);
+
+ // The intersection is on surface already
+ if (sIntersection.intersection.status == Intersection::Status::onSurface) {
+ // Release navigation step size
+ state.stepSize.release(ConstrainedStep::actor);
+ return Intersection::Status::onSurface;
+ } else if (sIntersection.intersection or sIntersection.alternative) {
+ // Path and overstep limit checking
+ double pLimit = state.stepSize.value(ConstrainedStep::aborter);
+ double oLimit = stepper.overstepLimit(state);
+ auto checkIntersection = [&](const Intersection& intersection) -> bool {
+ double cLimit = intersection.pathLength;
+ bool accept = (cLimit > oLimit and cLimit * cLimit < pLimit * pLimit);
+ if (accept) {
+ stepper.setStepSize(state, state.navDir * cLimit);
+ }
+ return accept;
+ };
+ // If either of the two intersections are viable return reachable
+ if (checkIntersection(sIntersection.intersection) or
+ (sIntersection.alternative and
+ checkIntersection(sIntersection.alternative))) {
+ return Intersection::Status::reachable;
+ }
+ }
+ return Intersection::Status::unreachable;
+}
+
+/// Update the Step size - single component
+///
+/// It takes a (valid) object intersection from the compatibleX(...)
+/// calls in the geometry and updates the step size
+///
+/// @param state [in,out] The stepping state (thread-local cache)
+/// @param oIntersection [in] The object that yielded this step size
+/// @param release [in] A release flag
+template <typename stepper_t, typename object_intersection_t>
+void updateSingleStepSize(typename stepper_t::State& state,
+ const object_intersection_t& oIntersection,
+ bool release = true) {
+ double stepSize = oIntersection.intersection.pathLength;
+ state.stepSize.update(stepSize, ConstrainedStep::actor, release);
+}
+
+} // namespace detail
+} // namespace Acts
\ No newline at end of file
diff --git a/Core/include/Acts/Propagator/detail/SteppingLogger.hpp b/Core/include/Acts/Propagator/detail/SteppingLogger.hpp
index ccaa3bd42754447da8f7d253fd992e1da4363753..f27c5bb6d7b73838045b31edc2f6e844fb285dca 100644
--- a/Core/include/Acts/Propagator/detail/SteppingLogger.hpp
+++ b/Core/include/Acts/Propagator/detail/SteppingLogger.hpp
@@ -9,7 +9,7 @@
#pragma once
#include "Acts/Geometry/GeometryID.hpp"
-#include "Acts/Propagator/detail/ConstrainedStep.hpp"
+#include "Acts/Propagator/ConstrainedStep.hpp"
#include "Acts/Utilities/Definitions.hpp"
namespace Acts {
diff --git a/Plugins/DD4hep/CMakeLists.txt b/Plugins/DD4hep/CMakeLists.txt
index 6069faab1de999432fff17ae45730ad5c97f41a8..813b57e319dba5d9fcae6f559274de8e22ae837b 100644
--- a/Plugins/DD4hep/CMakeLists.txt
+++ b/Plugins/DD4hep/CMakeLists.txt
@@ -4,7 +4,8 @@ add_library(
src/ConvertDD4hepDetector.cpp
src/ConvertDD4hepMaterial.cpp
src/DD4hepDetectorElement.cpp
- src/DD4hepLayerBuilder.cpp)
+ src/DD4hepLayerBuilder.cpp
+ src/DD4hepVolumeBuilder.cpp)
target_include_directories(
ActsDD4hepPlugin
PUBLIC
diff --git a/Tests/Core/Geometry/TrackingVolumeTests.cpp b/Tests/Core/Geometry/TrackingVolumeTests.cpp
index 42e37dabc3f1af1692ba67089e3c6094672e2b9e..b24a244d6f96442322b3d752ca84aa713f5ab9dd 100644
--- a/Tests/Core/Geometry/TrackingVolumeTests.cpp
+++ b/Tests/Core/Geometry/TrackingVolumeTests.cpp
@@ -26,7 +26,7 @@
#include "Acts/Propagator/ActionList.hpp"
#include "Acts/Propagator/StraightLineStepper.hpp"
#include "Acts/Propagator/Propagator.hpp"
-#include "Acts/Propagator/detail/ConstrainedStep.hpp"
+#include "Acts/Propagator/ConstrainedStep.hpp"
#include "Acts/Propagator/detail/DebugOutputActor.hpp"
#include "Acts/Propagator/detail/StandardAborters.hpp"
#include "Acts/Propagator/detail/SteppingLogger.hpp"
diff --git a/Tests/Core/Propagator/AbortListTests.cpp b/Tests/Core/Propagator/AbortListTests.cpp
index 1337255e9495e093c620c5d5e0471dced698fca1..4c613c4ebf7d0b041f93c0c242023630bf599950 100644
--- a/Tests/Core/Propagator/AbortListTests.cpp
+++ b/Tests/Core/Propagator/AbortListTests.cpp
@@ -20,7 +20,7 @@
// leave blank line
#include "Acts/Propagator/AbortList.hpp"
-#include "Acts/Propagator/detail/ConstrainedStep.hpp"
+#include "Acts/Propagator/ConstrainedStep.hpp"
#include "Acts/Propagator/detail/StandardAborters.hpp"
#include "Acts/Utilities/Definitions.hpp"
#include "Acts/Utilities/Units.hpp"
@@ -43,7 +43,6 @@ using PathLimit = detail::PathLimitReached;
using EndOfWorld = detail::EndOfWorldReached;
// the constrained step class
-using cstep = detail::ConstrainedStep;
/// This is a simple cache struct to mimic the
/// Navigator state
@@ -67,11 +66,12 @@ struct PropagatorState {
// This is a simple cache struct to mimic the
// Stepper cache in the propagation
struct StepperState {
- // accummulated path length cache
+ // Accummulated path length cache
double pathAccumulated = 0.;
-
+ // Navigation direction
+ NavigationDirection navDir = forward;
// adaptive sep size of the runge-kutta integration
- cstep stepSize = std::numeric_limits<double>::max();
+ ConstrainedStep stepSize = std::numeric_limits<double>::max();
};
/// emulate the options template
@@ -84,7 +84,7 @@ struct PropagatorState {
/// Debug output
/// the string where debug messages are stored (optionally)
- bool debug = false;
+ bool debug = true;
std::string debugString = "";
/// buffer & formatting for consistent output
size_t debugPfxWidth = 30;
@@ -117,6 +117,7 @@ BOOST_AUTO_TEST_CASE(AbortListTest_PathLimit) {
Result result;
AbortList<PathLimit> abortList;
+ abortList.get<PathLimit>().internalLimit = state.options.pathLimit;
// It should not abort yet
BOOST_CHECK(!abortList(result, state, stepper));
diff --git a/Tests/Core/Propagator/ActionListTests.cpp b/Tests/Core/Propagator/ActionListTests.cpp
index c32e2c0aa0cc840e8ce838dbaf8e913fc5f22cd5..e05d8a4cd8a71c585bb5f743e37e342c0c326f58 100644
--- a/Tests/Core/Propagator/ActionListTests.cpp
+++ b/Tests/Core/Propagator/ActionListTests.cpp
@@ -38,7 +38,6 @@ namespace Acts {
namespace Test {
// the constrained step class
-using cstep = detail::ConstrainedStep;
/// This is a simple cache struct to mimic the
/// Propagator cache
@@ -50,7 +49,7 @@ struct PropagatorState {
double pathAccumulated = 0.;
// adaptive sep size of the runge-kutta integration
- cstep stepSize = std::numeric_limits<double>::max();
+ ConstrainedStep stepSize = std::numeric_limits<double>::max();
};
/// emulate the options template
diff --git a/Tests/Core/Propagator/ConstrainedStepTests.cpp b/Tests/Core/Propagator/ConstrainedStepTests.cpp
index 0158333b547a391d4d64dba48be6f2e79be42e13..e977c7a1debac59621fbd359ea3d6061242424a4 100644
--- a/Tests/Core/Propagator/ConstrainedStepTests.cpp
+++ b/Tests/Core/Propagator/ConstrainedStepTests.cpp
@@ -19,7 +19,7 @@
#include <boost/test/output_test_stream.hpp>
// leave blank line
-#include "Acts/Propagator/detail/ConstrainedStep.hpp"
+#include "Acts/Propagator/ConstrainedStep.hpp"
namespace bdata = boost::unit_test::data;
namespace tt = boost::test_tools;
@@ -31,86 +31,84 @@ namespace Test {
// This tests the implementation of the AbortList
// and the standard aborters
BOOST_AUTO_TEST_CASE(ConstrainedStepTest) {
- using cstep = detail::ConstrainedStep;
-
// forward stepping test
- cstep stepSize_p = 0.25;
+ ConstrainedStep stepSize_p = 0.25;
// All of the types should be 0.25 now
- BOOST_CHECK_EQUAL(stepSize_p.values[cstep::accuracy],
+ BOOST_CHECK_EQUAL(stepSize_p.values[ConstrainedStep::accuracy],
std::numeric_limits<double>::max());
- BOOST_CHECK_EQUAL(stepSize_p.values[cstep::actor],
+ BOOST_CHECK_EQUAL(stepSize_p.values[ConstrainedStep::actor],
std::numeric_limits<double>::max());
- BOOST_CHECK_EQUAL(stepSize_p.values[cstep::aborter],
+ BOOST_CHECK_EQUAL(stepSize_p.values[ConstrainedStep::aborter],
std::numeric_limits<double>::max());
- BOOST_CHECK_EQUAL(stepSize_p.values[cstep::user], 0.25);
+ BOOST_CHECK_EQUAL(stepSize_p.values[ConstrainedStep::user], 0.25);
// Check the cast operation to double
BOOST_CHECK_EQUAL(stepSize_p, 0.25);
// now we update the accuracy
- stepSize_p.update(0.1, cstep::accuracy);
+ stepSize_p.update(0.1, ConstrainedStep::accuracy);
BOOST_CHECK_EQUAL(stepSize_p, 0.1);
// now we update the actor to smaller
- stepSize_p.update(0.05, cstep::actor);
+ stepSize_p.update(0.05, ConstrainedStep::actor);
BOOST_CHECK_EQUAL(stepSize_p, 0.05);
// we increase the actor, but do not release the step size
- stepSize_p.update(0.15, cstep::actor, false);
+ stepSize_p.update(0.15, ConstrainedStep::actor, false);
BOOST_CHECK_EQUAL(stepSize_p, 0.05);
// we increase the actor, but now DO release the step size
// it falls back to the accuracy
- stepSize_p.update(0.15, cstep::actor, true);
+ stepSize_p.update(0.15, ConstrainedStep::actor, true);
BOOST_CHECK_EQUAL(stepSize_p, 0.1);
// now set two and update them
- stepSize_p.update(0.05, cstep::user);
- stepSize_p.update(0.03, cstep::accuracy);
+ stepSize_p.update(0.05, ConstrainedStep::user);
+ stepSize_p.update(0.03, ConstrainedStep::accuracy);
BOOST_CHECK_EQUAL(stepSize_p, 0.03);
// now we release the accuracy - to the highest available value
- stepSize_p.release(cstep::accuracy);
- BOOST_CHECK_EQUAL(stepSize_p.values[cstep::accuracy],
+ stepSize_p.release(ConstrainedStep::accuracy);
+ BOOST_CHECK_EQUAL(stepSize_p.values[ConstrainedStep::accuracy],
std::numeric_limits<double>::max());
BOOST_CHECK_EQUAL(stepSize_p, 0.05);
// backward stepping test
- cstep stepSize_n = -0.25;
+ ConstrainedStep stepSize_n = -0.25;
// All of the types should be 0.25 now
- BOOST_CHECK_EQUAL(stepSize_n.values[cstep::accuracy],
+ BOOST_CHECK_EQUAL(stepSize_n.values[ConstrainedStep::accuracy],
-std::numeric_limits<double>::max());
- BOOST_CHECK_EQUAL(stepSize_n.values[cstep::actor],
+ BOOST_CHECK_EQUAL(stepSize_n.values[ConstrainedStep::actor],
-std::numeric_limits<double>::max());
- BOOST_CHECK_EQUAL(stepSize_n.values[cstep::aborter],
+ BOOST_CHECK_EQUAL(stepSize_n.values[ConstrainedStep::aborter],
-std::numeric_limits<double>::max());
- BOOST_CHECK_EQUAL(stepSize_n.values[cstep::user], -0.25);
+ BOOST_CHECK_EQUAL(stepSize_n.values[ConstrainedStep::user], -0.25);
// Check the cast operation to double
BOOST_CHECK_EQUAL(stepSize_n, -0.25);
// now we update the accuracy
- stepSize_n.update(-0.1, cstep::accuracy);
+ stepSize_n.update(-0.1, ConstrainedStep::accuracy);
BOOST_CHECK_EQUAL(stepSize_n, -0.1);
// now we update the actor to smaller
- stepSize_n.update(-0.05, cstep::actor);
+ stepSize_n.update(-0.05, ConstrainedStep::actor);
BOOST_CHECK_EQUAL(stepSize_n, -0.05);
// we increase the actor and accuracy is smaller again, without reset
- stepSize_n.update(-0.15, cstep::actor, false);
+ stepSize_n.update(-0.15, ConstrainedStep::actor, false);
BOOST_CHECK_EQUAL(stepSize_n, -0.05);
// we increase the actor and accuracy is smaller again, with reset
- stepSize_n.update(-0.15, cstep::actor, true);
+ stepSize_n.update(-0.15, ConstrainedStep::actor, true);
BOOST_CHECK_EQUAL(stepSize_n, -0.1);
// now set two and update them
- stepSize_n.update(-0.05, cstep::user);
- stepSize_n.update(-0.03, cstep::accuracy);
+ stepSize_n.update(-0.05, ConstrainedStep::user);
+ stepSize_n.update(-0.03, ConstrainedStep::accuracy);
BOOST_CHECK_EQUAL(stepSize_n, -0.03);
// now we release the accuracy - to the highest available value
- stepSize_n.release(cstep::accuracy);
- BOOST_CHECK_EQUAL(stepSize_n.values[cstep::accuracy],
+ stepSize_n.release(ConstrainedStep::accuracy);
+ BOOST_CHECK_EQUAL(stepSize_n.values[ConstrainedStep::accuracy],
-std::numeric_limits<double>::max());
BOOST_CHECK_EQUAL(stepSize_n, -0.05);
}
diff --git a/Tests/Core/Propagator/DirectNavigatorTests.cpp b/Tests/Core/Propagator/DirectNavigatorTests.cpp
index fb2fbfd43aa8077ee602fd2871db11483e5ef2d2..71e76b2321e7635c218ee4eed66cf8d1680d68f9 100644
--- a/Tests/Core/Propagator/DirectNavigatorTests.cpp
+++ b/Tests/Core/Propagator/DirectNavigatorTests.cpp
@@ -62,10 +62,12 @@ Stepper dstepper(bField);
ReferencePropagator rpropagator(std::move(estepper), std::move(navigator));
DirectPropagator dpropagator(std::move(dstepper), std::move(dnavigator));
-const int ntests = 500;
+const int ntests = 1000;
const int skip = 0;
bool debugMode = false;
bool referenceTiming = false;
+bool oversteppingTest = false;
+double oversteppingMaxStepSize = 1_mm;
/// The actual test nethod that runs the test
/// can be used with several propagator types
@@ -113,6 +115,9 @@ void runTest(const rpropagator_t& rprop, const dpropagator_t& dprop, double pT,
using Options = PropagatorOptions<RefereceActionList, ReferenceAbortList>;
Options pOptions(tgContext, mfContext);
pOptions.debug = debugMode;
+ if (oversteppingTest) {
+ pOptions.maxStepSize = oversteppingMaxStepSize;
+ }
// Surface collector configuration
auto& sCollector = pOptions.actionList.template get<SurfaceCollector<>>();
@@ -200,7 +205,7 @@ BOOST_DATA_TEST_CASE(
test_direct_navigator,
bdata::random((bdata::seed = 20,
bdata::distribution =
- std::uniform_real_distribution<>(0.5_GeV, 10_GeV))) ^
+ std::uniform_real_distribution<>(0.15_GeV, 10_GeV))) ^
bdata::random((bdata::seed = 21,
bdata::distribution =
std::uniform_real_distribution<>(-M_PI, M_PI))) ^
diff --git a/Tests/Core/Propagator/LoopProtectionTests.cpp b/Tests/Core/Propagator/LoopProtectionTests.cpp
index 2d2cf9b787887bf388f6da79f9c70becd2f0285f..3cc6a23168b04a592f2e634aed7073e0f82d21e7 100644
--- a/Tests/Core/Propagator/LoopProtectionTests.cpp
+++ b/Tests/Core/Propagator/LoopProtectionTests.cpp
@@ -20,6 +20,7 @@
#include "Acts/Propagator/Propagator.hpp"
#include "Acts/Propagator/detail/LoopProtection.hpp"
#include "Acts/Propagator/detail/StandardAborters.hpp"
+#include "Acts/Propagator/detail/DebugOutputActor.hpp"
#include "Acts/Tests/CommonHelpers/FloatComparisons.hpp"
#include "Acts/Geometry/GeometryContext.hpp"
#include "Acts/MagneticField/MagneticFieldContext.hpp"
@@ -46,6 +47,8 @@ struct SteppingState {
Vector3D pos = Vector3D(0., 0., 0.);
Vector3D dir = Vector3D(0., 0., 1);
double p = 100_MeV;
+
+ NavigationDirection navDir;
};
/// @brief mockup of stepping state
@@ -87,6 +90,11 @@ struct Options {
bool loopProtection = true;
double loopFraction = 0.5;
+ bool debug = false;
+ std::string debugString;
+ int debugMsgWidth = 60;
+ int debugPfxWidth = 30;
+
/// Contains: target aborters
AbortList<PathLimitReached> abortList;
};
@@ -183,10 +191,20 @@ BOOST_DATA_TEST_CASE(
Vector3D mom(px, py, pz);
CurvilinearParameters start(std::nullopt, pos, mom, q, 42.);
- PropagatorOptions<> options(tgContext, mfContext);
+ using DebugOutput = Acts::detail::DebugOutputActor;
+ using ProopagatorOptions =
+ PropagatorOptions<ActionList<DebugOutput>, AbortList<>>;
+ ProopagatorOptions options(tgContext, mfContext);
+ options.debug = false;
options.maxSteps = 1e6;
const auto& result = epropagator.propagate(start, options).value();
+ if (options.debug) {
+ const auto debugString =
+ result.template get<DebugOutput::result_type>().debugString;
+ std::cout << debugString << std::endl;
+ }
+
// this test assumes state.options.loopFraction = 0.5
CHECK_CLOSE_REL(px, -result.endParameters->momentum().x(), 1e-2);
CHECK_CLOSE_REL(py, -result.endParameters->momentum().y(), 1e-2);
diff --git a/Tests/Core/Propagator/NavigatorTests.cpp b/Tests/Core/Propagator/NavigatorTests.cpp
index b8ae317f7384144dafa92315a4675d3762479b9d..51d428ecf6c0d466d6ed33440d632c7a2220a75a 100644
--- a/Tests/Core/Propagator/NavigatorTests.cpp
+++ b/Tests/Core/Propagator/NavigatorTests.cpp
@@ -20,7 +20,8 @@
#include "Acts/Propagator/Navigator.hpp"
#include "Acts/Propagator/StraightLineStepper.hpp"
#include "Acts/Propagator/StepperConcept.hpp"
-#include "Acts/Propagator/detail/ConstrainedStep.hpp"
+#include "Acts/Propagator/ConstrainedStep.hpp"
+#include "Acts/Propagator/detail/SteppingHelper.hpp"
#include "Acts/Surfaces/CylinderSurface.hpp"
#include "Acts/Tests/CommonHelpers/CylindricalTrackingGeometry.hpp"
#include "Acts/Tests/CommonHelpers/FloatComparisons.hpp"
@@ -81,7 +82,15 @@ struct PropagatorState {
double pathAccumulated = 0.;
// adaptive sep size of the runge-kutta integration
- Cstep stepSize = Cstep(100_cm);
+ ConstrainedStep stepSize = ConstrainedStep(100_cm);
+
+ // Previous step size for overstep estimation (ignored here)
+ double previousStepSize = 0.;
+
+ /// The tolerance for the stepping
+ double tolerance = s_onSurfaceTolerance;
+
+ GeometryContext geoContext = GeometryContext();
};
/// Global particle position accessor
@@ -104,9 +113,33 @@ struct PropagatorState {
return s_onSurfaceTolerance;
}
- bool surfaceReached(const State& state, const Surface* surface) const {
- return surface->isOnSurface(tgContext, position(state), direction(state),
- true);
+ Intersection::Status updateSurfaceStatus(
+ State& state, const Surface& surface,
+ const BoundaryCheck& bcheck) const {
+ return detail::updateSingleSurfaceStatus<Stepper>(*this, state, surface,
+ bcheck);
+ }
+
+ template <typename object_intersection_t>
+ void updateStepSize(State& state,
+ const object_intersection_t& oIntersection,
+ bool release = true) const {
+ detail::updateSingleStepSize<Stepper>(state, oIntersection, release);
+ }
+
+ void setStepSize(
+ State& state, double stepSize,
+ ConstrainedStep::Type stype = ConstrainedStep::actor) const {
+ state.previousStepSize = state.stepSize;
+ state.stepSize.update(stepSize, stype, true);
+ }
+
+ void releaseStepSize(State& state) const {
+ state.stepSize.release(ConstrainedStep::actor);
+ }
+
+ std::string outputStepSize(const State& state) const {
+ return state.stepSize.toString();
}
BoundState boundState(State& state, const Surface& surface,
diff --git a/Tests/Core/Propagator/PropagatorTests.cpp b/Tests/Core/Propagator/PropagatorTests.cpp
index 5c5f18945f4316a3f33959184bfc04bc78f56195..1a408bca856bf3173644b0b0a11d9d220332f403 100644
--- a/Tests/Core/Propagator/PropagatorTests.cpp
+++ b/Tests/Core/Propagator/PropagatorTests.cpp
@@ -19,7 +19,7 @@
#include "Acts/Propagator/ActionList.hpp"
#include "Acts/Propagator/EigenStepper.hpp"
#include "Acts/Propagator/Propagator.hpp"
-#include "Acts/Propagator/detail/ConstrainedStep.hpp"
+#include "Acts/Propagator/ConstrainedStep.hpp"
#include "Acts/Propagator/detail/StandardAborters.hpp"
#include "Acts/Surfaces/CylinderSurface.hpp"
#include "Acts/Tests/CommonHelpers/FloatComparisons.hpp"
@@ -40,7 +40,6 @@ namespace Test {
GeometryContext tgContext = GeometryContext();
MagneticFieldContext mfContext = MagneticFieldContext();
-using cstep = detail::ConstrainedStep;
using Covariance = BoundSymMatrix;
/// An observer that measures the perpendicular distance
@@ -95,13 +94,13 @@ struct SurfaceObserver {
stepper.direction(state.stepping), true)
.pathLength;
// Adjust the step size so that we cannot cross the target surface
- state.stepping.stepSize.update(distance, cstep::actor);
+ state.stepping.stepSize.update(distance, ConstrainedStep::actor);
// return true if you fall below tolerance
if (std::abs(distance) <= tolerance) {
++result.surfaces_passed;
result.surface_passed_r = perp(stepper.position(state.stepping));
// release the step size, will be re-adjusted
- state.stepping.stepSize.release(cstep::actor);
+ state.stepping.stepSize.release(ConstrainedStep::actor);
}
}
}
diff --git a/Tests/Core/Propagator/StepperTests.cpp b/Tests/Core/Propagator/StepperTests.cpp
index 86a9e77cfab663da15114eedaf999d0952049641..1331eb1a06216edf0247167a033edf8578b81276 100644
--- a/Tests/Core/Propagator/StepperTests.cpp
+++ b/Tests/Core/Propagator/StepperTests.cpp
@@ -42,7 +42,6 @@ using namespace Acts::UnitLiterals;
namespace Acts {
namespace Test {
-using cstep = detail::ConstrainedStep;
using Covariance = BoundSymMatrix;
// Create a test context
diff --git a/Tests/Integration/BVHNavigationTest.cpp b/Tests/Integration/BVHNavigationTest.cpp
index 1e2f82b61a9c314c321b057b1a00e971de1400db..65097395725f15c9ac74ffc04dd3f2ffc29678bd 100644
--- a/Tests/Integration/BVHNavigationTest.cpp
+++ b/Tests/Integration/BVHNavigationTest.cpp
@@ -35,7 +35,7 @@
#include "Acts/Propagator/ActionList.hpp"
#include "Acts/Propagator/StraightLineStepper.hpp"
#include "Acts/Propagator/Propagator.hpp"
-#include "Acts/Propagator/detail/ConstrainedStep.hpp"
+#include "Acts/Propagator/ConstrainedStep.hpp"
#include "Acts/Propagator/detail/DebugOutputActor.hpp"
#include "Acts/Propagator/detail/StandardAborters.hpp"
#include "Acts/Propagator/detail/SteppingLogger.hpp"
diff --git a/Tests/Integration/PropagationTestBase.hpp b/Tests/Integration/PropagationTestBase.hpp
index 56034286855291cec81a5d3ee35ccdc011383284..98158db9b041950e4949dda2f496142b13617ae6 100644
--- a/Tests/Integration/PropagationTestBase.hpp
+++ b/Tests/Integration/PropagationTestBase.hpp
@@ -22,6 +22,8 @@
using namespace Acts::UnitLiterals;
+unsigned int itest = 0;
+
// datasets for Boost::Test data-driven test cases
namespace ds {
// track parameters
@@ -50,6 +52,9 @@ auto threeRandom = (rand1 ^ rand2 ^ rand2);
BOOST_DATA_TEST_CASE(forward_backward_propagation_,
ds::trackParameters* ds::propagationLimit, pT, phi, theta,
charge, plimit) {
+ std::cout << "Running (first patch) tests : " << itest << std::endl;
+ ++itest;
+
// foward backward check atlas stepper
foward_backward(apropagator, pT, phi, theta, charge, plimit, 1_um, 1_eV,
debug);
diff --git a/Tests/Integration/PropagationTestHelper.hpp b/Tests/Integration/PropagationTestHelper.hpp
index c764c27d29d0878c88aafaabce1caa5fb2eae7d5..29d5c08f1fc982000ec2a4b72be3cd4c954ae2b6 100644
--- a/Tests/Integration/PropagationTestHelper.hpp
+++ b/Tests/Integration/PropagationTestHelper.hpp
@@ -224,7 +224,7 @@ void foward_backward(const Propagator_type& propagator, double pT, double phi,
std::cout << " - resulted at position : "
<< fwdResult.endParameters->position() << std::endl;
- auto bwdOutput = fwdResult.template get<DebugOutput::result_type>();
+ auto bwdOutput = bwdResult.template get<DebugOutput::result_type>();
std::cout << ">>>>> Output for backward propagation " << std::endl;
std::cout << bwdOutput.debugString << std::endl;
std::cout << " - resulted at position : "