From 396e97d5282b96919345a069aef5778dc43a511a Mon Sep 17 00:00:00 2001
From: Tobias Boeckh <tobias.boeckh@cern.ch>
Date: Thu, 9 Nov 2023 15:37:03 +0100
Subject: [PATCH] store track parameters always from upstream to downstream in
 Trk::Track

---
 .../NtupleDumper/src/NtupleDumperAlg.cxx          | 15 ++-------------
 Tracking/Acts/FaserActsKalmanFilter/src/CKF2.cxx  |  8 ++++----
 .../src/CircleFitTrackSeedTool.cxx                |  4 ++--
 .../src/CreateTrkTrackTool.cxx                    |  7 ++++---
 .../src/CreateTrkTrackTool.h                      |  2 +-
 .../src/KalmanFitterTool.cxx                      | 13 +++++--------
 .../FaserActsKalmanFilter/src/KalmanFitterTool.h  |  2 +-
 7 files changed, 19 insertions(+), 32 deletions(-)

diff --git a/PhysicsAnalysis/NtupleDumper/src/NtupleDumperAlg.cxx b/PhysicsAnalysis/NtupleDumper/src/NtupleDumperAlg.cxx
index 8e540ca7..0e88d1e4 100644
--- a/PhysicsAnalysis/NtupleDumper/src/NtupleDumperAlg.cxx
+++ b/PhysicsAnalysis/NtupleDumper/src/NtupleDumperAlg.cxx
@@ -931,19 +931,8 @@ StatusCode NtupleDumperAlg::execute(const EventContext &ctx) const
     m_hitSet.push_back(hitSet.to_ulong());
     if (stationMap.count(1) == 0 || stationMap.count(2) == 0 || stationMap.count(3) == 0) continue;
 
-    const DataVector<const Trk::TrackParameters> *trackParameters = track->trackParameters();
-    const Trk::TrackParameters* upstreamParameters = *std::min_element(
-      trackParameters->begin(), trackParameters->end(),
-      [](const Trk::TrackParameters *lhs, const Trk::TrackParameters *rhs) {
-          return lhs->position().z() < rhs->position().z();
-      }
-    );
-    const Trk::TrackParameters* downstreamParameters = *std::max_element(
-      trackParameters->begin(), trackParameters->end(),
-      [](const Trk::TrackParameters *lhs, const Trk::TrackParameters *rhs) {
-          return lhs->position().z() < rhs->position().z();
-      }
-    );
+    const Trk::TrackParameters* upstreamParameters = track->trackParameters()->front();
+    const Trk::TrackParameters* downstreamParameters = track->trackParameters()->back();
 
     if ((upstreamParameters == nullptr) || (downstreamParameters == nullptr)) continue;
 
