feat: change the behaviour of enableCutsForSortedSP in seedFinder (#1213

) This PR splits `enableCutsForSortedSP` into two booleans: Instead of using `enableCutsForSortedSP` to enable both the cotTheta sorting of SPs and the cut based on the sorting, this PR makes `enableCutsForSortedSP ` enable only the sorting and `skipPreviousTopSP` is added for the additional cut to skip top SPs if cotThetaBB < cotThetaT For ITk strips the cuts are not applied.
acts-project · Apr 28, 2022 · b2c7017 · b2c7017
1 parent fc2ee7f
commit b2c7017
Show file tree

Hide file tree

Showing 8 changed files with 43 additions and 68 deletions.
diff --git a/Core/include/Acts/Seeding/SeedFinderOrthogonal.ipp b/Core/include/Acts/Seeding/SeedFinderOrthogonal.ipp
@@ -219,8 +219,8 @@ void SeedFinderOrthogonal<external_spacepoint_t>::filterCandidates(
   std::vector<LinCircle> linCircleTop;
   linCircleTop.reserve(top.size());
 
-  transformCoordinates(bottom, middle, true, false, linCircleBottom);
-  transformCoordinates(top, middle, false, false, linCircleTop);
+  transformCoordinates(bottom, middle, true, linCircleBottom);
+  transformCoordinates(top, middle, false, linCircleTop);
 
   std::vector<float> tanLM;
   std::vector<float> tanMT;

diff --git a/Core/include/Acts/Seeding/SeedFinderUtils.hpp b/Core/include/Acts/Seeding/SeedFinderUtils.hpp
@@ -53,18 +53,16 @@ LinCircle transformCoordinates(external_spacepoint_t& sp,
 /// @param[in] vec The list of bottom or top spacepoints
 /// @param[in] spM The middle spacepoint.
 /// @param[in] bottom Should be true if vec are bottom spacepoints.
-/// @param[in] enableCutsForSortedSP enables sorting of cotTheta.
 /// @param[out] linCircleVec The output vector to write to.
 template <typename external_spacepoint_t>
 void transformCoordinates(
     std::vector<InternalSpacePoint<external_spacepoint_t>*>& vec,
     InternalSpacePoint<external_spacepoint_t>& spM, bool bottom,
-    bool enableCutsForSortedSP, std::vector<LinCircle>& linCircleVec);
+    std::vector<LinCircle>& linCircleVec);
 
 template <typename external_spacepoint_t, typename callable_t>
 void transformCoordinates(std::vector<external_spacepoint_t*>& vec,
                           external_spacepoint_t& spM, bool bottom,
-                          bool enableCutsForSortedSP,
                           std::vector<LinCircle>& linCircleVec,
                           callable_t&& extractFunction);
 }  // namespace Acts

diff --git a/Core/include/Acts/Seeding/SeedFinderUtils.ipp b/Core/include/Acts/Seeding/SeedFinderUtils.ipp
@@ -62,7 +62,7 @@ template <typename external_spacepoint_t>
 void transformCoordinates(
     std::vector<InternalSpacePoint<external_spacepoint_t>*>& vec,
     InternalSpacePoint<external_spacepoint_t>& spM, bool bottom,
-    bool enableCutsForSortedSP, std::vector<LinCircle>& linCircleVec) {
+    std::vector<LinCircle>& linCircleVec) {
   auto extractFunction =
       [](const InternalSpacePoint<external_spacepoint_t>& obj)
       -> std::array<float, 6> {
@@ -72,13 +72,12 @@ void transformCoordinates(
   };
 
   return transformCoordinates<InternalSpacePoint<external_spacepoint_t>>(
-      vec, spM, bottom, enableCutsForSortedSP, linCircleVec, extractFunction);
+      vec, spM, bottom, linCircleVec, extractFunction);
 }
 
 template <typename external_spacepoint_t, typename callable_t>
 void transformCoordinates(std::vector<external_spacepoint_t*>& vec,
                           external_spacepoint_t& spM, bool bottom,
-                          bool enableCutsForSortedSP,
                           std::vector<LinCircle>& linCircleVec,
                           callable_t&& extractFunction) {
   auto [xM, yM, zM, rM, varianceRM, varianceZM] = extractFunction(spM);
@@ -132,15 +131,13 @@ void transformCoordinates(std::vector<external_spacepoint_t*>& vec,
     sp->setCotTheta(cot_theta);
   }
   // sort the SP in order of cotTheta
-  if (enableCutsForSortedSP) {
-    std::sort(vec.begin(), vec.end(),
-              [](external_spacepoint_t* a, external_spacepoint_t* b) -> bool {
-                return (a->cotTheta() < b->cotTheta());
-              });
-    std::sort(linCircleVec.begin(), linCircleVec.end(),
-              [](const LinCircle& a, const LinCircle& b) -> bool {
-                return (a.cotTheta < b.cotTheta);
-              });
-  }
+  std::sort(vec.begin(), vec.end(),
+            [](external_spacepoint_t* a, external_spacepoint_t* b) -> bool {
+              return (a->cotTheta() < b->cotTheta());
+            });
+  std::sort(linCircleVec.begin(), linCircleVec.end(),
+            [](const LinCircle& a, const LinCircle& b) -> bool {
+              return (a.cotTheta < b.cotTheta);
+            });
 }
 }  // namespace Acts
