refactor: Align initial track param covariance (#2491)

Currently the initial track param covariance will depend on the chosen seeder. Truth smeared will use the configured sigmas to cook up a cov while to other seeder take a directly configured one.

Here I align the truth smeared to the others which also allows to experiment with the covariance separately from the smearing.

Examples/Algorithms/TrackFinding/include/ActsExamples/TrackFinding/TrackParamsEstimationAlgorithm.hpp

@@ -69,18 +69,11 @@ class TrackParamsEstimationAlgorithm final : public IAlgorithm {
     std::shared_ptr<const Acts::MagneticFieldProvider> magneticField;
     /// The minimum magnetic field to trigger the track parameters estimation
     double bFieldMin = 0.1 * Acts::UnitConstants::T;
-    /// Constant term of the loc0 resolution.
-    double sigmaLoc0 = 25 * Acts::UnitConstants::um;
-    /// Constant term of the loc1 resolution.
-    double sigmaLoc1 = 100 * Acts::UnitConstants::um;
-    /// Phi angular resolution.
-    double sigmaPhi = 0.02 * Acts::UnitConstants::degree;
-    /// Theta angular resolution.
-    double sigmaTheta = 0.02 * Acts::UnitConstants::degree;
-    /// q/p resolution.
-    double sigmaQOverP = 0.1 / Acts::UnitConstants::GeV;
-    /// Time resolution.
-    double sigmaT0 = 10 * Acts::UnitConstants::ns;
+    /// Initial covariance matrix diagonal.
+    std::array<double, 6> initialSigmas = {
+        25 * Acts::UnitConstants::um,       100 * Acts::UnitConstants::um,
+        0.02 * Acts::UnitConstants::degree, 0.02 * Acts::UnitConstants::degree,
+        0.1 / Acts::UnitConstants::GeV,     10 * Acts::UnitConstants::ns};
     /// Inflate initial covariance.
     std::array<double, 6> initialVarInflation = {1., 1., 1., 1., 1., 1.};
     /// Particle hypothesis.

Examples/Algorithms/TrackFinding/src/TrackParamsEstimationAlgorithm.cpp

@@ -58,24 +58,10 @@ ActsExamples::TrackParamsEstimationAlgorithm::TrackParamsEstimationAlgorithm(
   m_outputTracks.maybeInitialize(m_cfg.outputProtoTracks);
 
   // Set up the track parameters covariance (the same for all tracks)
-  m_covariance(Acts::eBoundLoc0, Acts::eBoundLoc0) =
-      m_cfg.initialVarInflation[Acts::eBoundLoc0] * cfg.sigmaLoc0 *
-      m_cfg.sigmaLoc0;
-  m_covariance(Acts::eBoundLoc1, Acts::eBoundLoc1) =
-      m_cfg.initialVarInflation[Acts::eBoundLoc1] * cfg.sigmaLoc1 *
-      m_cfg.sigmaLoc1;
-  m_covariance(Acts::eBoundPhi, Acts::eBoundPhi) =
-      m_cfg.initialVarInflation[Acts::eBoundPhi] * cfg.sigmaPhi *
-      m_cfg.sigmaPhi;
-  m_covariance(Acts::eBoundTheta, Acts::eBoundTheta) =
-      m_cfg.initialVarInflation[Acts::eBoundTheta] * cfg.sigmaTheta *
-      m_cfg.sigmaTheta;
-  m_covariance(Acts::eBoundQOverP, Acts::eBoundQOverP) =
-      m_cfg.initialVarInflation[Acts::eBoundQOverP] * cfg.sigmaQOverP *
-      m_cfg.sigmaQOverP;
-  m_covariance(Acts::eBoundTime, Acts::eBoundTime) =
-      m_cfg.initialVarInflation[Acts::eBoundTime] * m_cfg.sigmaT0 *
-      m_cfg.sigmaT0;
+  for (std::size_t i = Acts::eBoundLoc0; i < Acts::eBoundSize; ++i) {
+    m_covariance(i, i) = m_cfg.initialVarInflation[i] * m_cfg.initialSigmas[i] *
+                         m_cfg.initialSigmas[i];
+  }
 }
 
 ActsExamples::ProcessCode ActsExamples::TrackParamsEstimationAlgorithm::execute(

Examples/Algorithms/TruthTracking/ActsExamples/TruthTracking/ParticleSmearing.cpp

@@ -127,20 +127,29 @@ ActsExamples::ProcessCode ActsExamples::ParticleSmearing::execute(
 
       // build the track covariance matrix using the smearing sigmas
       Acts::BoundSquareMatrix cov = Acts::BoundSquareMatrix::Zero();
-      cov(Acts::eBoundLoc0, Acts::eBoundLoc0) =
-          m_cfg.initialVarInflation[Acts::eBoundLoc0] * sigmaD0 * sigmaD0;
-      cov(Acts::eBoundLoc1, Acts::eBoundLoc1) =
-          m_cfg.initialVarInflation[Acts::eBoundLoc1] * sigmaZ0 * sigmaZ0;
-      cov(Acts::eBoundTime, Acts::eBoundTime) =
-          m_cfg.initialVarInflation[Acts::eBoundTime] * sigmaT0 * sigmaT0;
-      cov(Acts::eBoundPhi, Acts::eBoundPhi) =
-          m_cfg.initialVarInflation[Acts::eBoundPhi] * sigmaPhi * sigmaPhi;
-      cov(Acts::eBoundTheta, Acts::eBoundTheta) =
-          m_cfg.initialVarInflation[Acts::eBoundTheta] * sigmaTheta *
-          sigmaTheta;
-      cov(Acts::eBoundQOverP, Acts::eBoundQOverP) =
-          m_cfg.initialVarInflation[Acts::eBoundQOverP] * sigmaQOverP *
-          sigmaQOverP;
+      if (m_cfg.initialSigmas) {
+        // use the initial sigmas if set
+        for (std::size_t i = Acts::eBoundLoc0; i < Acts::eBoundSize; ++i) {
+          cov(i, i) = m_cfg.initialVarInflation[i] * (*m_cfg.initialSigmas)[i] *
+                      (*m_cfg.initialSigmas)[i];
+        }
+      } else {
+        // otherwise use the smearing sigmas
+        cov(Acts::eBoundLoc0, Acts::eBoundLoc0) =
+            m_cfg.initialVarInflation[Acts::eBoundLoc0] * sigmaD0 * sigmaD0;
+        cov(Acts::eBoundLoc1, Acts::eBoundLoc1) =
+            m_cfg.initialVarInflation[Acts::eBoundLoc1] * sigmaZ0 * sigmaZ0;
+        cov(Acts::eBoundTime, Acts::eBoundTime) =
+            m_cfg.initialVarInflation[Acts::eBoundTime] * sigmaT0 * sigmaT0;
+        cov(Acts::eBoundPhi, Acts::eBoundPhi) =
+            m_cfg.initialVarInflation[Acts::eBoundPhi] * sigmaPhi * sigmaPhi;
+        cov(Acts::eBoundTheta, Acts::eBoundTheta) =
+            m_cfg.initialVarInflation[Acts::eBoundTheta] * sigmaTheta *
+            sigmaTheta;
+        cov(Acts::eBoundQOverP, Acts::eBoundQOverP) =
+            m_cfg.initialVarInflation[Acts::eBoundQOverP] * sigmaQOverP *
+            sigmaQOverP;
+      }
 
       parameters.emplace_back(perigee, params, cov, particleHypothesis);
     }

Examples/Algorithms/TruthTracking/ActsExamples/TruthTracking/ParticleSmearing.hpp

@@ -21,6 +21,7 @@
 #include <array>
 #include <limits>
 #include <memory>
+#include <optional>
 #include <string>
 
 namespace ActsExamples {
@@ -58,6 +59,9 @@ class ParticleSmearing final : public IAlgorithm {
     double sigmaTheta = 1 * Acts::UnitConstants::degree;
     /// Relative momentum resolution.
     double sigmaPRel = 0.05;
+    /// Optional. Initial covariance matrix diagonal. Overwrites the default if
+    /// set.
+    std::optional<std::array<double, 6>> initialSigmas = std::nullopt;
     /// Inflate the initial covariance matrix
     std::array<double, 6> initialVarInflation = {1., 1., 1., 1., 1., 1.};
     /// Random numbers service.

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

@@ -258,6 +258,7 @@ def addSeeding(
             rnd,
             selectedParticles,
             particleSmearingSigmas,
+            initialSigmas,
             initialVarInflation,
             particleHypothesis,
             logLevel,
@@ -317,12 +318,7 @@ def addSeeding(
             trackingGeometry=trackingGeometry,
             magneticField=field,
             **acts.examples.defaultKWArgs(
-                sigmaLoc0=initialSigmas[0] if initialSigmas is not None else None,
-                sigmaLoc1=initialSigmas[1] if initialSigmas is not None else None,
-                sigmaPhi=initialSigmas[2] if initialSigmas is not None else None,
-                sigmaTheta=initialSigmas[3] if initialSigmas is not None else None,
-                sigmaQOverP=initialSigmas[4] if initialSigmas is not None else None,
-                sigmaT0=initialSigmas[5] if initialSigmas is not None else None,
+                initialSigmas=initialSigmas,
                 initialVarInflation=initialVarInflation,
                 particleHypothesis=particleHypothesis,
             ),
@@ -393,6 +389,7 @@ def addTruthSmearedSeeding(
     rnd: Optional[acts.examples.RandomNumbers],
     selectedParticles: str,
     particleSmearingSigmas: ParticleSmearingSigmas,
+    initialSigmas: Optional[List[float]],
     initialVarInflation: List[float],
     particleHypothesis: Optional[acts.ParticleHypothesis],
     logLevel: acts.logging.Level = None,
@@ -420,6 +417,7 @@ def addTruthSmearedSeeding(
             sigmaPhi=particleSmearingSigmas.phi,
             sigmaTheta=particleSmearingSigmas.theta,
             sigmaPRel=particleSmearingSigmas.pRel,
+            initialSigmas=initialSigmas,
             initialVarInflation=initialVarInflation,
             particleHypothesis=particleHypothesis,
         ),

Examples/Python/src/TrackFinding.cpp

@@ -261,8 +261,8 @@ void addTrackFinding(Context& ctx) {
       ActsExamples::TrackParamsEstimationAlgorithm, mex,
       "TrackParamsEstimationAlgorithm", inputSeeds, inputProtoTracks,
       outputTrackParameters, outputSeeds, outputProtoTracks, trackingGeometry,
-      magneticField, bFieldMin, sigmaLoc0, sigmaLoc1, sigmaPhi, sigmaTheta,
-      sigmaQOverP, sigmaT0, initialVarInflation, particleHypothesis);
+      magneticField, bFieldMin, initialSigmas, initialVarInflation,
+      particleHypothesis);
 
   {
     using Alg = ActsExamples::TrackFindingAlgorithm;

Examples/Python/src/TruthTracking.cpp

@@ -90,7 +90,7 @@ void addTruthTracking(Context& ctx) {
       ActsExamples::ParticleSmearing, mex, "ParticleSmearing", inputParticles,
       outputTrackParameters, sigmaD0, sigmaD0PtA, sigmaD0PtB, sigmaZ0,
       sigmaZ0PtA, sigmaZ0PtB, sigmaT0, sigmaPhi, sigmaTheta, sigmaPRel,
-      initialVarInflation, particleHypothesis, randomNumbers);
+      initialSigmas, initialVarInflation, particleHypothesis, randomNumbers);
 
   {
     using Alg = ActsExamples::ParticleSelector;

Examples/Run/Reconstruction/Common/RecCKFTracks.cpp

@@ -245,12 +245,8 @@ int runRecCKFTracks(
     paramsEstimationCfg.trackingGeometry = trackingGeometry;
     paramsEstimationCfg.magneticField = magneticField;
     paramsEstimationCfg.bFieldMin = 0.1_T;
-    paramsEstimationCfg.sigmaLoc0 = 25._um;
-    paramsEstimationCfg.sigmaLoc1 = 100._um;
-    paramsEstimationCfg.sigmaPhi = 0.02_degree;
-    paramsEstimationCfg.sigmaTheta = 0.02_degree;
-    paramsEstimationCfg.sigmaQOverP = 0.1 / 1._GeV;
-    paramsEstimationCfg.sigmaT0 = 1400._s;
+    paramsEstimationCfg.initialSigmas = {25._um,      100._um,      0.02_degree,
+                                         0.02_degree, 0.1 / 1._GeV, 1400._s};
     paramsEstimationCfg.initialVarInflation =
         vm["ckf-initial-variance-inflation"].template as<Options::Reals<6>>();