Secvertex amvf to be merged into main EICrecon repo

src/algorithms/tracking/SecondaryVertexFinder.h

@@ -0,0 +1,91 @@
+// SPDX-License-Identifier: LGPL-3.0-or-later
+// Copyright (C) 2025 Dongwi H. Dongwi (Bishoy)
+
+#pragma once
+
+#include <Acts/Geometry/GeometryContext.hpp>
+#include <Acts/MagneticField/MagneticFieldContext.hpp>
+#include <Acts/Propagator/EigenStepper.hpp>
+#include <Acts/Vertexing/AdaptiveGridDensityVertexFinder.hpp>
+#include <Acts/Vertexing/AdaptiveMultiVertexFinder.hpp>
+#include <Acts/Vertexing/AdaptiveMultiVertexFitter.hpp>
+#include <Acts/Vertexing/HelicalTrackLinearizer.hpp>
+#include <Acts/Vertexing/ImpactPointEstimator.hpp>
+#include <Acts/Vertexing/Vertex.hpp>
+#include <Acts/Vertexing/VertexingOptions.hpp>
+#include <ActsExamples/EventData/Trajectories.hpp>
+#include <algorithms/algorithm.h>
+#include <edm4eic/ReconstructedParticleCollection.h>
+#include <edm4eic/VertexCollection.h>
+#include <spdlog/logger.h>
+#include <gsl/pointers>
+#include <memory>
+#include <string>
+#include <string_view>
+#include <vector>
+
+#include "ActsGeometryProvider.h"
+#include "DD4hepBField.h"
+#include "SecondaryVertexFinderConfig.h"
+#include "algorithms/interfaces/WithPodConfig.h"
+
+namespace eicrecon {
+
+using SecondaryVertexFinderAlgorithm =
+    algorithms::Algorithm<algorithms::Input<edm4eic::ReconstructedParticleCollection,
+                                            std::vector<const ActsExamples::Trajectories>>,
+                          algorithms::Output<edm4eic::VertexCollection, edm4eic::VertexCollection>>;
+
+class SecondaryVertexFinder : public SecondaryVertexFinderAlgorithm,
+                              public WithPodConfig<eicrecon::SecondaryVertexFinderConfig> {
+public:
+  SecondaryVertexFinder(std::string_view name)
+      : SecondaryVertexFinderAlgorithm{name,
+                                       {"inputReconstructedParticles", "inputActsTrajectories"},
+                                       {"outputPrimaryVertices", "outputSecondaryVertices"},
+                                       ""} {}
+
+  void init() final;
+
+  void process(const Input&, const Output&) const final;
+
+  // FIXME this is not compliant with algorithms interface
+  void applyLogger(std::shared_ptr<spdlog::logger> log) { m_log = log; };
+
+private:
+  // Calculate an initial Primary Vertex
+  void calculatePrimaryVertex(
+      const edm4eic::ReconstructedParticleCollection&,
+      const std::vector<gsl::not_null<const ActsExamples::Trajectories*>>& trajectories,
+      Acts::EigenStepper<>, edm4eic::VertexCollection&) const;
+
+  //Calculate secondary vertex and store secVertex container
+  void calculateSecondaryVertex(
+      const edm4eic::ReconstructedParticleCollection&,
+      const std::vector<gsl::not_null<const ActsExamples::Trajectories*>>& trajectories,
+      Acts::EigenStepper<>, edm4eic::VertexCollection&) const;
+
+  // Functions to be used to check efficacy of sec. vertex
+  void setVertexContainer(std::vector<Acts::Vertex> inputcontainer) {
+    vtx_container = inputcontainer;
+  };
+
+  //set up Impact estimator
+  using ImpactPointEstimator   = Acts::ImpactPointEstimator;
+  using LinearizerSec          = Acts::HelicalTrackLinearizer;
+  using VertexFitterSec        = Acts::AdaptiveMultiVertexFitter;
+  using VertexFinderSec        = Acts::AdaptiveMultiVertexFinder;
+  using VertexFinderOptionsSec = Acts::VertexingOptions;
+  using seedFinder             = Acts::AdaptiveGridDensityVertexFinder;
+
+  std::shared_ptr<spdlog::logger> m_log;
+  std::shared_ptr<const ActsGeometryProvider> m_geoSvc;
+
+  std::shared_ptr<const eicrecon::BField::DD4hepBField> m_BField = nullptr;
+  Acts::GeometryContext m_geoctx;
+  Acts::MagneticFieldContext m_fieldctx;
+  SecondaryVertexFinderConfig m_cfg;
+  std::vector<Acts::Vertex> vtx_container;
+};
+
+} // namespace eicrecon

