From ae5640603e5117f3168c606da490a9b23e9bbc28 Mon Sep 17 00:00:00 2001
From: "Alexander J. Pfleger" <70842573+AJPfleger@users.noreply.github.com>
Date: Thu, 22 Feb 2024 16:52:24 +0100
Subject: [PATCH] fix(gx2f): propagate final covariance for trackstates (#2949)

## Issue
I figured out, that the trackstates in the final track contain the wrong covariance. This occurred, while investigating the unbiased pulls. There, some dimensions had a negative covariance.

https://github.com/acts-project/acts/blob/19464fa0d24f031083ef3ad3108b74f2b5fa5454/Examples/Io/Root/src/RootTrackStatesWriter.cpp#L456

## Reason
After the last update we calculate the covariance of the initial parameters and add this covariance to the final track. However, the last propagation was done with the initial guess for the covariance. Therefore, the wrong covariance got propagated and assigned to each trackstate.

## Solution
Propagate an additional time, with the final parameters+covariance.

## Future Plans
Since an additional propagation is quite expensive, we should look into the following two ideas:
- calculate covariance after each update (matrix inverse)
- toggle the re-propagation if we are not interested in the propagated covariance

## Notes
This also removes `nUpdate` from the Actor since it wasn't needed anymore.

## Blocked by
- https://github.com/acts-project/acts/pull/2972
- https://github.com/acts-project/acts/pull/2966
---
 .../TrackFitting/GlobalChiSquareFitter.hpp    | 39 +++++++++++++++----
 1 file changed, 31 insertions(+), 8 deletions(-)

diff --git a/Core/include/Acts/TrackFitting/GlobalChiSquareFitter.hpp b/Core/include/Acts/TrackFitting/GlobalChiSquareFitter.hpp
index ef7b6b559c4..a3a5cd637e5 100644
--- a/Core/include/Acts/TrackFitting/GlobalChiSquareFitter.hpp
+++ b/Core/include/Acts/TrackFitting/GlobalChiSquareFitter.hpp
@@ -369,13 +369,6 @@ class Gx2Fitter {
     /// Calibration context for the fit
     const CalibrationContext* calibrationContext{nullptr};
 
-    /// The current iteration of the fitter.
-    /// The variable is updated in fit().
-    /// The actor needs to know the current iteration for adding new
-    /// trackStates. During the first iteration, each measurement surfaces will
-    /// be added to the track.
-    std::size_t nUpdate = Acts::MultiTrajectoryTraits::kInvalid;
-
     /// @brief Gx2f actor operation
     ///
     /// @tparam propagator_state_t is the type of Propagator state
@@ -648,7 +641,6 @@ class Gx2Fitter {
       gx2fActor.extensions = gx2fOptions.extensions;
       gx2fActor.calibrationContext = &gx2fOptions.calibrationContext.get();
       gx2fActor.actorLogger = m_actorLogger.get();
-      gx2fActor.nUpdate = nUpdate;
 
       auto propagatorState = m_propagator.makeState(params, propagatorOptions);
 
@@ -797,6 +789,37 @@ class Gx2Fitter {
 
     ACTS_VERBOSE("final covariance:\n" << fullCovariancePredicted);
 
+    // Propagate again with the final covariance matrix. This is necessary to
+    // obtain the propagated covariance for each state.
+    if (gx2fOptions.nUpdateMax > 0) {
+      ACTS_VERBOSE("Propagate with the final covariance.");
+      // update covariance
+      ACTS_VERBOSE("finaldeltaParams:\n" << deltaParams);
+      params.covariance() = fullCovariancePredicted;
+
+      // set up propagator and co
+      Acts::GeometryContext geoCtx = gx2fOptions.geoContext;
+      Acts::MagneticFieldContext magCtx = gx2fOptions.magFieldContext;
+      // Set options for propagator
+      PropagatorOptions propagatorOptions(geoCtx, magCtx);
+      auto& gx2fActor = propagatorOptions.actionList.template get<GX2FActor>();
+      gx2fActor.inputMeasurements = &inputMeasurements;
+      gx2fActor.extensions = gx2fOptions.extensions;
+      gx2fActor.calibrationContext = &gx2fOptions.calibrationContext.get();
+      gx2fActor.actorLogger = m_actorLogger.get();
+
+      auto propagatorState = m_propagator.makeState(params, propagatorOptions);
+
+      auto& r = propagatorState.template get<Gx2FitterResult<traj_t>>();
+      r.fittedStates = &trackContainer.trackStateContainer();
+
+      // Clear the track container. It could be more performant to update the
+      // existing states, but this needs some more thinking.
+      trackContainer.clear();
+
+      m_propagator.template propagate(propagatorState);
+    }
+
     if (!trackContainer.hasColumn(
             Acts::hashString(Gx2fConstants::gx2fnUpdateColumn))) {
       trackContainer.template addColumn<std::size_t>("Gx2fnUpdateColumn");