refactor: Align IVF and AMVF in examples (#2144)

IVF and AMVF use different seeders and the code diverged quite a bit in the examples. In this PR I try to realign them

Core/include/Acts/Vertexing/IterativeVertexFinder.hpp

@@ -94,7 +94,7 @@ class IterativeVertexFinder {
     /// ImpactPointEstimator
     IPEstimator ipEst;
 
-    /// Vertex finder configuration variables
+    // Vertex finder configuration variables
     bool useBeamConstraint = false;
     double significanceCutSeeding = 10;
     double maximumChi2cutForSeeding = 36.;
@@ -105,6 +105,9 @@ class IterativeVertexFinder {
     bool doMaxTracksCut = false;
     int maxTracks = 5000;
     double cutOffTrackWeight = 0.01;
+    /// If `reassignTracksAfterFirstFit` is set this threshold will be used to
+    /// decide if a track should be checked for reassignment to other vertices
+    double cutOffTrackWeightReassign = 1;
   };
 
   /// State struct

Core/include/Acts/Vertexing/IterativeVertexFinder.ipp

@@ -22,19 +22,18 @@ auto Acts::IterativeVertexFinder<vfitter_t, sfinder_t>::find(
   int nInterations = 0;
   // begin iterating
   while (seedTracks.size() > 1 && nInterations < m_cfg.maxVertices) {
-    /// Begin seeding
+    /// Do seeding
     auto seedRes = getVertexSeed(seedTracks, vertexingOptions);
+
     if (!seedRes.ok()) {
       return seedRes.error();
     }
-    // retrieve the seed vertex as the last element in
-    // the seed vertexCollection
-    Vertex<InputTrack_t>& seedVertex = *seedRes;
+
+    const auto& seedVertex = *seedRes;
 
     if (seedVertex.fullPosition()[eZ] ==
         vertexingOptions.vertexConstraint.position().z()) {
-      ACTS_DEBUG(
-          "No seed found anymore. Break and stop primary vertex finding.");
+      ACTS_DEBUG("No more seed found. Break and stop primary vertex finding.");
       break;
     }
 
@@ -87,7 +86,8 @@ auto Acts::IterativeVertexFinder<vfitter_t, sfinder_t>::find(
       }
     }
     /// End vertex fit
-    ACTS_DEBUG("Vertex position after fit: " << currentVertex.fullPosition());
+    ACTS_DEBUG("Vertex position after fit: "
+               << currentVertex.fullPosition().transpose());
 
     // Number degrees of freedom
     double ndf = currentVertex.fitQuality().second;
@@ -150,7 +150,6 @@ auto Acts::IterativeVertexFinder<vfitter_t, sfinder_t>::find(
     }
 
     nInterations++;
-
   }  // end while loop
 
   return vertexCollection;
@@ -161,27 +160,36 @@ auto Acts::IterativeVertexFinder<vfitter_t, sfinder_t>::getVertexSeed(
     const std::vector<const InputTrack_t*>& seedTracks,
     const VertexingOptions<InputTrack_t>& vertexingOptions) const
     -> Result<Vertex<InputTrack_t>> {
-  typename sfinder_t::State state;
-  auto res = m_cfg.seedFinder.find(seedTracks, vertexingOptions, state);
-  if (res.ok()) {
-    auto vertexCollection = *res;
-    if (vertexCollection.empty()) {
-      ACTS_DEBUG(
-          "No seed found. Number of input tracks: " << seedTracks.size());
-      return VertexingError::SeedingError;
-    }
-    // current seed is last element in collection
-    Vertex<InputTrack_t> seedVertex = vertexCollection.back();
-    ACTS_DEBUG("Seed found at position: ("
-               << seedVertex.fullPosition()[eX] << ", "
-               << seedVertex.fullPosition()[eY] << ", "
-               << seedVertex.fullPosition()[eZ] << ", " << seedVertex.time()
-               << "). Number of input tracks: " << seedTracks.size());
-    return seedVertex;
-  } else {
-    ACTS_DEBUG("No seed found. Number of input tracks: " << seedTracks.size());
+  typename sfinder_t::State finderState;
+  auto res = m_cfg.seedFinder.find(seedTracks, vertexingOptions, finderState);
+
+  if (!res.ok()) {
+    ACTS_DEBUG("Seeding error: internal. Number of input tracks: "
+               << seedTracks.size());
     return VertexingError::SeedingError;
   }
+
+  const auto& vertexCollection = *res;
+
+  if (vertexCollection.empty()) {
+    ACTS_DEBUG("Seeding error: no seeds. Number of input tracks: "
+               << seedTracks.size());
+    return VertexingError::SeedingError;
+  }
+
+  ACTS_DEBUG("Found " << vertexCollection.size() << " seeds");
+
+  // retrieve the seed vertex as the last element in
+  // the seed vertexCollection
+  Vertex<InputTrack_t> seedVertex = vertexCollection.back();
+
+  ACTS_DEBUG("Considering seed at position: ("
+             << seedVertex.fullPosition()[eX] << ", "
+             << seedVertex.fullPosition()[eY] << ", "
+             << seedVertex.fullPosition()[eZ] << ", " << seedVertex.time()
+             << "). Number of input tracks: " << seedTracks.size());
+
+  return seedVertex;
 }
 
 template <typename vfitter_t, typename sfinder_t>
@@ -231,7 +239,7 @@ Acts::IterativeVertexFinder<vfitter_t, sfinder_t>::getCompatibility(
        (vertex.fullCovariance() * linTrack.positionJacobian.transpose()))
           .template block<2, 2>(0, 0);
   weightReduced += errorVertexReduced;
-  weightReduced = weightReduced.inverse();
+  weightReduced = weightReduced.inverse().eval();
 
   // Calculate compatibility / chi2
   Vector2 trackParameters2D =
@@ -392,8 +400,8 @@ Acts::IterativeVertexFinder<vfitter_t, sfinder_t>::fillPerigeesToFit(
       }
 
       if (*distanceRes / error < m_cfg.significanceCutSeeding) {
-        if (count % m_cfg.splitVerticesTrkInvFraction == 0 ||
-            !m_cfg.createSplitVertices) {
+        if (!m_cfg.createSplitVertices ||
+            count % m_cfg.splitVerticesTrkInvFraction == 0) {
           perigeesToFitOut.push_back(sTrack);
           ++count;
         } else {
@@ -429,7 +437,7 @@ Acts::IterativeVertexFinder<vfitter_t, sfinder_t>::reassignTracksToNewVertex(
     for (auto tracksIter = tracksBegin; tracksIter != tracksEnd;) {
       // consider only tracks that are not too tightly assigned to other
       // vertex
-      if (tracksIter->trackWeight > m_cfg.cutOffTrackWeight) {
+      if (tracksIter->trackWeight > m_cfg.cutOffTrackWeightReassign) {
         tracksIter++;
         continue;
       }
@@ -452,8 +460,8 @@ Acts::IterativeVertexFinder<vfitter_t, sfinder_t>::reassignTracksToNewVertex(
       }
       double chi2OldVtx = *resultOld;
 
-      ACTS_DEBUG("Compatibility to new vertex: " << chi2NewVtx);
-      ACTS_DEBUG("Compatibility to old vertex: " << chi2OldVtx);
+      ACTS_DEBUG("Compatibility to new vs old vertex: " << chi2NewVtx << " vs "
+                                                        << chi2OldVtx);
 
       if (chi2NewVtx < chi2OldVtx) {
         perigeesToFit.push_back(tracksIter->originalParams);
@@ -525,7 +533,7 @@ Acts::IterativeVertexFinder<vfitter_t, sfinder_t>::reassignTracksToNewVertex(
   if (!isGoodVertex) {
     removeAllTracks(perigeesToFit, seedTracks);
 
-    ACTS_DEBUG("Going to new  iteration with "
+    ACTS_DEBUG("Going to new iteration with "
                << seedTracks.size() << "seed tracks after BAD vertex.");
   }
 

Core/include/Acts/Vertexing/ZScanVertexFinder.hpp

@@ -86,7 +86,6 @@ class ZScanVertexFinder {
                         getDefaultLogger("ZScanVertexFinder", Logging::INFO))
       : m_cfg(std::move(cfg)),
         m_extractParameters([](T params) { return params; }),
-
         m_logger(std::move(logger)) {}
 
   /// @brief Constructor for user-defined InputTrack_t type =!

Examples/Algorithms/Vertexing/include/ActsExamples/Vertexing/AdaptiveMultiVertexFinderAlgorithm.hpp

@@ -36,11 +36,13 @@ class AdaptiveMultiVertexFinderAlgorithm final : public IAlgorithm {
   using Linearizer = Acts::HelicalTrackLinearizer<Propagator>;
   using Fitter =
       Acts::AdaptiveMultiVertexFitter<Acts::BoundTrackParameters, Linearizer>;
-  using SeedFinder = Acts::TrackDensityVertexFinder<
+  using Seeder = Acts::TrackDensityVertexFinder<
       Fitter, Acts::GaussianTrackDensity<Acts::BoundTrackParameters>>;
-  using Finder = Acts::AdaptiveMultiVertexFinder<Fitter, SeedFinder>;
+  using Finder = Acts::AdaptiveMultiVertexFinder<Fitter, Seeder>;
+  using Options = Acts::VertexingOptions<Acts::BoundTrackParameters>;
 
-  using VertexCollection = std::vector<Acts::Vertex<Fitter::InputTrack_t>>;
+  using VertexCollection =
+      std::vector<Acts::Vertex<Acts::BoundTrackParameters>>;
 
   struct Config {
     /// Optional. Input track parameters collection

Examples/Algorithms/Vertexing/include/ActsExamples/Vertexing/IterativeVertexFinderAlgorithm.hpp

@@ -16,7 +16,7 @@
 #include "Acts/Vertexing/HelicalTrackLinearizer.hpp"
 #include "Acts/Vertexing/ImpactPointEstimator.hpp"
 #include "Acts/Vertexing/IterativeVertexFinder.hpp"
-#include "Acts/Vertexing/ZScanVertexFinder.hpp"
+#include "Acts/Vertexing/TrackDensityVertexFinder.hpp"
 #include "ActsExamples/EventData/ProtoVertex.hpp"
 #include "ActsExamples/EventData/Track.hpp"
 #include "ActsExamples/EventData/Trajectories.hpp"
@@ -30,16 +30,15 @@ namespace ActsExamples {
 class IterativeVertexFinderAlgorithm final : public IAlgorithm {
  public:
   using Propagator = Acts::Propagator<Acts::EigenStepper<>>;
-  using PropagatorOptions = Acts::PropagatorOptions<>;
+  using IPEstimator =
+      Acts::ImpactPointEstimator<Acts::BoundTrackParameters, Propagator>;
   using Linearizer = Acts::HelicalTrackLinearizer<Propagator>;
-  using VertexFitter =
+  using Fitter =
       Acts::FullBilloirVertexFitter<Acts::BoundTrackParameters, Linearizer>;
-  using ImpactPointEstimator =
-      Acts::ImpactPointEstimator<Acts::BoundTrackParameters, Propagator>;
-  using VertexSeeder = Acts::ZScanVertexFinder<VertexFitter>;
-  using VertexFinder = Acts::IterativeVertexFinder<VertexFitter, VertexSeeder>;
-  using VertexFinderOptions =
-      Acts::VertexingOptions<Acts::BoundTrackParameters>;
+  using Seeder = Acts::TrackDensityVertexFinder<
+      Fitter, Acts::GaussianTrackDensity<Acts::BoundTrackParameters>>;
+  using Finder = Acts::IterativeVertexFinder<Fitter, Seeder>;
+  using Options = Acts::VertexingOptions<Acts::BoundTrackParameters>;
 
   using VertexCollection =
       std::vector<Acts::Vertex<Acts::BoundTrackParameters>>;

Examples/Algorithms/Vertexing/src/AdaptiveMultiVertexFinderAlgorithm.cpp

@@ -93,7 +93,7 @@ ActsExamples::AdaptiveMultiVertexFinderAlgorithm::execute(
   Fitter fitter(fitterCfg, logger().cloneWithSuffix("AMVFitter"));
 
   // Set up the vertex seed finder
-  SeedFinder seedFinder;
+  Seeder seedFinder;
 
   Finder::Config finderConfig(std::move(fitter), seedFinder, ipEstimator,
                               std::move(linearizer), m_cfg.bField);
@@ -133,8 +133,7 @@ ActsExamples::AdaptiveMultiVertexFinderAlgorithm::execute(
   Finder::State state;
 
   // Default vertexing options, this is where e.g. a constraint could be set
-  using VertexingOptions = Acts::VertexingOptions<Acts::BoundTrackParameters>;
-  VertexingOptions finderOpts(ctx.geoContext, ctx.magFieldContext);
+  Options finderOpts(ctx.geoContext, ctx.magFieldContext);
 
   VertexCollection vertices;
 

Examples/Algorithms/Vertexing/src/IterativeVertexFinderAlgorithm.cpp

@@ -90,26 +90,24 @@ ActsExamples::ProcessCode ActsExamples::IterativeVertexFinderAlgorithm::execute(
   // Set up propagator with void navigator
   auto propagator = std::make_shared<Propagator>(
       stepper, Acts::detail::VoidNavigator{}, logger().cloneWithSuffix("Prop"));
-  PropagatorOptions propagatorOpts(ctx.geoContext, ctx.magFieldContext);
   // Setup the vertex fitter
-  VertexFitter::Config vertexFitterCfg;
-  VertexFitter vertexFitter(vertexFitterCfg);
+  Fitter::Config vertexFitterCfg;
+  Fitter vertexFitter(vertexFitterCfg);
   // Setup the track linearizer
   Linearizer::Config linearizerCfg(m_cfg.bField, propagator);
   Linearizer linearizer(linearizerCfg, logger().cloneWithSuffix("HelLin"));
   // Setup the seed finder
-  ImpactPointEstimator::Config ipEstCfg(m_cfg.bField, propagator);
-  ImpactPointEstimator ipEst(ipEstCfg);
-  VertexSeeder::Config seederCfg(ipEst);
-  VertexSeeder seeder(seederCfg);
+  IPEstimator::Config ipEstCfg(m_cfg.bField, propagator);
+  IPEstimator ipEst(ipEstCfg);
+  Seeder seeder;
   // Set up the actual vertex finder
-  VertexFinder::Config finderCfg(vertexFitter, std::move(linearizer),
-                                 std::move(seeder), ipEst);
+  Finder::Config finderCfg(vertexFitter, std::move(linearizer),
+                           std::move(seeder), ipEst);
   finderCfg.maxVertices = 200;
-  finderCfg.reassignTracksAfterFirstFit = true;
-  VertexFinder finder(finderCfg);
-  VertexFinder::State state(*m_cfg.bField, ctx.magFieldContext);
-  VertexFinderOptions finderOpts(ctx.geoContext, ctx.magFieldContext);
+  finderCfg.reassignTracksAfterFirstFit = false;
+  Finder finder(finderCfg, logger().cloneWithSuffix("Finder"));
+  Finder::State state(*m_cfg.bField, ctx.magFieldContext);
+  Options finderOpts(ctx.geoContext, ctx.magFieldContext);
 
   // find vertices and measure elapsed time
   auto t1 = std::chrono::high_resolution_clock::now();