src/algorithms/tracking/SecondaryVertexFinderConfig.h

@@ -0,0 +1,39 @@
+#pragma once
+
+#include <Acts/Definitions/Units.hpp>
+
+namespace eicrecon {
+
+struct SecondaryVertexFinderConfig {
+  unsigned int maxVertices      = 20;
+  unsigned int maxIterations    = 500;
+  unsigned int maxSecIterations = 1000;
+  //Max z interval used for adding tracks to fit: when adding a new vertex to the multi vertex fit,
+  //only the tracks whose z at PCA is closer to the seeded vertex than tracksMaxZinterval
+  //are added to this new vertex
+  float tracksMaxZinterval    = 3 * Acts::UnitConstants::mm;
+  float tracksMaxZintervalSec = 10 * Acts::UnitConstants::mm;
+
+  // Max chi2 value for which tracks are considered compatible with
+  // the fitted vertex. These tracks are removed from the seedTracks
+  // after the fit has been performed.
+  double maxVertexChi2 = 18.42;
+  bool doFullSplitting = false;
+  // 5 corresponds to a p-value of ~0.92 using `chi2(x=5,ndf=2)`
+  float tracksMaxSignificance      = 6.7;
+  float maxMergeVertexSignificance = 3;
+  float minWeight                  = 1e-04;
+  float maxDistToLinPoint          = 5.5 * Acts::UnitConstants::mm;
+  // Bin extent in z-direction
+  float spatialBinExtent         = 25 * Acts::UnitConstants::um;
+  Acts::Vector4 initialVariances = Acts::Vector4{1e+2, 1e+2, 1e+2, 1e+8};
+  // Bin extent in t-direction
+  float temporalBinExtent          = 19. * Acts::UnitConstants::mm;
+  bool doSmoothing                 = true;
+  bool reassignTracksAfterFirstFit = true;
+  bool useTime                     = false;
+  // Use seed vertex as a constraint for the fit
+  bool useSeedConstraint = true;
+};
+
+} // namespace eicrecon

src/algorithms/tracking/TrackProjector.cc