diff --git a/Core/include/Acts/Seeding/Seedfinder.ipp b/Core/include/Acts/Seeding/Seedfinder.ipp
@@ -239,9 +239,8 @@ void Seedfinder<external_spacepoint_t, platform_t>::createSeedsForGroup(
     state.linCircleTop.clear();
 
     transformCoordinates(state.compatBottomSP, *spM, true,
-                         m_config.enableCutsForSortedSP, state.linCircleBottom);
-    transformCoordinates(state.compatTopSP, *spM, false,
-                         m_config.enableCutsForSortedSP, state.linCircleTop);
+                         state.linCircleBottom);
+    transformCoordinates(state.compatTopSP, *spM, false, state.linCircleTop);
 
     state.topSpVec.clear();
     state.curvatures.clear();
@@ -293,39 +292,27 @@ void Seedfinder<external_spacepoint_t, platform_t>::createSeedsForGroup(
 
         float deltaCotTheta = cotThetaB - lt.cotTheta;
         float deltaCotTheta2 = deltaCotTheta * deltaCotTheta;
-        float error;
-        float dCotThetaMinusError2;
-        if (m_config.enableCutsForSortedSP) {
-          // if the error is larger than the difference in theta, no need to
-          // compare with scattering
-          if (deltaCotTheta2 - error2 > scatteringInRegion2) {
+        // Apply a cut on the compatibility between the r-z slope of the two
+        // seed segments. This is done by comparing the squared difference
+        // between slopes, and comparing to the squared uncertainty in this
+        // difference - we keep a seed if the difference is compatible within
+        // the assumed uncertainties. The uncertainties get contribution from
+        // the  space-point-related squared error (error2) and a scattering term
+        // calculated assuming the minimum pt we expect to reconstruct
+        // (scatteringInRegion2). This assumes gaussian error propagation which
+        // allows just adding the two errors if they are uncorrelated (which is
+        // fair for scattering and measurement uncertainties)
+        if (deltaCotTheta2 > (error2 + scatteringInRegion2)) {
+          // additional cut to skip top SPs when producing triplets
+          if (m_config.skipPreviousTopSP) {
             // break if cotThetaB < lt.cotTheta because the SP are sorted by
             // cotTheta
             if (cotThetaB - lt.cotTheta < 0) {
               break;
             }
-            // since cotThetaB > lt.cotTheta and the SP are sorted by cotTheta,
-            // the next bottom SP is expected to have cotThetaB > lt.cotTheta as
-            // well and deltaCotTheta2 - error2 > sigmaSquaredScatteringMinPt
             t0 = t + 1;
-            continue;
-          }
-        } else {
-          if (deltaCotTheta2 - error2 > 0) {
-            deltaCotTheta = std::abs(deltaCotTheta);
-            // if deltaTheta larger than the scattering for the lower pT cut,
-            // skip
-            error = std::sqrt(error2);
-            dCotThetaMinusError2 =
-                deltaCotTheta2 + error2 - 2 * deltaCotTheta * error;
-            // avoid taking root of scatteringInRegion
-            // if left side of ">" is positive, both sides of unequality can be
-            // squared
-            // (scattering is always positive)
-            if (dCotThetaMinusError2 > scatteringInRegion2) {
-              continue;
-            }
           }
+          continue;
         }
 
         // protects against division by 0
@@ -347,7 +334,7 @@ void Seedfinder<external_spacepoint_t, platform_t>::createSeedsForGroup(
 
         // refinement of the cut on the compatibility between the r-z slope of
         // the two seed segments using a scattering term scaled by the actual
-        // measured pT
+        // measured pT (p2scatterSigma)
         float iHelixDiameter2 = B2 / S2;
         // calculate scattering for p(T) calculated from seed curvature
         float pT2scatterSigma = iHelixDiameter2 * m_config.sigmapT2perRadius;
@@ -364,21 +351,15 @@ void Seedfinder<external_spacepoint_t, platform_t>::createSeedsForGroup(
         // from rad to deltaCotTheta
         float p2scatterSigma = pT2scatterSigma * iSinTheta2;
         // if deltaTheta larger than allowed scattering for calculated pT, skip
-        if (m_config.enableCutsForSortedSP) {
-          if (deltaCotTheta2 - error2 > p2scatterSigma) {
+        if (deltaCotTheta2 > (error2 + p2scatterSigma)) {
+          if (m_config.skipPreviousTopSP) {
             if (cotThetaB - lt.cotTheta < 0) {
               break;
             }
             t0 = t;
-            continue;
-          }
-        } else {
-          if ((deltaCotTheta2 - error2 > 0) &&
-              (dCotThetaMinusError2 > p2scatterSigma)) {
-            continue;
           }
+          continue;
         }
-
         // A and B allow calculation of impact params in U/V plane with linear
         // function
         // (in contrast to having to solve a quadratic function in x/y plane)

diff --git a/Core/include/Acts/Seeding/SeedfinderConfig.hpp b/Core/include/Acts/Seeding/SeedfinderConfig.hpp
@@ -75,9 +75,8 @@ struct SeedfinderConfig {
   // non equidistant binning in z
   std::vector<float> zBinEdges;
 
-  // sort the SP in transformCoordinates method and enables compatibility cuts
-  // based on the sorting of cotTheta
-  bool enableCutsForSortedSP = false;
+  // additional cut to skip top SPs when producing triplets
+  bool skipPreviousTopSP = false;
 
   // FIXME: this is not used yet
   //        float upperPtResolutionPerSeed = 20* Acts::GeV;

diff --git a/Examples/Python/src/TrackFinding.cpp b/Examples/Python/src/TrackFinding.cpp
@@ -108,7 +108,7 @@ void addTrackFinding(Context& ctx) {
     ACTS_PYTHON_MEMBER(impactMax);
     ACTS_PYTHON_MEMBER(deltaZMax);
     ACTS_PYTHON_MEMBER(zBinEdges);
-    ACTS_PYTHON_MEMBER(enableCutsForSortedSP);
+    ACTS_PYTHON_MEMBER(skipPreviousTopSP);
     ACTS_PYTHON_MEMBER(interactionPointCut);
     ACTS_PYTHON_MEMBER(zBinEdges);
     ACTS_PYTHON_MEMBER(zBinsCustomLooping);

diff --git a/Examples/Python/tests/root_file_hashes.txt b/Examples/Python/tests/root_file_hashes.txt
@@ -1,9 +1,9 @@
-test_ckf_tracks_example[full_seeding-generic]__trackstates_ckf.root: d9f9e9a6cf4e123d1ef48070acf5e04bfae85f8d36ed6aeb171c9557de26f5e4
-test_ckf_tracks_example[full_seeding-generic]__tracksummary_ckf.root: ae2b0adff43681516be6c05a55ca586601c36a90f7e76ef861959fd921b37767
-test_ckf_tracks_example[full_seeding-generic]__performance_seeding_trees.root: 48f4a7c4a3aea7188eea2417d7cde50d61d7e7adb39a4c093a24a87c7576c2f5
-test_ckf_tracks_example[full_seeding-odd]__trackstates_ckf.root: 50b9bd069f3360f208defe4694d59dfd69e773888479a84edc43be1509cfdd87
-test_ckf_tracks_example[full_seeding-odd]__tracksummary_ckf.root: cb91fcc0a9ecc78aa361796c2393a13a0a9b47f22843f07e7f93e55d590473f8
-test_ckf_tracks_example[full_seeding-odd]__performance_seeding_trees.root: fb9769a523e3082b672358670715e2469e9544cd8554b8605981b003c0495231
+test_ckf_tracks_example[full_seeding-generic]__trackstates_ckf.root: 67e0899d0891a0b330555c07dbd31952698d3888262b3fb827331b3a43b8e8ea
+test_ckf_tracks_example[full_seeding-generic]__tracksummary_ckf.root: 9b1a1b432d6c82a60f08fdfd84f90ee566bdd2a86bee1c755cb064e71222712e
+test_ckf_tracks_example[full_seeding-generic]__performance_seeding_trees.root: 4890c1f1fb8618f59d4e5c2c0338af0595f2b92847b8eaced8a24e0fc1eb2934
+test_ckf_tracks_example[full_seeding-odd]__trackstates_ckf.root: 18f1328a34920b2485bc129aa7ac0ab2bd13432fbd6e379295e3c4ee3a62006f
+test_ckf_tracks_example[full_seeding-odd]__tracksummary_ckf.root: f323c0d6bfa66429e74f103269489e9b70a85f92e5cf319d6b5a91834a339226
+test_ckf_tracks_example[full_seeding-odd]__performance_seeding_trees.root: 74a1754b03812343dfbec54276152294ce92da3e35d6ad43a393580101b9a678
 test_ckf_tracks_example[truth_estimated-generic]__trackstates_ckf.root: 9e2afb3defe3f4f32d7fe5b4ca342654e4ca97d585569985998e229dfa8f0da6
 test_ckf_tracks_example[truth_estimated-generic]__tracksummary_ckf.root: 94267dea26347dea32a235c5bd034f2fc52bf5ab40be42be311aa7ef89bd7f09
 test_ckf_tracks_example[truth_estimated-generic]__performance_seeding_trees.root: 471373aceff30b77fc0844eef16fd423911f46981ebe7e78a6cbb239c315d04f

diff --git a/Examples/Scripts/Python/itk_seeding.py b/Examples/Scripts/Python/itk_seeding.py
@@ -72,7 +72,7 @@ def runITkSeeding(field, csvInputDir, outputDir, s=None):
         deltaZMax=900 * u.mm,
         interactionPointCut=True,
         zBinEdges=gridConfig.zBinEdges,
-        enableCutsForSortedSP=True,  # enable cotTheta sorting in SeedFinder
+        skipPreviousTopSP=True,
         rRangeMiddleSP=[
             [40.0, 90.0],
             [40.0, 200.0],