introduce option to use full vertex covariance for IP estimation and remove...

introduce option to use full vertex covariance for IP estimation and remove unique_ptr for TrackToVertexIPEstimator return value

Core/include/Acts/Vertexing/AdaptiveMultiVertexFinder.hpp

@@ -149,6 +149,12 @@ class AdaptiveMultiVertexFinder {
     // performance. To be further investigated.
     bool refitAfterBadVertex = true;
 
+    // Use the full available vertex covariance information after
+    // seeding for the IP estimation. In original implementation
+    // this is not (!) done, however, this is probably not correct.
+    // So definitely consider setting this to true.
+    bool useVertexCovForIPEstimation = false;
+
   };  // Config struct
 
   /// @brief Constructor used if InputTrack_t type == BoundParameters

Core/include/Acts/Vertexing/AdaptiveMultiVertexFinder.ipp

@@ -192,20 +192,21 @@ auto Acts::AdaptiveMultiVertexFinder<vfitter_t, sfinder_t>::getIPSignificance(
   // After all, the vertex seed does have a non-zero convariance in general and
   // it probably should be used.
   Vertex<InputTrack_t> newVtx = vtx;
-  newVtx.setFullCovariance(SpacePointSymMatrix::Zero());
-
-  double significance = 0.;
+  if (not m_cfg.useVertexCovForIPEstimation) {
+    newVtx.setFullCovariance(SpacePointSymMatrix::Zero());
+  }
 
   auto estRes = m_cfg.ipEstimator.estimate(*track, newVtx);
   if (!estRes.ok()) {
     return estRes.error();
   }
 
-  auto ipas = std::move(*estRes);
+  ImpactParametersAndSigma ipas = *estRes;
 
-  if (ipas->sigmad0 > 0 && ipas->sigmaz0 > 0) {
-    significance = std::sqrt(std::pow(ipas->IPd0 / ipas->sigmad0, 2) +
-                             std::pow(ipas->IPz0 / ipas->sigmaz0, 2));
+  double significance = 0.;
+  if (ipas.sigmad0 > 0 && ipas.sigmaz0 > 0) {
+    significance = std::sqrt(std::pow(ipas.IPd0 / ipas.sigmad0, 2) +
+                             std::pow(ipas.IPz0 / ipas.sigmaz0, 2));
   }
 
   return significance;
@@ -267,8 +268,6 @@ auto Acts::AdaptiveMultiVertexFinder<vfitter_t, sfinder_t>::
       }
     }
     if (nearTrackFound) {
-      // TODO: check athena actualcandidate position here (has not changed?)
-      // TODO: so do I want to change the vtx position here?
       vtx.setFullPosition(SpacePointVector(0., 0., newZ, 0.));
 
       // Update vertex info for current vertex

Core/include/Acts/Vertexing/IterativeVertexFinder.hpp

@@ -16,7 +16,6 @@
 #include "Acts/Vertexing/FullBilloirVertexFitter.hpp"
 #include "Acts/Vertexing/HelicalTrackLinearizer.hpp"
 #include "Acts/Vertexing/ImpactPoint3dEstimator.hpp"
-#include "Acts/Vertexing/TrackToVertexIPEstimator.hpp"
 #include "Acts/Vertexing/Vertex.hpp"
 #include "Acts/Vertexing/VertexFinderOptions.hpp"
 #include "Acts/Vertexing/VertexFitterConcept.hpp"

Core/include/Acts/Vertexing/TrackToVertexIPEstimator.hpp

@@ -77,7 +77,7 @@ class TrackToVertexIPEstimator {
   ///
   /// @param track Track to estimate IP from
   /// @param vtx Vertex the track belongs to
-  Result<std::unique_ptr<ImpactParametersAndSigma>> estimate(
+  Result<ImpactParametersAndSigma> estimate(
       const BoundParameters& track, const Vertex<input_track_t>& vtx) const;
 
  private:

Core/include/Acts/Vertexing/TrackToVertexIPEstimator.ipp

@@ -10,8 +10,7 @@
 
 template <typename input_track_t, typename propagator_t,
           typename propagator_options_t>
-Acts::Result<std::unique_ptr<Acts::ImpactParametersAndSigma>>
-Acts::TrackToVertexIPEstimator<
+Acts::Result<Acts::ImpactParametersAndSigma> Acts::TrackToVertexIPEstimator<
     input_track_t, propagator_t,
     propagator_options_t>::estimate(const BoundParameters& track,
                                     const Vertex<input_track_t>& vtx) const {
@@ -19,11 +18,8 @@ Acts::TrackToVertexIPEstimator<
   // towards
   // the vertex position. By this time the vertex should NOT contain this
   // trajectory anymore
-
-  const Vector3D& lp = vtx.position();
-
   const std::shared_ptr<PerigeeSurface> perigeeSurface =
-      Surface::makeShared<PerigeeSurface>(lp);
+      Surface::makeShared<PerigeeSurface>(vtx.position());
 
   // Do the propagation to linPoint
   auto result =
@@ -47,18 +43,18 @@ Acts::TrackToVertexIPEstimator<
 
   ActsVectorD<2> d0JacXY(-std::sin(phi), std::cos(phi));
 
-  std::unique_ptr<ImpactParametersAndSigma> newIPandSigma =
-      std::make_unique<ImpactParametersAndSigma>();
-  newIPandSigma->IPd0 = d0;
+  ImpactParametersAndSigma newIPandSigma;
+
+  newIPandSigma.IPd0 = d0;
   double d0_PVcontrib = d0JacXY.transpose() * (vrtXYCov * d0JacXY);
   if (d0_PVcontrib >= 0) {
-    newIPandSigma->sigmad0 = std::sqrt(
+    newIPandSigma.sigmad0 = std::sqrt(
         d0_PVcontrib + perigeeCov(ParID_t::eLOC_D0, ParID_t::eLOC_D0));
-    newIPandSigma->PVsigmad0 = std::sqrt(d0_PVcontrib);
+    newIPandSigma.PVsigmad0 = std::sqrt(d0_PVcontrib);
   } else {
-    newIPandSigma->sigmad0 =
+    newIPandSigma.sigmad0 =
         std::sqrt(perigeeCov(ParID_t::eLOC_D0, ParID_t::eLOC_D0));
-    newIPandSigma->PVsigmad0 = 0;
+    newIPandSigma.PVsigmad0 = 0;
   }
 
   ActsSymMatrixD<2> covPerigeeZ0Theta;
@@ -72,32 +68,32 @@ Acts::TrackToVertexIPEstimator<
   ActsVectorD<2> z0JacZ0Theta(std::sin(theta), z0 * std::cos(theta));
 
   if (vtxZZCov >= 0) {
-    newIPandSigma->IPz0SinTheta = z0 * std::sin(theta);
-    newIPandSigma->sigmaz0SinTheta = std::sqrt(
+    newIPandSigma.IPz0SinTheta = z0 * std::sin(theta);
+    newIPandSigma.sigmaz0SinTheta = std::sqrt(
         z0JacZ0Theta.transpose() * (covPerigeeZ0Theta * z0JacZ0Theta) +
         std::sin(theta) * vtxZZCov * std::sin(theta));
 
-    newIPandSigma->PVsigmaz0SinTheta =
+    newIPandSigma.PVsigmaz0SinTheta =
         std::sqrt(std::sin(theta) * vtxZZCov * std::sin(theta));
-    newIPandSigma->IPz0 = z0;
-    newIPandSigma->sigmaz0 =
+    newIPandSigma.IPz0 = z0;
+    newIPandSigma.sigmaz0 =
         std::sqrt(vtxZZCov + perigeeCov(ParID_t::eLOC_Z0, ParID_t::eLOC_Z0));
-    newIPandSigma->PVsigmaz0 = std::sqrt(vtxZZCov);
+    newIPandSigma.PVsigmaz0 = std::sqrt(vtxZZCov);
   } else {
     ACTS_WARNING(
         "Contribution to z0_err from PV is negative: Error in PV "
         "error matrix! Removing contribution from PV");
-    newIPandSigma->IPz0SinTheta = z0 * std::sin(theta);
+    newIPandSigma.IPz0SinTheta = z0 * std::sin(theta);
     double sigma2z0sinTheta =
         (z0JacZ0Theta.transpose() * (covPerigeeZ0Theta * z0JacZ0Theta));
-    newIPandSigma->sigmaz0SinTheta = std::sqrt(sigma2z0sinTheta);
-    newIPandSigma->PVsigmaz0SinTheta = 0;
+    newIPandSigma.sigmaz0SinTheta = std::sqrt(sigma2z0sinTheta);
+    newIPandSigma.PVsigmaz0SinTheta = 0;
 
-    newIPandSigma->IPz0 = z0;
-    newIPandSigma->sigmaz0 =
+    newIPandSigma.IPz0 = z0;
+    newIPandSigma.sigmaz0 =
         std::sqrt(perigeeCov(ParID_t::eLOC_Z0, ParID_t::eLOC_Z0));
-    newIPandSigma->PVsigmaz0 = 0;
+    newIPandSigma.PVsigmaz0 = 0;
   }
 
-  return std::move(newIPandSigma);
+  return newIPandSigma;
 }

Core/include/Acts/Vertexing/ZScanVertexFinder.ipp

@@ -27,25 +27,25 @@ auto Acts::ZScanVertexFinder<vfitter_t>::find(
     const BoundParameters& params = m_extractParameters(*iTrk);
 
     std::pair<double, double> z0AndWeight;
-    std::unique_ptr<ImpactParametersAndSigma> ipas = nullptr;
+    ImpactParametersAndSigma ipas;
     if (useConstraint &&
         vFinderOptions.vertexConstraint.covariance()(0, 0) != 0) {
       auto estRes =
           m_cfg.ipEstimator.estimate(params, vFinderOptions.vertexConstraint);
       if (estRes.ok()) {
-        ipas = std::move(*estRes);
+        ipas = *estRes;
       } else {
         return estRes.error();
       }
     }
 
-    if (ipas != nullptr && ipas->sigmad0 > 0) {
+    if (ipas.sigmad0 > 0) {
       // calculate z0
       z0AndWeight.first =
-          ipas->IPz0 + vFinderOptions.vertexConstraint.position().z();
+          ipas.IPz0 + vFinderOptions.vertexConstraint.position().z();
 
       // calculate chi2 of IP
-      double chi2IP = std::pow(ipas->IPd0 / ipas->sigmad0, 2);
+      double chi2IP = std::pow(ipas.IPd0 / ipas.sigmad0, 2);
 
       if (!m_cfg.disableAllWeights) {
         z0AndWeight.second =

Tests/UnitTests/Core/Vertexing/TrackToVertexIPEstimatorTests.cpp

@@ -144,10 +144,10 @@ BOOST_AUTO_TEST_CASE(track_to_vertex_ip_estimator_test) {
         BoundParameters(tgContext, std::move(covMat), paramVec, perigeeSurface);
 
     // Check if IP are retrieved
-    std::unique_ptr<ImpactParametersAndSigma> output =
+    ImpactParametersAndSigma output =
         ipEst.estimate(track, myConstraint).value();
-    BOOST_CHECK_NE(output->IPd0, 0.);
-    BOOST_CHECK_NE(output->IPz0, 0.);
+    BOOST_CHECK_NE(output.IPd0, 0.);
+    BOOST_CHECK_NE(output.IPz0, 0.);
   }
 }
 }  // namespace Test