@@ -43,7 +43,7 @@ void TrackProjector::init() {
 
 void TrackProjector::process(const Input& input, const Output& output) const {
   const auto [acts_trajectories, tracks] = input;
-  auto [track_segments]                  = output;
+  const auto [track_segments]            = output;
 
   debug("Track projector event process. Num of input trajectories: {}",
         std::size(acts_trajectories));

src/factories/tracking/SecondaryVertexFinder_factory.h

@@ -0,0 +1,73 @@
+// SPDX-License-Identifier: LGPL-3.0-or-later
+// Copyright (C) 2025 Dongwi H. Dongwi (Bishoy)
+
+#pragma once
+
+#include <ActsExamples/EventData/Trajectories.hpp>
+#include <JANA/JEvent.h>
+#include <edm4eic/Vertex.h>
+#include <edm4eic/TrackParameters.h>
+#include <memory>
+#include <string>
+#include <utility>
+#include <vector>
+
+#include "algorithms/tracking/SecondaryVertexFinderConfig.h"
+#include "algorithms/tracking/SecondaryVertexFinder.h"
+#include "extensions/jana/JOmniFactory.h"
+
+namespace eicrecon {
+
+class SecondaryVertexFinder_factory
+    : public JOmniFactory<SecondaryVertexFinder_factory, SecondaryVertexFinderConfig> {
+
+private:
+  using AlgoT = eicrecon::SecondaryVertexFinder;
+  std::unique_ptr<AlgoT> m_algo;
+
+  PodioInput<edm4eic::ReconstructedParticle> m_reco_input{this};
+  Input<ActsExamples::Trajectories> m_acts_trajectories_input{this};
+  PodioOutput<edm4eic::Vertex> m_prm_vertices_output{this};
+  PodioOutput<edm4eic::Vertex> m_sec_vertices_output{this};
+
+  ParameterRef<unsigned int> m_maxVertices{this, "maxVertices", config().maxVertices,
+                                           "Maximum num vertices that can be found"};
+  ParameterRef<bool> m_reassignTracksAfterFirstFit{
+      this, "reassignTracksAfterFirstFit", config().reassignTracksAfterFirstFit,
+      "Whether or not to reassign tracks after first fit"};
+  ParameterRef<float> m_tracksMaxZinterval{this, "tracksMaxZinterval", config().tracksMaxZinterval,
+                                           "Max z interval for Acts::AdaptiveMultiVertexFinder."};
+  ParameterRef<unsigned int> m_maxIterations{this, "maxIterations", config().maxIterations,
+                                             "Max iterations for Acts::AdaptiveMultivertexFinder"};
+  ParameterRef<float> m_maxDistToLinPoint{
+      this, "maxDistToLinPoint", config().maxDistToLinPoint,
+      "Max disttance to line point (pca) for Acts::AdaptiveMultivertexFinder"};
+
+  Service<ACTSGeo_service> m_ACTSGeoSvc{this};
+
+public:
+  void Configure() {
+    m_algo = std::make_unique<AlgoT>(this->GetPrefix());
+    m_algo->applyConfig(config());
+    m_algo->applyLogger(logger());
+    m_algo->init();
+  }
+
+  void Process(int32_t, uint64_t) {
+    std::vector<gsl::not_null<const ActsExamples::Trajectories*>> acts_trajectories_input;
+    for (auto acts_traj : m_acts_trajectories_input()) {
+      acts_trajectories_input.push_back(acts_traj);
+    }
+    m_algo->process(
+        {
+            m_reco_input(),
+            acts_trajectories_input,
+        },
+        {
+            m_prm_vertices_output().get(),
+            m_sec_vertices_output().get(),
+        });
+  }
+};
+
+} // namespace eicrecon

src/global/tracking/tracking.cc

@@ -30,6 +30,7 @@
 #include "factories/tracking/AmbiguitySolver_factory.h"
 #include "factories/tracking/CKFTracking_factory.h"
 #include "factories/tracking/IterativeVertexFinder_factory.h"
+#include "factories/tracking/SecondaryVertexFinder_factory.h"
 #include "factories/tracking/TrackParamTruthInit_factory.h"
 #include "factories/tracking/TrackProjector_factory.h"
 #include "factories/tracking/TrackPropagation_factory.h"
@@ -392,6 +393,14 @@ void InitPlugin(JApplication* app) {
   app->Add(new JOmniFactoryGeneratorT<CollectionCollector_factory<edm4eic::Track, true>>(
       "CombinedTracks", {"CentralCKFTracks", "B0TrackerCKFTracks"}, {"CombinedTracks"}, app));
 
+  app->Add(new JOmniFactoryGeneratorT<SecondaryVertexFinder_factory>(
+      "SecondaryTrackVerticesAMVF", {"ReconstructedParticles", "CentralCKFActsTrajectories"},
+      {
+          "PrimaryTrackVerticesAMVF",
+          "SecondaryTrackVerticesAMVF",
+      },
+      {}, app));
+
   app->Add(new JOmniFactoryGeneratorT<
            CollectionCollector_factory<edm4eic::MCRecoTrackParticleAssociation, true>>(
       "CombinedTrackAssociations", {"CentralCKFTrackAssociations", "B0TrackerCKFTrackAssociations"},

src/services/io/podio/JEventProcessorPODIO.cc

@@ -210,6 +210,8 @@ JEventProcessorPODIO::JEventProcessorPODIO() {
       "CentralCKFTracks",
       "CentralCKFTrackAssociations",
       "CentralCKFTrackParameters",
+      "PrimaryTrackVerticesAMVF",
+      "SecondaryTrackVerticesAMVF",
       // tracking properties - true seeding
       "CentralCKFTruthSeededTrajectoriesUnfiltered",
       "CentralCKFTruthSeededTracksUnfiltered",