diff --git a/Tracking/Acts/FaserActsKalmanFilter/src/CKF2.cxx b/Tracking/Acts/FaserActsKalmanFilter/src/CKF2.cxx
index 81f7a978..57d78d9d 100644
--- a/Tracking/Acts/FaserActsKalmanFilter/src/CKF2.cxx
+++ b/Tracking/Acts/FaserActsKalmanFilter/src/CKF2.cxx
@@ -213,10 +213,10 @@ StatusCode CKF2::execute() {
   //   ATH_MSG_DEBUG("  position: " << params.position(gctx).transpose());
   //   ATH_MSG_DEBUG("  momentum: " << params.momentum().transpose());
   //   ATH_MSG_DEBUG("  charge:   " << params.charge());
-  //   std::unique_ptr<Trk::Track> track = m_createTrkTrackTool->createTrack(gctx, traj);
+  //   std::unique_ptr<Trk::Track> track = m_createTrkTrackTool->createTrack(gctx, traj, m_backwardPropagation);
   //   if (track != nullptr) {
   //     m_numberOfSelectedTracks++;
-  //     std::unique_ptr<Trk::Track> track2 = m_kalmanFitterTool1->fit(ctx, gctx, track.get(), Acts::BoundVector::Zero(), m_isMC, origin);
+  //     std::unique_ptr<Trk::Track> track2 = m_kalmanFitterTool1->fit(ctx, gctx, track.get(), Acts::BoundVector::Zero(), m_isMC);
   //     if (track2) {
   //       outputAllTracks->push_back(std::move(track2));
   //     } else {
@@ -244,11 +244,11 @@ StatusCode CKF2::execute() {
     ATH_MSG_DEBUG("  position: " << params.position(gctx).transpose());
     ATH_MSG_DEBUG("  momentum: " << params.momentum().transpose());
     ATH_MSG_DEBUG("  charge:   " << params.charge());
-    std::unique_ptr<Trk::Track> track = m_createTrkTrackTool->createTrack(gctx, traj);
+    std::unique_ptr<Trk::Track> track = m_createTrkTrackTool->createTrack(gctx, traj, m_backwardPropagation);
     if (track != nullptr) {
       m_numberOfSelectedTracks++;
       std::unique_ptr<Trk::Track> track2 = m_kalmanFitterTool1->fit(
-        ctx, gctx, track.get(), Acts::BoundVector::Zero(), m_isMC, origin, m_backwardPropagation);
+        ctx, gctx, track.get(), Acts::BoundVector::Zero(), m_isMC);
       if (track2) {
         outputTracks->push_back(std::move(track2));
       } else {
diff --git a/Tracking/Acts/FaserActsKalmanFilter/src/CircleFitTrackSeedTool.cxx b/Tracking/Acts/FaserActsKalmanFilter/src/CircleFitTrackSeedTool.cxx
index 4fe45004..51703a12 100644
--- a/Tracking/Acts/FaserActsKalmanFilter/src/CircleFitTrackSeedTool.cxx
+++ b/Tracking/Acts/FaserActsKalmanFilter/src/CircleFitTrackSeedTool.cxx
@@ -334,12 +334,12 @@ double CircleFitTrackSeedTool::Seed::getX(double z) {
 
 void CircleFitTrackSeedTool::Seed::getChi2() {
   chi2 = 0;
-  for (const Acts::Vector3 &pos : positions) {
+  for (const Acts::Vector3 &pos : fakePositions) {
     m_dy = pos.y() - getY(pos.z());
     chi2 += (m_dy * m_dy) / (m_sigma_y * m_sigma_y);
   }
 
-  for (const Acts::Vector3 &pos : positions) {
+  for (const Acts::Vector3 &pos : fakePositions) {
     m_dx = pos.x() - getX(pos.z());
     chi2 += (m_dx * m_dx) / (m_sigma_x * m_sigma_x);
   }
diff --git a/Tracking/Acts/FaserActsKalmanFilter/src/CreateTrkTrackTool.cxx b/Tracking/Acts/FaserActsKalmanFilter/src/CreateTrkTrackTool.cxx
index 99468901..9a9bf115 100644
--- a/Tracking/Acts/FaserActsKalmanFilter/src/CreateTrkTrackTool.cxx
+++ b/Tracking/Acts/FaserActsKalmanFilter/src/CreateTrkTrackTool.cxx
@@ -17,7 +17,7 @@ StatusCode CreateTrkTrackTool::finalize() {
 }
 
 std::unique_ptr<Trk::Track>
-CreateTrkTrackTool::createTrack(const Acts::GeometryContext &gctx, const FaserActsRecMultiTrajectory &traj) const {
+CreateTrkTrackTool::createTrack(const Acts::GeometryContext &gctx, const FaserActsRecMultiTrajectory &traj, bool backwardPropagation) const {
   std::unique_ptr<Trk::Track> newtrack = nullptr;
   DataVector<const Trk::TrackStateOnSurface>* finalTrajectory = new DataVector<const Trk::TrackStateOnSurface>{};
   using ConstTrackStateProxy = Acts::detail_lt::TrackStateProxy<IndexSourceLink, 6, true>;
@@ -69,9 +69,10 @@ CreateTrkTrackTool::createTrack(const Acts::GeometryContext &gctx, const FaserAc
 									      quality,
 									      nullptr,
 									      typePattern);
-      if (perState)
-      {
+      if ((perState) && (!backwardPropagation)) {
 	        finalTrajectory->insert(finalTrajectory->begin(), perState);
+      } else if ((perState) && (backwardPropagation)) {
+	        finalTrajectory->push_back(perState);
       }
     }
     return;
diff --git a/Tracking/Acts/FaserActsKalmanFilter/src/CreateTrkTrackTool.h b/Tracking/Acts/FaserActsKalmanFilter/src/CreateTrkTrackTool.h
index d0edfd5c..5fcf1f8c 100644
--- a/Tracking/Acts/FaserActsKalmanFilter/src/CreateTrkTrackTool.h
+++ b/Tracking/Acts/FaserActsKalmanFilter/src/CreateTrkTrackTool.h
@@ -18,7 +18,7 @@ public:
   virtual ~CreateTrkTrackTool() = default;
   virtual StatusCode initialize() override;
   virtual StatusCode finalize() override;
-  std::unique_ptr<Trk::Track> createTrack(const Acts::GeometryContext &gctx, const FaserActsRecMultiTrajectory &traj) const;
+  std::unique_ptr<Trk::Track> createTrack(const Acts::GeometryContext &gctx, const FaserActsRecMultiTrajectory &traj, bool backwardPropagation=false) const;
   const Trk::TrackParameters* ConvertActsTrackParameterToATLAS(const Acts::BoundTrackParameters &actsParameter, const Acts::GeometryContext& gctx) const;
 private:
   const FaserSCT_ID* m_idHelper {nullptr};
diff --git a/Tracking/Acts/FaserActsKalmanFilter/src/KalmanFitterTool.cxx b/Tracking/Acts/FaserActsKalmanFilter/src/KalmanFitterTool.cxx
index 942d649b..5315404c 100644
--- a/Tracking/Acts/FaserActsKalmanFilter/src/KalmanFitterTool.cxx
+++ b/Tracking/Acts/FaserActsKalmanFilter/src/KalmanFitterTool.cxx
@@ -365,7 +365,7 @@ resi.push_back({local.x(),local.y(),parameter.position(gctx).x(), parameter.posi
 std::unique_ptr<Trk::Track>
 KalmanFitterTool::fit(const EventContext &ctx, const Acts::GeometryContext &gctx,
                       Trk::Track* inputTrack, const Acts::BoundVector& inputVector,
-                      bool isMC, double origin, bool backwardPropagation) const {
+                      bool isMC) const {
   std::unique_ptr<Trk::Track> newTrack = nullptr;
   std::vector<FaserActsRecMultiTrajectory> myTrajectories;
 
@@ -380,6 +380,9 @@ KalmanFitterTool::fit(const EventContext &ctx, const Acts::GeometryContext &gctx
   }
 
 
+  // set the start position 5 mm in front of the first track measurement
+  double origin = inputTrack->trackParameters()->front()->position().z() - 10;
+
   auto pSurface = Acts::Surface::makeShared<Acts::PlaneSurface>(
       Acts::Vector3 {0, 0, origin}, Acts::Vector3{0, 0, -1});
 
@@ -396,16 +399,10 @@ KalmanFitterTool::fit(const EventContext &ctx, const Acts::GeometryContext &gctx
   FaserActsOutlierFinder faserActsOutlierFinder{0};
   faserActsOutlierFinder.cluster_z=-1000000.;
   faserActsOutlierFinder.StateChiSquaredPerNumberDoFCut=10000.;
-  Acts::PropagatorPlainOptions pOptions;
-  if (backwardPropagation) {
-    pOptions.direction = Acts::backward;
-  } else {
-    pOptions.direction = Acts::forward;
-  }
   Acts::KalmanFitterOptions<MeasurementCalibrator, FaserActsOutlierFinder, Acts::VoidReverseFilteringLogic>
       kfOptions(gctx, mfContext, calibContext, MeasurementCalibrator(measurements),
                 faserActsOutlierFinder, Acts::VoidReverseFilteringLogic(), Acts::LoggerWrapper{*m_logger},
-                pOptions, &(*pSurface));
+                Acts::PropagatorPlainOptions(), &(*pSurface));
   kfOptions.multipleScattering = false;
   kfOptions.energyLoss = false;
 
diff --git a/Tracking/Acts/FaserActsKalmanFilter/src/KalmanFitterTool.h b/Tracking/Acts/FaserActsKalmanFilter/src/KalmanFitterTool.h
index 759a456b..4901d324 100644
--- a/Tracking/Acts/FaserActsKalmanFilter/src/KalmanFitterTool.h
+++ b/Tracking/Acts/FaserActsKalmanFilter/src/KalmanFitterTool.h
@@ -90,7 +90,7 @@ public:
   std::unique_ptr<Trk::Track> fit(const EventContext &ctx, const Acts::GeometryContext &gctx,
                                   Trk::Track *inputTrack,
                                   const Acts::BoundVector& inputVector = Acts::BoundVector::Zero(),
-                                  bool isMC=false, double origin=0, bool backwardPropagation=false) const;
+                                  bool isMC=false) const;
   std::vector<TSOS4Residual> getUnbiasedResidual(const EventContext &ctx, const Acts::GeometryContext &gctx,
                                   Trk::Track *inputTrack,
                                   const Acts::BoundVector& inputVector = Acts::BoundVector::Zero(),
-- 
GitLab