Skip to content

Commit

Permalink
refactor: Refactor target direction logic in KF/CKF (#2539)
Browse files Browse the repository at this point in the history
Not sure what the use case for this was but to me it seems like the target propagation logic is flawed in KF/CKF. If we propagate to target from the last state the smoothing seems unnecessary. Also since the target surface is a perigee surface for us there are edge cases where the last state has a smaller "distance" to the perigee than the first one.

In general I'd say the user should take such a decision but this is too big of a change at the moment.

discovered with #2518
  • Loading branch information
andiwand committed Oct 19, 2023
1 parent 25331d0 commit fd80ca7
Show file tree
Hide file tree
Showing 4 changed files with 97 additions and 24 deletions.
57 changes: 47 additions & 10 deletions Core/include/Acts/TrackFinding/CombinatorialKalmanFilter.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -47,6 +47,16 @@

namespace Acts {

enum class CombinatorialKalmanFilterTargetSurfaceStrategy {
/// Use the first trackstate to reach target surface
first,
/// Use the last trackstate to reach target surface
last,
/// Use the first or last trackstate to reach target surface depending on the
/// distance
firstOrLast,
};

/// Track quality summary for one trajectory.
///
/// This could be used to decide if a track is to be recorded when the
Expand Down Expand Up @@ -200,6 +210,11 @@ struct CombinatorialKalmanFilterOptions {
/// The smoothing target surface
const Surface* smoothingTargetSurface = nullptr;

/// Strategy to propagate to reference surface
CombinatorialKalmanFilterTargetSurfaceStrategy
smoothingTargetSurfaceStrategy =
CombinatorialKalmanFilterTargetSurfaceStrategy::firstOrLast;

/// Whether to consider multiple scattering.
bool multipleScattering = true;

Expand Down Expand Up @@ -334,6 +349,11 @@ class CombinatorialKalmanFilter {
/// The smoothing target surface
const Surface* smoothingTargetSurface = nullptr;

/// Strategy to propagate to reference surface
CombinatorialKalmanFilterTargetSurfaceStrategy
smoothingTargetSurfaceStrategy =
CombinatorialKalmanFilterTargetSurfaceStrategy::firstOrLast;

/// Whether to consider multiple scattering.
bool multipleScattering = true;

Expand Down Expand Up @@ -1183,6 +1203,7 @@ class CombinatorialKalmanFilter {
// Screen output for debugging
ACTS_VERBOSE("Apply smoothing on " << nStates
<< " filtered track states.");

// Smooth the track states
auto smoothRes =
m_extensions.smoother(state.geoContext, *result.fittedStates,
Expand All @@ -1197,8 +1218,7 @@ class CombinatorialKalmanFilter {
return Result<void>::success();
}

// Obtain the smoothed parameters at first/last measurement state.
// The first state can also be a material state
// Obtain the smoothed parameters at first/last measurement state
auto firstCreatedState =
result.fittedStates->getTrackState(firstStateIndex);
auto lastCreatedMeasurement =
Expand All @@ -1214,7 +1234,8 @@ class CombinatorialKalmanFilter {
.closest();
};

// The smoothed free params at the first/last measurement state
// The smoothed free params at the first/last measurement state.
// (the first state can also be a material state)
auto firstParams = MultiTrajectoryHelpers::freeSmoothed(
state.options.geoContext, firstCreatedState);
auto lastParams = MultiTrajectoryHelpers::freeSmoothed(
Expand All @@ -1225,18 +1246,31 @@ class CombinatorialKalmanFilter {
const auto lastIntersection = target(lastParams);

// Update the stepping parameters - in order to progress to destination.
// At the same time, reverse navigation direction for further
// stepping if necessary.
// At the same time, reverse navigation direction for further stepping if
// necessary.
// @note The stepping parameters is updated to the smoothed parameters at
// either the first measurement state or the last measurement state. It
// assumes the target surface is not within the first and the last
// smoothed measurement state. Also, whether the intersection is on
// surface is not checked here.
bool closerToFirstCreatedState =
std::abs(firstIntersection.pathLength()) <=
std::abs(lastIntersection.pathLength());
bool useFirstTrackState = true;
switch (smoothingTargetSurfaceStrategy) {
case CombinatorialKalmanFilterTargetSurfaceStrategy::first:
useFirstTrackState = true;
break;
case CombinatorialKalmanFilterTargetSurfaceStrategy::last:
useFirstTrackState = false;
break;
case CombinatorialKalmanFilterTargetSurfaceStrategy::firstOrLast:
useFirstTrackState = std::abs(firstIntersection.pathLength()) <=
std::abs(lastIntersection.pathLength());
break;
default:
ACTS_ERROR("Unknown target surface strategy");
return KalmanFitterError::SmoothFailed;
}
bool reverseDirection = false;
if (closerToFirstCreatedState) {
if (useFirstTrackState) {
stepper.resetState(state.stepping, firstCreatedState.smoothed(),
firstCreatedState.smoothedCovariance(),
firstCreatedState.referenceSurface(),
Expand All @@ -1256,9 +1290,10 @@ class CombinatorialKalmanFilter {
"target surface");
state.options.direction = state.options.direction.invert();
}
const auto& surface = closerToFirstCreatedState
const auto& surface = useFirstTrackState
? firstCreatedState.referenceSurface()
: lastCreatedMeasurement.referenceSurface();

ACTS_VERBOSE(
"Smoothing successful, updating stepping state to smoothed "
"parameters at surface "
Expand Down Expand Up @@ -1374,6 +1409,8 @@ class CombinatorialKalmanFilter {
propOptions.actionList.template get<CombinatorialKalmanFilterActor>();
combKalmanActor.filterTargetSurface = tfOptions.filterTargetSurface;
combKalmanActor.smoothingTargetSurface = tfOptions.smoothingTargetSurface;
combKalmanActor.smoothingTargetSurfaceStrategy =
tfOptions.smoothingTargetSurfaceStrategy;
combKalmanActor.multipleScattering = tfOptions.multipleScattering;
combKalmanActor.energyLoss = tfOptions.energyLoss;
combKalmanActor.smoothing = tfOptions.smoothing;
Expand Down
53 changes: 42 additions & 11 deletions Core/include/Acts/TrackFitting/KalmanFitter.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -46,6 +46,16 @@

namespace Acts {

enum class KalmanFitterTargetSurfaceStrategy {
/// Use the first trackstate to reach target surface
first,
/// Use the last trackstate to reach target surface
last,
/// Use the first or last trackstate to reach target surface depending on the
/// distance
firstOrLast,
};

/// Extension struct which holds delegates to customise the KF behavior
template <typename traj_t>
struct KalmanFitterExtensions {
Expand Down Expand Up @@ -157,6 +167,10 @@ struct KalmanFitterOptions {
/// The reference Surface
const Surface* referenceSurface = nullptr;

/// Strategy to propagate to reference surface
KalmanFitterTargetSurfaceStrategy referenceSurfaceStrategy =
KalmanFitterTargetSurfaceStrategy::firstOrLast;

/// Whether to consider multiple scattering
bool multipleScattering = true;

Expand Down Expand Up @@ -293,6 +307,10 @@ class KalmanFitter {
/// The target surface
const Surface* targetSurface = nullptr;

/// Strategy to propagate to target surface
KalmanFitterTargetSurfaceStrategy targetSurfaceStrategy =
KalmanFitterTargetSurfaceStrategy::firstOrLast;

/// Allows retrieving measurements for a surface
const std::map<GeometryIdentifier, SourceLink>* inputMeasurements = nullptr;

Expand Down Expand Up @@ -916,10 +934,8 @@ class KalmanFitter {
return KalmanFitterError::SmoothFailed;
}
// Screen output for debugging
if (logger().doPrint(Logging::VERBOSE)) {
ACTS_VERBOSE("Apply smoothing on " << nStates
<< " filtered track states.");
}
ACTS_VERBOSE("Apply smoothing on " << nStates
<< " filtered track states.");

// Smooth the track states
auto smoothRes =
Expand Down Expand Up @@ -963,18 +979,31 @@ class KalmanFitter {
const auto lastIntersection = target(lastParams);

// Update the stepping parameters - in order to progress to destination.
// At the same time, reverse navigation direction for further
// stepping if necessary.
// At the same time, reverse navigation direction for further stepping if
// necessary.
// @note The stepping parameters is updated to the smoothed parameters at
// either the first measurement state or the last measurement state. It
// assumes the target surface is not within the first and the last
// smoothed measurement state. Also, whether the intersection is on
// surface is not checked here.
bool closerToFirstCreatedState =
std::abs(firstIntersection.pathLength()) <=
std::abs(lastIntersection.pathLength());
bool useFirstTrackState = true;
switch (targetSurfaceStrategy) {
case KalmanFitterTargetSurfaceStrategy::first:
useFirstTrackState = true;
break;
case KalmanFitterTargetSurfaceStrategy::last:
useFirstTrackState = false;
break;
case KalmanFitterTargetSurfaceStrategy::firstOrLast:
useFirstTrackState = std::abs(firstIntersection.pathLength()) <=
std::abs(lastIntersection.pathLength());
break;
default:
ACTS_ERROR("Unknown target surface strategy");
return KalmanFitterError::SmoothFailed;
}
bool reverseDirection = false;
if (closerToFirstCreatedState) {
if (useFirstTrackState) {
stepper.resetState(state.stepping, firstCreatedState.smoothed(),
firstCreatedState.smoothedCovariance(),
firstCreatedState.referenceSurface(),
Expand All @@ -994,9 +1023,10 @@ class KalmanFitter {
"target surface");
state.options.direction = state.options.direction.invert();
}
const auto& surface = closerToFirstCreatedState
const auto& surface = useFirstTrackState
? firstCreatedState.referenceSurface()
: lastCreatedMeasurement.referenceSurface();

ACTS_VERBOSE(
"Smoothing successful, updating stepping state to smoothed "
"parameters at surface "
Expand Down Expand Up @@ -1230,6 +1260,7 @@ class KalmanFitter {
auto& kalmanActor = kalmanOptions.actionList.template get<KalmanActor>();
kalmanActor.inputMeasurements = &inputMeasurements;
kalmanActor.targetSurface = kfOptions.referenceSurface;
kalmanActor.targetSurfaceStrategy = kfOptions.referenceSurfaceStrategy;
kalmanActor.multipleScattering = kfOptions.multipleScattering;
kalmanActor.energyLoss = kfOptions.energyLoss;
kalmanActor.reversedFiltering = kfOptions.reversedFiltering;
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -17,6 +17,7 @@
#include "Acts/Propagator/Propagator.hpp"
#include "Acts/Surfaces/PerigeeSurface.hpp"
#include "Acts/Surfaces/Surface.hpp"
#include "Acts/TrackFinding/CombinatorialKalmanFilter.hpp"
#include "Acts/TrackFitting/GainMatrixSmoother.hpp"
#include "Acts/TrackFitting/GainMatrixUpdater.hpp"
#include "Acts/TrackFitting/KalmanFitter.hpp"
Expand Down Expand Up @@ -115,6 +116,8 @@ ActsExamples::ProcessCode ActsExamples::TrackFindingAlgorithm::execute(
ActsExamples::TrackFindingAlgorithm::TrackFinderOptions options(
ctx.geoContext, ctx.magFieldContext, ctx.calibContext, slAccessorDelegate,
extensions, pOptions, pSurface.get());
options.smoothingTargetSurfaceStrategy =
Acts::CombinatorialKalmanFilterTargetSurfaceStrategy::first;

// Perform the track finding for all initial parameters
ACTS_DEBUG("Invoke track finding with " << initialParameters.size()
Expand Down
8 changes: 5 additions & 3 deletions Examples/Algorithms/TrackFitting/src/KalmanFitterFunction.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -81,13 +81,13 @@ struct KalmanFitterFunctionImpl final : public TrackFitterFunction {
bool energyLoss = false;
Acts::FreeToBoundCorrection freeToBoundCorrection;

IndexSourceLink::SurfaceAccessor m_slSurfaceAccessor;
IndexSourceLink::SurfaceAccessor slSurfaceAccessor;

KalmanFitterFunctionImpl(Fitter&& f, DirectFitter&& df,
const Acts::TrackingGeometry& trkGeo)
: fitter(std::move(f)),
directFitter(std::move(df)),
m_slSurfaceAccessor{trkGeo} {}
slSurfaceAccessor{trkGeo} {}

template <typename calibrator_t>
auto makeKfOptions(const GeneralFitterOptions& options,
Expand All @@ -107,14 +107,16 @@ struct KalmanFitterFunctionImpl final : public TrackFitterFunction {
options.geoContext, options.magFieldContext, options.calibrationContext,
extensions, options.propOptions, &(*options.referenceSurface));

kfOptions.referenceSurfaceStrategy =
Acts::KalmanFitterTargetSurfaceStrategy::first;
kfOptions.multipleScattering = multipleScattering;
kfOptions.energyLoss = energyLoss;
kfOptions.freeToBoundCorrection = freeToBoundCorrection;
kfOptions.extensions.calibrator.connect<&calibrator_t::calibrate>(
&calibrator);
kfOptions.extensions.surfaceAccessor
.connect<&IndexSourceLink::SurfaceAccessor::operator()>(
&m_slSurfaceAccessor);
&slSurfaceAccessor);

return kfOptions;
}
Expand Down

0 comments on commit fd80ca7

Please sign in to comment.