Merge branch 'main' into inequality

Core/include/Acts/TrackFinding/CombinatorialKalmanFilter.hpp

@@ -863,79 +863,57 @@ class CombinatorialKalmanFilter {
         auto currentBranch = result.activeBranches.back();
         TrackIndexType prevTip = currentBranch.tipIndex();
 
-        // Add a state only if there already is measurement on this branch
-        if (currentBranch.nMeasurements() > 0) {
-          ACTS_VERBOSE(
-              "Record hole or passive material state on surface since we "
-              "already have measurements on the branch");
-
-          // No source links on surface, add either hole or passive material
-          // TrackState. No storage allocation for uncalibrated/calibrated
-          // measurement and filtered parameter
-          auto stateMask = PM::Predicted | PM::Jacobian;
-
-          // Transport the covariance to a curvilinear surface
-          stepper.transportCovarianceToCurvilinear(state.stepping);
+        // No source links on surface, add either hole or passive material
+        // TrackState. No storage allocation for uncalibrated/calibrated
+        // measurement and filtered parameter
+        auto stateMask = PM::Predicted | PM::Jacobian;
 
-          // Update state and stepper with pre material effects
-          materialInteractor(surface, state, stepper, navigator,
-                             MaterialUpdateStage::PreUpdate);
+        // Transport the covariance to a curvilinear surface
+        stepper.transportCovarianceToCurvilinear(state.stepping);
 
-          // Transport & bind the state to the current surface
-          auto boundStateRes =
-              stepper.boundState(state.stepping, *surface, false);
-          if (!boundStateRes.ok()) {
-            return boundStateRes.error();
-          }
-          auto& boundState = *boundStateRes;
-          auto& [boundParams, jacobian, pathLength] = boundState;
-          boundParams.covariance() = state.stepping.cov;
-
-          // Add a hole or material track state
-          TrackIndexType currentTip = addNonSourcelinkState(
-              stateMask, boundState, result, isSensitive, prevTip);
-          auto nonSourcelinkState =
-              result.trackStates->getTrackState(currentTip);
-          currentBranch.tipIndex() = currentTip;
-          if (isSensitive) {
-            currentBranch.nHoles()++;
-          }
+        // Update state and stepper with pre material effects
+        materialInteractor(surface, state, stepper, navigator,
+                           MaterialUpdateStage::PreUpdate);
 
-          BranchStopperResult branchStopperResult =
-              m_extensions.branchStopper(currentBranch, nonSourcelinkState);
+        // Transport & bind the state to the current surface
+        auto boundStateRes =
+            stepper.boundState(state.stepping, *surface, false);
+        if (!boundStateRes.ok()) {
+          return boundStateRes.error();
+        }
+        auto& boundState = *boundStateRes;
+        auto& [boundParams, jacobian, pathLength] = boundState;
+        boundParams.covariance() = state.stepping.cov;
 
-          // Check the branch
-          if (branchStopperResult == BranchStopperResult::Continue) {
-            // Remembered the active branch and its state
-          } else {
-            // No branch on this surface
-            nBranchesOnSurface = 0;
-            if (branchStopperResult == BranchStopperResult::StopAndKeep) {
-              storeLastActiveBranch(result);
-            }
-            // Remove the branch from list
-            result.activeBranches.pop_back();
-          }
+        // Add a hole or material track state to the multitrajectory
+        TrackIndexType currentTip = addNonSourcelinkState(
+            stateMask, boundState, result, isSensitive, prevTip);
+        auto nonSourcelinkState = result.trackStates->getTrackState(currentTip);
+        currentBranch.tipIndex() = currentTip;
 
-          // Update state and stepper with post material effects
-          materialInteractor(surface, state, stepper, navigator,
-                             MaterialUpdateStage::PostUpdate);
-        } else if (isMaterial) {
-          // Since we did not find a measurement yet we don't create a track
-          // state but still need to do the material update
-          ACTS_VERBOSE("Handle material without creating a track state");
+        if (isSensitive) {
+          currentBranch.nHoles()++;
+        }
 
-          // Transport the covariance to a curvilinear surface
-          stepper.transportCovarianceToCurvilinear(state.stepping);
+        BranchStopperResult branchStopperResult =
+            m_extensions.branchStopper(currentBranch, nonSourcelinkState);
 
-          // Update state and stepper with full material effects
-          materialInteractor(surface, state, stepper, navigator,
-                             MaterialUpdateStage::FullUpdate);
+        // Check the branch
+        if (branchStopperResult == BranchStopperResult::Continue) {
+          // Remembered the active branch and its state
         } else {
-          ACTS_VERBOSE(
-              "Skipping surface as there are no measurements on the branch yet "
-              "and it does not have material");
+          // No branch on this surface
+          nBranchesOnSurface = 0;
+          if (branchStopperResult == BranchStopperResult::StopAndKeep) {
+            storeLastActiveBranch(result);
+          }
+          // Remove the branch from list
+          result.activeBranches.pop_back();
         }
+
+        // Update state and stepper with post material effects
+        materialInteractor(surface, state, stepper, navigator,
+                           MaterialUpdateStage::PostUpdate);
       } else {
         // Neither measurement nor material on surface, this branch is still
         // valid. Count the branch on current surface
@@ -1166,31 +1144,13 @@ class CombinatorialKalmanFilter {
       auto currentBranch = result.activeBranches.back();
       TrackIndexType currentTip = currentBranch.tipIndex();
 
-      ACTS_VERBOSE("Find track with entry index = "
-                   << currentTip << " and there are nMeasurements = "
-                   << currentBranch.nMeasurements()
+      ACTS_VERBOSE("Storing track "
+                   << currentBranch.index() << " with tip index " << currentTip
+                   << ". nMeasurements = " << currentBranch.nMeasurements()
                    << ", nOutliers = " << currentBranch.nOutliers()
-                   << ", nHoles = " << currentBranch.nHoles() << " on track");
-
-      std::optional<TrackStateProxy> lastMeasurement;
-      for (const auto& trackState : currentBranch.trackStatesReversed()) {
-        if (trackState.typeFlags().test(TrackStateFlag::MeasurementFlag) &&
-            !trackState.typeFlags().test(TrackStateFlag::OutlierFlag)) {
-          lastMeasurement = trackState;
-          break;
-        }
-      }
+                   << ", nHoles = " << currentBranch.nHoles());
 
-      if (lastMeasurement.has_value()) {
-        currentBranch.tipIndex() = lastMeasurement->index();
-        result.collectedTracks.push_back(currentBranch);
-        ACTS_VERBOSE("Last measurement found on track with entry index = "
-                     << currentTip << " and measurement index = "
-                     << lastMeasurement->index());
-      } else {
-        ACTS_VERBOSE(
-            "No measurement found on track with entry index = " << currentTip);
-      }
+      result.collectedTracks.push_back(currentBranch);
     }
 
     CombinatorialKalmanFilterExtensions<track_container_t> m_extensions;

Core/include/Acts/Utilities/TrackHelpers.hpp

@@ -368,13 +368,14 @@ Result<void> extrapolateTracksToReferenceSurface(
 
 /// Helper function to calculate a number of track level quantities and store
 /// them on the track itself
-/// @note The input track needs to be mutable, so @c ReadOnly=false
-/// @tparam track_container_t the track container backend
-/// @tparam track_state_container_t the track state container backend
-/// @tparam holder_t the holder type for the track container backends
+/// @tparam track_proxy_t The track proxy type
 /// @param track A mutable track proxy to operate on
 template <TrackProxyConcept track_proxy_t>
-void calculateTrackQuantities(track_proxy_t track) {
+void calculateTrackQuantities(track_proxy_t track)
+  requires(!track_proxy_t::ReadOnly)
+{
+  using ConstTrackStateProxy = typename track_proxy_t::ConstTrackStateProxy;
+
   track.chi2() = 0;
   track.nDoF() = 0;
 
@@ -383,8 +384,8 @@ void calculateTrackQuantities(track_proxy_t track) {
   track.nSharedHits() = 0;
   track.nOutliers() = 0;
 
-  for (const auto &trackState : track.trackStatesReversed()) {
-    auto typeFlags = trackState.typeFlags();
+  for (ConstTrackStateProxy trackState : track.trackStatesReversed()) {
+    ConstTrackStateType typeFlags = trackState.typeFlags();
 
     if (typeFlags.test(Acts::TrackStateFlag::HoleFlag)) {
       track.nHoles()++;
@@ -401,6 +402,97 @@ void calculateTrackQuantities(track_proxy_t track) {
   }
 }
 
+/// Helper function to trim track states from the front of a track
+/// @tparam track_proxy_t the track proxy type
+/// @param track the track to trim
+/// @param trimHoles whether to trim holes
+/// @param trimOutliers whether to trim outliers
+/// @param trimMaterial whether to trim pure material states
+template <TrackProxyConcept track_proxy_t>
+void trimTrackFront(track_proxy_t track, bool trimHoles, bool trimOutliers,
+                    bool trimMaterial)
+  requires(!track_proxy_t::ReadOnly)
+{
+  using TrackStateProxy = typename track_proxy_t::TrackStateProxy;
+
+  // TODO specialize if track is forward linked
+
+  std::optional<TrackStateProxy> front;
+
+  for (TrackStateProxy trackState : track.trackStatesReversed()) {
+    TrackStateType typeFlags = trackState.typeFlags();
+    if (trimHoles && typeFlags.test(TrackStateFlag::HoleFlag)) {
+      continue;
+    }
+    if (trimOutliers && typeFlags.test(TrackStateFlag::OutlierFlag)) {
+      continue;
+    }
+    if (trimMaterial && typeFlags.test(TrackStateFlag::MaterialFlag) &&
+        !typeFlags.test(TrackStateFlag::MeasurementFlag)) {
+      continue;
+    }
+
+    front = trackState;
+  }
+
+  if (front.has_value()) {
+    front.value().previous() = TrackStateProxy::kInvalid;
+  }
+}
+
+/// Helper function to trim track states from the back of a track
+/// @tparam track_proxy_t the track proxy type
+/// @param track the track to trim
+/// @param trimHoles whether to trim holes
+/// @param trimOutliers whether to trim outliers
+/// @param trimMaterial whether to trim pure material states
+template <TrackProxyConcept track_proxy_t>
+void trimTrackBack(track_proxy_t track, bool trimHoles, bool trimOutliers,
+                   bool trimMaterial)
+  requires(!track_proxy_t::ReadOnly)
+{
+  using TrackStateProxy = typename track_proxy_t::TrackStateProxy;
+
+  std::optional<TrackStateProxy> back;
+
+  for (TrackStateProxy trackState : track.trackStatesReversed()) {
+    back = trackState;
+
+    TrackStateType typeFlags = trackState.typeFlags();
+    if (trimHoles && typeFlags.test(TrackStateFlag::HoleFlag)) {
+      continue;
+    }
+    if (trimOutliers && typeFlags.test(TrackStateFlag::OutlierFlag)) {
+      continue;
+    }
+    if (trimMaterial && typeFlags.test(TrackStateFlag::MaterialFlag) &&
+        !typeFlags.test(TrackStateFlag::MeasurementFlag)) {
+      continue;
+    }
+
+    break;
+  }
+
+  if (back.has_value()) {
+    track.tipIndex() = back.value().index();
+  }
+}
+
+/// Helper function to trim track states from the front and back of a track
+/// @tparam track_proxy_t the track proxy type
+/// @param track the track to trim
+/// @param trimHoles whether to trim holes
+/// @param trimOutliers whether to trim outliers
+/// @param trimMaterial whether to trim pure material states
+template <TrackProxyConcept track_proxy_t>
+void trimTrack(track_proxy_t track, bool trimHoles, bool trimOutliers,
+               bool trimMaterial)
+  requires(!track_proxy_t::ReadOnly)
+{
+  trimTrackFront(track, trimHoles, trimOutliers, trimMaterial);
+  trimTrackBack(track, trimHoles, trimOutliers, trimMaterial);
+}
+
 }  // namespace Acts
 
 namespace std {

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

@@ -126,6 +126,8 @@ class TrackFindingAlgorithm final : public IAlgorithm {
     bool stayOnSeed = false;
     /// Compute shared hit information
     bool computeSharedHits = false;
+    /// Whether to trim the tracks
+    bool trimTracks = true;
 
     // Pixel and strip volume ids to be used for maxPixel/StripHoles cuts
     std::set<Acts::GeometryIdentifier::Value> pixelVolumes;

Examples/Algorithms/TrackFinding/src/TrackFindingAlgorithm.cpp

@@ -419,6 +419,10 @@ ProcessCode TrackFindingAlgorithm::execute(const AlgorithmContext& ctx) const {
       }
     });
 
+    // trim the track if requested
+    if (m_cfg.trimTracks) {
+      Acts::trimTrack(track, true, true, true);
+    }
     Acts::calculateTrackQuantities(track);
 
     if (m_trackSelector.has_value() && !m_trackSelector->isValidTrack(track)) {

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

@@ -146,8 +146,9 @@
         "stripVolumes",
         "maxPixelHoles",
         "maxStripHoles",
+        "trimTracks",
     ],
-    defaults=[15.0, 25.0, 10, None, None, None, None, None, None, None],
+    defaults=[15.0, 25.0, 10, None, None, None, None, None, None, None, None],
 )
 
 AmbiguityResolutionConfig = namedtuple(
@@ -1534,6 +1535,7 @@ def addCKFTracks(
             stripVolumes=ckfConfig.stripVolumes,
             maxPixelHoles=ckfConfig.maxPixelHoles,
             maxStripHoles=ckfConfig.maxStripHoles,
+            trimTracks=ckfConfig.trimTracks,
         ),
     )
     s.addAlgorithm(trackFinder)

Examples/Python/src/TrackFinding.cpp

@@ -341,6 +341,7 @@ void addTrackFinding(Context& ctx) {
     ACTS_PYTHON_MEMBER(stripVolumes);
     ACTS_PYTHON_MEMBER(maxPixelHoles);
     ACTS_PYTHON_MEMBER(maxStripHoles);
+    ACTS_PYTHON_MEMBER(trimTracks);
     ACTS_PYTHON_STRUCT_END();
   }
 

Examples/Python/tests/root_file_hashes.txt

@@ -41,11 +41,11 @@ test_ckf_tracks_example[generic-truth_estimated]__tracksummary_ckf.root: 5a973ec
 test_ckf_tracks_example[generic-truth_estimated]__performance_seeding.root: 1facb05c066221f6361b61f015cdf0918e94d9f3fce2269ec7b6a4dffeb2bc7e
 test_ckf_tracks_example[generic-truth_smeared]__trackstates_ckf.root: b10f61d3b68ecc3d1910a17aeadd01a1e23b31b6418935809d1f12e95eac607a
 test_ckf_tracks_example[generic-truth_smeared]__tracksummary_ckf.root: 69733ede1fc08370b5c0d0535f274b59bb51217239bd6645ff62ec9dadaa1f41
-test_ckf_tracks_example[odd-full_seeding]__trackstates_ckf.root: 6b9b97545e9c4c833547adcd54c78235165e7e184d52548ecb5f16c421f41b37
-test_ckf_tracks_example[odd-full_seeding]__tracksummary_ckf.root: 5144c31f9f5281f1c287de978035c2584e9c7b1d46aece6afaecd5c1c86b0020
+test_ckf_tracks_example[odd-full_seeding]__trackstates_ckf.root: 33a68b20005a92b88bf7d01179a4b4fe3c1476a2019176904356924c0ad68f4b
+test_ckf_tracks_example[odd-full_seeding]__tracksummary_ckf.root: 749a53533fed3d550c2e37f504bac5a8ab5c5cdf60d5055dcd545aa44613601a
 test_ckf_tracks_example[odd-full_seeding]__performance_seeding_trees.root: 43c58577aafe07645e5660c4f43904efadf91d8cda45c5c04c248bbe0f59814f
 test_ckf_tracks_example[odd-truth_estimated]__trackstates_ckf.root: 7750f58b970c79dc3c937482e2a6e4576bc5d495fd6e4576d63ac2882d8283d4
-test_ckf_tracks_example[odd-truth_estimated]__tracksummary_ckf.root: 3754817f9cf0972b92147e80bf36b2dc79a9c96e161b236a483a2104d2e33cd6
+test_ckf_tracks_example[odd-truth_estimated]__tracksummary_ckf.root: f5d5d5521e367dd2e26365b896b079279a5a8b97b7e11d38eb0eb317740ad4dd
 test_ckf_tracks_example[odd-truth_estimated]__performance_seeding.root: 1a36b7017e59f1c08602ef3c2cb0483c51df248f112e3780c66594110719c575
 test_ckf_tracks_example[odd-truth_smeared]__trackstates_ckf.root: 33398059bf968f7279d4cc706f3b914fcf6f010ae82c43df8875785bda6f7cbe
 test_ckf_tracks_example[odd-truth_smeared]__tracksummary_ckf.root: 3d01335a51fb03c78a174c8b11d473b38ba9c4ed3d6cba201b3925f463411708