feat: Initial sigma qop coefficients for Examples `TrackParamsEstim…

…ationAlgorithm` (acts-project#3111)

Adds a `qop` coefficient term to the initial covariance of the estimated track parameters. Since the resolution of `qop` or `p` usually depends on the momentum itself this should allow us to model the initial curvature uncertainty more effectively and improve the track finding procedure.

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

@@ -71,9 +71,25 @@ class TrackParamsEstimationAlgorithm final : public IAlgorithm {
     double bFieldMin = 0.1 * Acts::UnitConstants::T;
     /// Initial covariance matrix diagonal.
     std::array<double, 6> initialSigmas = {
-        1 * Acts::UnitConstants::mm,     1 * Acts::UnitConstants::mm,
-        1 * Acts::UnitConstants::degree, 1 * Acts::UnitConstants::degree,
-        0.1 / Acts::UnitConstants::GeV,  1 * Acts::UnitConstants::ns};
+        1 * Acts::UnitConstants::mm,
+        1 * Acts::UnitConstants::mm,
+        1 * Acts::UnitConstants::degree,
+        1 * Acts::UnitConstants::degree,
+        0 * Acts::UnitConstants::e / Acts::UnitConstants::GeV,
+        1 * Acts::UnitConstants::ns};
+    /// Initial q/p coefficient covariance matrix diagonal.
+    std::array<double, 6> initialSimgaQoverPCoefficients = {
+        0 * Acts::UnitConstants::mm /
+            (Acts::UnitConstants::e * Acts::UnitConstants::GeV),
+        0 * Acts::UnitConstants::mm /
+            (Acts::UnitConstants::e * Acts::UnitConstants::GeV),
+        0 * Acts::UnitConstants::degree /
+            (Acts::UnitConstants::e * Acts::UnitConstants::GeV),
+        0 * Acts::UnitConstants::degree /
+            (Acts::UnitConstants::e * Acts::UnitConstants::GeV),
+        0.1,
+        0 * Acts::UnitConstants::ns /
+            (Acts::UnitConstants::e * Acts::UnitConstants::GeV)};
     /// Inflate initial covariance.
     std::array<double, 6> initialVarInflation = {1., 1., 1., 1., 1., 1.};
     /// Inflate time covariance if no time measurement is available.
@@ -101,10 +117,6 @@ class TrackParamsEstimationAlgorithm final : public IAlgorithm {
  private:
   Config m_cfg;
 
-  /// The track parameters covariance (assumed to be the same for all estimated
-  /// track parameters for the moment)
-  Acts::BoundSquareMatrix m_covariance = Acts::BoundSquareMatrix::Zero();
-
   ReadDataHandle<SimSeedContainer> m_inputSeeds{this, "InputSeeds"};
   ReadDataHandle<ProtoTrackContainer> m_inputTracks{this, "InputTracks"};
 

Examples/Algorithms/TrackFinding/src/TrackParamsEstimationAlgorithm.cpp

@@ -9,6 +9,7 @@
 #include "ActsExamples/TrackFinding/TrackParamsEstimationAlgorithm.hpp"
 
 #include "Acts/Definitions/Algebra.hpp"
+#include "Acts/Definitions/TrackParametrization.hpp"
 #include "Acts/EventData/ParticleHypothesis.hpp"
 #include "Acts/EventData/SourceLink.hpp"
 #include "Acts/Geometry/GeometryIdentifier.hpp"
@@ -32,6 +33,40 @@
 #include <utility>
 #include <vector>
 
+namespace ActsExamples {
+
+namespace {
+
+Acts::BoundSquareMatrix makeInitialCovariance(
+    const TrackParamsEstimationAlgorithm::Config& config,
+    const Acts::BoundVector& params, const SimSpacePoint& sp) {
+  Acts::BoundSquareMatrix result = Acts::BoundSquareMatrix::Zero();
+
+  for (std::size_t i = Acts::eBoundLoc0; i < Acts::eBoundSize; ++i) {
+    double sigma = config.initialSigmas[i];
+
+    // Add momentum dependent uncertainties
+    sigma +=
+        config.initialSimgaQoverPCoefficients[i] * params[Acts::eBoundQOverP];
+
+    double var = sigma * sigma;
+
+    // Inflate the time uncertainty if no time measurement is available
+    if (i == Acts::eBoundTime && !sp.t().has_value()) {
+      var *= config.noTimeVarInflation;
+    }
+
+    // Inflate the initial covariance
+    var *= config.initialVarInflation[i];
+
+    result(i, i) = var;
+  }
+
+  return result;
+}
+
+}  // namespace
+
 ActsExamples::TrackParamsEstimationAlgorithm::TrackParamsEstimationAlgorithm(
     ActsExamples::TrackParamsEstimationAlgorithm::Config cfg,
     Acts::Logging::Level lvl)
@@ -56,12 +91,6 @@ ActsExamples::TrackParamsEstimationAlgorithm::TrackParamsEstimationAlgorithm(
   m_outputTrackParameters.initialize(m_cfg.outputTrackParameters);
   m_outputSeeds.maybeInitialize(m_cfg.outputSeeds);
   m_outputTracks.maybeInitialize(m_cfg.outputProtoTracks);
-
-  // Set up the track parameters covariance (the same for all tracks)
-  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(
@@ -132,11 +161,8 @@ ActsExamples::ProcessCode ActsExamples::TrackParamsEstimationAlgorithm::execute(
 
     const auto& params = optParams.value();
 
-    Acts::BoundSquareMatrix cov = m_covariance;
-    if (!bottomSP->t().has_value()) {
-      // Inflate the time uncertainty if no time measurement is available
-      cov(Acts::eBoundTime, Acts::eBoundTime) *= m_cfg.noTimeVarInflation;
-    }
+    Acts::BoundSquareMatrix cov =
+        makeInitialCovariance(m_cfg, params, *bottomSP);
 
     trackParameters.emplace_back(surface->getSharedPtr(), params, cov,
                                  m_cfg.particleHypothesis);
@@ -161,3 +187,4 @@ ActsExamples::ProcessCode ActsExamples::TrackParamsEstimationAlgorithm::execute(
 
   return ProcessCode::SUCCESS;
 }
+}  // namespace ActsExamples

Examples/Python/src/TrackFinding.cpp

@@ -292,8 +292,8 @@ void addTrackFinding(Context& ctx) {
       ActsExamples::TrackParamsEstimationAlgorithm, mex,
       "TrackParamsEstimationAlgorithm", inputSeeds, inputProtoTracks,
       outputTrackParameters, outputSeeds, outputProtoTracks, trackingGeometry,
-      magneticField, bFieldMin, initialSigmas, initialVarInflation,
-      particleHypothesis);
+      magneticField, bFieldMin, initialSigmas, initialSimgaQoverPCoefficients,
+      initialVarInflation, noTimeVarInflation, particleHypothesis);
 
   {
     using Alg = ActsExamples::TrackFindingAlgorithm;

Examples/Python/tests/root_file_hashes.txt

@@ -44,19 +44,19 @@ test_digitization_example_input[smeared]__particles.root: 5fe7dda2933ee6b9615b06
 test_digitization_example_input[smeared]__measurements.root: 95ef20bcdc349da5bc025faf41cb76855ac087ae75df4f4bf953db5b23927aa5
 test_digitization_example_input[geometric]__particles.root: 5fe7dda2933ee6b9615b064d192322fe07831133cd998e5ed99a3b992b713a10
 test_digitization_example_input[geometric]__measurements.root: fa4729e28fdbbc459400dc3b7cc896c7e0cca047ab889c55c2f148d13b361637
-test_ckf_tracks_example[generic-full_seeding]__trackstates_ckf.root: 6d632ce48408da81e79d494f9cd38d69457046bcb437b6a5df5a166e13e1f3db
-test_ckf_tracks_example[generic-full_seeding]__tracksummary_ckf.root: b37e882b3760706c40893cb2df3de9ca2350907ba9122298cf3527c151426c51
+test_ckf_tracks_example[generic-full_seeding]__trackstates_ckf.root: 2b6778011c20e77ddd3337a1d23dbe19df7e98733ca6f54fd82d7c4361f9eaae
+test_ckf_tracks_example[generic-full_seeding]__tracksummary_ckf.root: 6713eeba4bd96216a7385928ad2f634378e95c4dc6fa42e8b575ac626afa1828
 test_ckf_tracks_example[generic-full_seeding]__performance_seeding_trees.root: 0e0676ffafdb27112fbda50d1cf627859fa745760f98073261dcf6db3f2f991e
-test_ckf_tracks_example[generic-truth_estimated]__trackstates_ckf.root: c03ce2da814c3d85fb288693a62c3903a1d932cf6936044a15fbcc61572df3dd
-test_ckf_tracks_example[generic-truth_estimated]__tracksummary_ckf.root: 269fbbedc2f7aa163d29bfc0bda18312cef953b89e97a638ebb94799a93411d5
+test_ckf_tracks_example[generic-truth_estimated]__trackstates_ckf.root: c5e3e74aa94795fa991cf6e3c4922604fa0ff946c5a4a51a16a299f53534eda6
+test_ckf_tracks_example[generic-truth_estimated]__tracksummary_ckf.root: bd0c382f85a95549fbde5ca977a7311b0405ec669f21131e5d3629e9e24e0bea
 test_ckf_tracks_example[generic-truth_estimated]__performance_seeding.root: 1facb05c066221f6361b61f015cdf0918e94d9f3fce2269ec7b6a4dffeb2bc7e
 test_ckf_tracks_example[generic-truth_smeared]__trackstates_ckf.root: 2a1569f867a155ed0034c6436d5b90550cbc0ae3d3dfa0be17a8a571dd338d7f
 test_ckf_tracks_example[generic-truth_smeared]__tracksummary_ckf.root: f3d8b236e144e139a264964bfa2ba6807cbe977099448e2a7d4159d43d497228
-test_ckf_tracks_example[odd-full_seeding]__trackstates_ckf.root: 49df09e94671a042b8ecb3de394eba273cca7be1d5fcdab6f1127a36f1d6cf44
-test_ckf_tracks_example[odd-full_seeding]__tracksummary_ckf.root: 86b718b20553274933e90599918c0ebd3e97777d1c17dc20ad5facfe011916f2
+test_ckf_tracks_example[odd-full_seeding]__trackstates_ckf.root: 4d3234a5050dfa02c7b699a376ba6b6052e8929b1edc597de42335074be41472
+test_ckf_tracks_example[odd-full_seeding]__tracksummary_ckf.root: e4afdc4d0a90adf60101a036d215735d44627dfa1526b675d6bba600d5998e61
 test_ckf_tracks_example[odd-full_seeding]__performance_seeding_trees.root: 43c58577aafe07645e5660c4f43904efadf91d8cda45c5c04c248bbe0f59814f
-test_ckf_tracks_example[odd-truth_estimated]__trackstates_ckf.root: a278382e1dfe4b3f037f12a62d0bd1ff06a005320ae7691e339b09269b3e6165
-test_ckf_tracks_example[odd-truth_estimated]__tracksummary_ckf.root: 6485e1daf43a2f40e2696c016a0e3dca61a462da7bd6b2a54305fb34a4377ed8
+test_ckf_tracks_example[odd-truth_estimated]__trackstates_ckf.root: 097539b396b5bc7ed43f1e3005354bcbcfba4d24d0244379516dd4d327a1ae95
+test_ckf_tracks_example[odd-truth_estimated]__tracksummary_ckf.root: 9c1c0d99c1211c5505ac822bde7d9a4351a46f15e648a53192a5340630c4779a
 test_ckf_tracks_example[odd-truth_estimated]__performance_seeding.root: 1a36b7017e59f1c08602ef3c2cb0483c51df248f112e3780c66594110719c575
 test_ckf_tracks_example[odd-truth_smeared]__trackstates_ckf.root: fbcad6efe99cc8bda41c49f1ee81fe48c817cbc01940f0577d0f1b93dbc7af55
 test_ckf_tracks_example[odd-truth_smeared]__tracksummary_ckf.root: 3bad4c3191fe62044df6c57dbe2c461943387b1cebba938c8400ef03cff83b72

Examples/Python/tests/test_examples.py

@@ -1223,7 +1223,7 @@ def test_ML_Ambiguity_Solver(tmp_path, assert_root_hash):
     env["ACTS_LOG_FAILURE_THRESHOLD"] = "ERROR"
     try:
         subprocess.check_call(
-            [sys.executable, str(script), "-n5", "--MLSolver"],
+            [sys.executable, str(script), "-n1", "--MLSolver"],
             cwd=tmp_path,
             env=env,
             stderr=subprocess.STDOUT,