feat: Add/activate reverse filter cov inflation for KF and GSF

CI/physmon/workflows/physmon_trackfitting_kf.py

@@ -24,6 +24,8 @@
         field=setup.field,
         digiConfigFile=setup.digiConfig,
         outputDir=tp,
+        reverseFilteringMomThreshold=float("inf"),
+        reverseFilteringCovarianceScaling=100.0,
         s=s,
     )
 

CI/physmon/workflows/physmon_trackrefitting_kf.py

@@ -24,6 +24,8 @@
         field=setup.field,
         digiConfigFile=setup.digiConfig,
         outputDir=tp,
+        reverseFilteringMomThreshold=float("inf"),
+        reverseFilteringCovarianceScaling=100.0,
         s=s,
     )
 

Core/include/Acts/TrackFitting/GaussianSumFitter.hpp

@@ -321,11 +321,23 @@ struct GaussianSumFitter {
                                   : sParameters.referenceSurface();
 
       const auto& params = *fwdGsfResult.lastMeasurementState;
+      std::vector<
+          std::tuple<double, BoundVector, std::optional<BoundSquareMatrix>>>
+          inflatedParamVector;
+      for (auto& [w, p, cov] : params.components()) {
+        inflatedParamVector.emplace_back(
+            w, p, cov.value() * options.reverseFilteringCovarianceScaling);
+      }
+
+      MultiComponentBoundTrackParameters inflatedParams(
+          params.referenceSurface().getSharedPtr(),
+          std::move(inflatedParamVector), params.particleHypothesis());
+
       auto state =
           m_propagator.template makeState<MultiComponentBoundTrackParameters,
                                           decltype(bwdPropOptions),
                                           MultiStepperSurfaceReached>(
-              params, target, bwdPropOptions);
+              inflatedParams, target, bwdPropOptions);
 
       assert(
           (fwdGsfResult.lastMeasurementTip != MultiTrajectoryTraits::kInvalid &&

Core/include/Acts/TrackFitting/GsfOptions.hpp

@@ -112,6 +112,8 @@ struct GsfOptions {
 
   bool disableAllMaterialHandling = false;
 
+  double reverseFilteringCovarianceScaling = 1.0;
+
   std::string_view finalMultiComponentStateColumn = "";
 
   ComponentMergeMethod componentMergeMethod = ComponentMergeMethod::eMaxWeight;

Core/include/Acts/TrackFitting/KalmanFitter.hpp

@@ -484,7 +484,9 @@ class KalmanFitter {
             result.fittedStates->applyBackwards(
                 result.lastMeasurementIndex, [&](auto trackState) {
                   auto fSurface = &trackState.referenceSurface();
-                  if (!rangeContainsValue(result.passedAgainSurfaces,
+                  if (trackState.typeFlags().test(
+                          TrackStateFlag::MeasurementFlag) &&
+                      !rangeContainsValue(result.passedAgainSurfaces,
                                           fSurface)) {
                     // If reversed filtering missed this surface, then there is
                     // no smoothed parameter

Examples/Algorithms/TrackFitting/include/ActsExamples/TrackFitting/TrackFitterFunction.hpp

@@ -66,6 +66,7 @@ std::shared_ptr<TrackFitterFunction> makeKalmanFitterFunction(
     std::shared_ptr<const Acts::MagneticFieldProvider> magneticField,
     bool multipleScattering = true, bool energyLoss = true,
     double reverseFilteringMomThreshold = 0.0,
+    double reverseFilteringCovarianceScaling = 1.0,
     Acts::FreeToBoundCorrection freeToBoundCorrection =
         Acts::FreeToBoundCorrection(),
     const Acts::Logger& logger = *Acts::getDefaultLogger("Kalman",
@@ -89,14 +90,16 @@ enum class MixtureReductionAlgorithm { weightCut, KLDistance };
 /// parameters and covariance
 /// @param mixtureReductionAlgorithm How to reduce the number of components
 /// in a mixture
+/// @param reverseFilteringCovarianceScaling How the covariance matrices are
+/// inflated before the reverse filtering pass
 /// @param logger a logger instance
 std::shared_ptr<TrackFitterFunction> makeGsfFitterFunction(
     std::shared_ptr<const Acts::TrackingGeometry> trackingGeometry,
     std::shared_ptr<const Acts::MagneticFieldProvider> magneticField,
     BetheHeitlerApprox betheHeitlerApprox, std::size_t maxComponents,
     double weightCutoff, Acts::ComponentMergeMethod componentMergeMethod,
     MixtureReductionAlgorithm mixtureReductionAlgorithm,
-    const Acts::Logger& logger);
+    double reverseFilteringCovarianceScaling, const Acts::Logger& logger);
 
 /// Makes a fitter function object for the Global Chi Square Fitter (GX2F)
 ///

Examples/Algorithms/TrackFitting/src/GsfFitterFunction.cpp

@@ -87,6 +87,7 @@ struct GsfFitterFunctionImpl final : public ActsExamples::TrackFitterFunction {
       MixtureReductionAlgorithm::KLDistance;
   Acts::ComponentMergeMethod mergeMethod =
       Acts::ComponentMergeMethod::eMaxWeight;
+  double reverseFilteringCovarianceScaling = 1.0;
 
   IndexSourceLink::SurfaceAccessor m_slSurfaceAccessor;
 
@@ -115,6 +116,8 @@ struct GsfFitterFunctionImpl final : public ActsExamples::TrackFitterFunction {
     gsfOptions.abortOnError = abortOnError;
     gsfOptions.disableAllMaterialHandling = disableAllMaterialHandling;
     gsfOptions.componentMergeMethod = mergeMethod;
+    gsfOptions.reverseFilteringCovarianceScaling =
+        reverseFilteringCovarianceScaling;
 
     gsfOptions.extensions.calibrator.connect<&calibrator_t::calibrate>(
         &calibrator);
@@ -195,7 +198,7 @@ std::shared_ptr<TrackFitterFunction> ActsExamples::makeGsfFitterFunction(
     BetheHeitlerApprox betheHeitlerApprox, std::size_t maxComponents,
     double weightCutoff, Acts::ComponentMergeMethod componentMergeMethod,
     MixtureReductionAlgorithm mixtureReductionAlgorithm,
-    const Acts::Logger& logger) {
+    double reverseFilteringCovarianceScaling, const Acts::Logger& logger) {
   // Standard fitter
   MultiStepper stepper(magneticField, logger.cloneWithSuffix("Step"));
   const auto& geo = *trackingGeometry;
@@ -229,6 +232,8 @@ std::shared_ptr<TrackFitterFunction> ActsExamples::makeGsfFitterFunction(
   fitterFunction->weightCutoff = weightCutoff;
   fitterFunction->mergeMethod = componentMergeMethod;
   fitterFunction->reductionAlg = mixtureReductionAlgorithm;
+  fitterFunction->reverseFilteringCovarianceScaling =
+      reverseFilteringCovarianceScaling;
 
   return fitterFunction;
 }

Examples/Algorithms/TrackFitting/src/KalmanFitterFunction.cpp

@@ -76,7 +76,7 @@ struct KalmanFitterFunctionImpl final : public TrackFitterFunction {
   Acts::GainMatrixUpdater kfUpdater;
   Acts::GainMatrixSmoother kfSmoother;
   SimpleReverseFilteringLogic reverseFilteringLogic;
-
+  double reverseFilteringCovarianceScaling;
   bool multipleScattering = false;
   bool energyLoss = false;
   Acts::FreeToBoundCorrection freeToBoundCorrection;
@@ -114,6 +114,8 @@ struct KalmanFitterFunctionImpl final : public TrackFitterFunction {
     kfOptions.freeToBoundCorrection = freeToBoundCorrection;
     kfOptions.extensions.calibrator.connect<&calibrator_t::calibrate>(
         &calibrator);
+    kfOptions.reversedFilteringCovarianceScaling =
+        reverseFilteringCovarianceScaling;
 
     if (options.doRefit) {
       kfOptions.extensions.surfaceAccessor
@@ -159,6 +161,7 @@ ActsExamples::makeKalmanFitterFunction(
     std::shared_ptr<const Acts::MagneticFieldProvider> magneticField,
     bool multipleScattering, bool energyLoss,
     double reverseFilteringMomThreshold,
+    double reverseFilteringCovarianceScaling,
     Acts::FreeToBoundCorrection freeToBoundCorrection,
     const Acts::Logger& logger) {
   // Stepper should be copied into the fitters
@@ -191,6 +194,8 @@ ActsExamples::makeKalmanFitterFunction(
   fitterFunction->reverseFilteringLogic.momentumThreshold =
       reverseFilteringMomThreshold;
   fitterFunction->freeToBoundCorrection = freeToBoundCorrection;
+  fitterFunction->reverseFilteringCovarianceScaling =
+      reverseFilteringCovarianceScaling;
 
   return fitterFunction;
 }

Examples/Python/python/acts/examples/reconstruction.py

@@ -1289,6 +1289,7 @@ def addKalmanTracks(
     trackingGeometry: acts.TrackingGeometry,
     field: acts.MagneticFieldProvider,
     reverseFilteringMomThreshold: float = 0 * u.GeV,
+    reverseFilteringCovarianceScaling: float = 1.0,
     inputProtoTracks: str = "truth_particle_tracks",
     multipleScattering: bool = True,
     energyLoss: bool = True,
@@ -1302,6 +1303,7 @@ def addKalmanTracks(
         "multipleScattering": multipleScattering,
         "energyLoss": energyLoss,
         "reverseFilteringMomThreshold": reverseFilteringMomThreshold,
+        "reverseFilteringCovarianceScaling": reverseFilteringCovarianceScaling,
         "freeToBoundCorrection": acts.examples.FreeToBoundCorrection(False),
         "level": customLogLevel(),
     }
@@ -1362,6 +1364,7 @@ def addTruthTrackingGsf(
         "componentMergeMethod": acts.examples.ComponentMergeMethod.maxWeight,
         "mixtureReductionAlgorithm": acts.examples.MixtureReductionAlgorithm.KLDistance,
         "weightCutoff": 1.0e-4,
+        "reverseFilteringCovarianceScaling": 100.0,
         "level": customLogLevel(),
     }
 

Examples/Python/src/TrackFitting.cpp

@@ -53,17 +53,18 @@ void addTrackFitting(Context& ctx) {
            std::shared_ptr<const Acts::MagneticFieldProvider> magneticField,
            bool multipleScattering, bool energyLoss,
            double reverseFilteringMomThreshold,
+           double reverseFilteringCovarianceScaling,
            Acts::FreeToBoundCorrection freeToBoundCorrection,
            Logging::Level level) {
           return ActsExamples::makeKalmanFitterFunction(
               trackingGeometry, magneticField, multipleScattering, energyLoss,
-              reverseFilteringMomThreshold, freeToBoundCorrection,
-              *Acts::getDefaultLogger("Kalman", level));
+              reverseFilteringMomThreshold, reverseFilteringCovarianceScaling,
+              freeToBoundCorrection, *Acts::getDefaultLogger("Kalman", level));
         },
-        py::arg("trackingGeometry"), py::arg("magneticField"),
-        py::arg("multipleScattering"), py::arg("energyLoss"),
-        py::arg("reverseFilteringMomThreshold"),
-        py::arg("freeToBoundCorrection"), py::arg("level"));
+        "trackingGeometry"_a, "magneticField"_a, "multipleScattering"_a,
+        "energyLoss"_a, "reverseFilteringMomThreshold"_a,
+        "reverseFilteringCovarianceScaling"_a, "freeToBoundCorrection"_a,
+        "level"_a);
 
     py::class_<MeasurementCalibrator, std::shared_ptr<MeasurementCalibrator>>(
         mex, "MeasurementCalibrator");
@@ -108,17 +109,17 @@ void addTrackFitting(Context& ctx) {
            BetheHeitlerApprox betheHeitlerApprox, std::size_t maxComponents,
            double weightCutoff, Acts::ComponentMergeMethod componentMergeMethod,
            ActsExamples::MixtureReductionAlgorithm mixtureReductionAlgorithm,
-           Logging::Level level) {
+           double reverseFilteringCovarianceScaling, Logging::Level level) {
           return ActsExamples::makeGsfFitterFunction(
               trackingGeometry, magneticField, betheHeitlerApprox,
               maxComponents, weightCutoff, componentMergeMethod,
-              mixtureReductionAlgorithm,
+              mixtureReductionAlgorithm, reverseFilteringCovarianceScaling,
               *Acts::getDefaultLogger("GSFFunc", level));
         },
-        py::arg("trackingGeometry"), py::arg("magneticField"),
-        py::arg("betheHeitlerApprox"), py::arg("maxComponents"),
-        py::arg("weightCutoff"), py::arg("componentMergeMethod"),
-        py::arg("mixtureReductionAlgorithm"), py::arg("level"));
+        "trackingGeometry"_a, "magneticField"_a, "betheHeitlerApprox"_a,
+        "maxComponents"_a, "weightCutoff"_a, "componentMergeMethod"_a,
+        "mixtureReductionAlgorithm"_a, "reverseFilteringCovarianceScaling"_a,
+        "level"_a);
 
     mex.def(
         "makeGlobalChiSquareFitterFunction",

Examples/Scripts/Python/truth_tracking_gsf_refitting.py

@@ -22,6 +22,8 @@ def runRefittingGsf(
         field,
         digiConfigFile=digiConfigFile,
         outputDir=outputDir,
+        reverseFilteringMomThreshold=0.0,
+        reverseFilteringCovarianceScaling=1.0,
         s=s,
     )
 
@@ -36,6 +38,7 @@ def runRefittingGsf(
         "componentMergeMethod": acts.examples.ComponentMergeMethod.maxWeight,
         "mixtureReductionAlgorithm": acts.examples.MixtureReductionAlgorithm.KLDistance,
         "weightCutoff": 1.0e-4,
+        "reverseFilteringCovarianceScaling": 100.0,
         "level": acts.logging.INFO,
     }
 

Examples/Scripts/Python/truth_tracking_kalman.py

@@ -18,6 +18,7 @@ def runTruthTrackingKalman(
     inputHitsPath: Optional[Path] = None,
     decorators=[],
     reverseFilteringMomThreshold=0 * u.GeV,
+    reverseFilteringCovarianceScaling=1,
     s: acts.examples.Sequencer = None,
 ):
     from acts.examples.simulation import (
@@ -122,6 +123,7 @@ def runTruthTrackingKalman(
         trackingGeometry,
         field,
         reverseFilteringMomThreshold,
+        reverseFilteringCovarianceScaling,
     )
 
     s.addAlgorithm(

Examples/Scripts/Python/truth_tracking_kalman_refitting.py

@@ -18,20 +18,24 @@ def runRefittingKf(
     multipleScattering: bool = True,
     energyLoss: bool = True,
     reverseFilteringMomThreshold=0 * u.GeV,
+    reverseFilteringCovarianceScaling=1.0,
     s: acts.examples.Sequencer = None,
 ):
     s = runTruthTrackingKalman(
         trackingGeometry,
         field,
         digiConfigFile=digiConfigFile,
         outputDir=outputDir,
+        reverseFilteringMomThreshold=reverseFilteringMomThreshold,
+        reverseFilteringCovarianceScaling=reverseFilteringCovarianceScaling,
         s=s,
     )
 
     kalmanOptions = {
         "multipleScattering": multipleScattering,
         "energyLoss": energyLoss,
         "reverseFilteringMomThreshold": reverseFilteringMomThreshold,
+        "reverseFilteringCovarianceScaling": reverseFilteringCovarianceScaling,
         "freeToBoundCorrection": acts.examples.FreeToBoundCorrection(False),
         "level": acts.logging.INFO,
     }