src/seeding/triplet_finding.hpp

0001 /** TRACCC library, part of the ACTS project (R&D line)
0002  *
0003  * (c) 2021-2025 CERN for the benefit of the ACTS project
0004  *
0005  * Mozilla Public License Version 2.0
0006  */
0007
0008 #pragma once
0009
0010 // Local include(s).
0011 #include "traccc/edm/spacepoint_collection.hpp"
0012 #include "traccc/seeding/detail/doublet.hpp"
0013 #include "traccc/seeding/detail/lin_circle.hpp"
0014 #include "traccc/seeding/detail/seeding_config.hpp"
0015 #include "traccc/seeding/detail/triplet.hpp"
0016 #include "traccc/seeding/triplet_finding_helper.hpp"
0017 #include "traccc/utils/messaging.hpp"
0018
0019 // VecMem include(s).
0020 #include <vecmem/memory/memory_resource.hpp>
0021
0022 // System include(s).
0023 #include <cassert>
0024 #include <functional>
0025
0026 namespace traccc::host::details {
0027
0028 /// Triplet finding to search the compatible combintations of two doublets which
0029 /// share same middle spacepoint
0030 struct triplet_finding : public messaging {
0031   /// Constructor for the triplet finding
0032   ///
0033   /// @param finder_config Seed finding configuration parameters
0034   /// @param filter_config Seed filtering configuration parameters
0035   /// @param mr The memory resource to use
0036   ///
0037   triplet_finding(
0038       const seedfinder_config& finding_config,
0039       const seedfilter_config& filter_config, vecmem::memory_resource& mr,
0040       std::unique_ptr<const Logger> logger = getDummyLogger().clone())
0041       : messaging(std::move(logger)),
0042         m_finding_config{finding_config},
0043         m_filter_config{filter_config},
0044         m_mr{mr} {}
0045
0046   /// Callable operator for triplet finding per middle-bottom doublet
0047   ///
0048   /// @param spacepoints All the spacepoints in the event
0049   /// @param sp_grid The spacepoint grid to use
0050   /// @param mid_bot_doublet is the current middle-bottom doublets
0051   /// @param mid_bot_lc is transformed coordinate of @c mid_bot_doublet
0052   /// @param mid_top_doublets is the vector of middle-top doublets which share
0053   ///                         same middle spacepoint with current
0054   ///                         middle-bottom doublet
0055   /// @param mid_top_lcs is transformed coordinates of
0056   ///                    @c doublets_mid_top
0057   ///
0058   /// @return a vector of triplets
0059   triplet_collection_types::host operator()(
0060       const edm::spacepoint_collection::const_device& spacepoints,
0061       const traccc::details::spacepoint_grid_types::host& sp_grid,
0062       const doublet& mid_bot_doublet, const lin_circle& mid_bot_lc,
0063       const doublet_collection_types::host& mid_top_doublets,
0064       const lin_circle_collection_types::host& mid_top_lcs) const {
0065     // Create the output.
0066     triplet_collection_types::host result{&(m_mr.get())};
0067
0068     // Access the middle spacepoint that all the doublets share.
0069     const edm::spacepoint_collection::const_device::const_proxy_type spM =
0070         spacepoints.at(sp_grid.bin(
0071             mid_bot_doublet.sp1.bin_idx)[mid_bot_doublet.sp1.sp_idx]);
0072
0073     // Calculate quantities that help deciding if two doublets are
0074     // compatible.
0075     const scalar iSinTheta2 =
0076         1.f + mid_bot_lc.cotTheta() * mid_bot_lc.cotTheta();
0077     scalar scatteringInRegion2 =
0078         m_finding_config.maxScatteringAngle2 * iSinTheta2;
0079     scatteringInRegion2 *=
0080         m_finding_config.sigmaScattering * m_finding_config.sigmaScattering;
0081     scalar curvature, impact_parameter;
0082
0083     // Compare this mid-bottom doublet with all the mid-top doublets.
0084     assert(mid_top_doublets.size() == mid_top_lcs.size());
0085     for (std::size_t i = 0; i < mid_top_doublets.size(); ++i) {
0086       const doublet& mid_top_doublet = mid_top_doublets.at(i);
0087       const lin_circle& mid_top_lc = mid_top_lcs.at(i);
0088
0089       if (!triplet_finding_helper::isCompatible(
0090               spM, mid_bot_lc, mid_top_lc, m_finding_config, iSinTheta2,
0091               scatteringInRegion2, curvature, impact_parameter)) {
0092         continue;
0093       }
0094
0095       result.push_back({mid_bot_doublet.sp2,  // bottom
0096                         mid_bot_doublet.sp1,  // middle
0097                         mid_top_doublet.sp2,  // top
0098                         curvature,            // curvature
0099                         -impact_parameter * m_filter_config.impactWeightFactor,
0100                         mid_bot_lc.Zo()});
0101     }
0102
0103     // Set the triplet weights in a super complicated double loop over the
0104     // triplets found in the previous step.
0105     for (std::size_t i = 0; i < result.size(); ++i) {
0106       triplet& current_triplet = result[i];
0107       const sp_location& current_spT_idx = current_triplet.sp3;
0108       const edm::spacepoint_collection::const_device::const_proxy_type
0109           current_spT = spacepoints.at(
0110               sp_grid.bin(current_spT_idx.bin_idx)[current_spT_idx.sp_idx]);
0111       const scalar currentTop_r = current_spT.radius();
0112
0113       // if two compatible seeds with high distance in r are found,
0114       // compatible seeds span 5 layers
0115       // -> very good seed
0116       std::vector<scalar> compatibleSeedR;
0117       scalar lowerLimitCurv =
0118           current_triplet.curvature - m_filter_config.deltaInvHelixDiameter;
0119       scalar upperLimitCurv =
0120           current_triplet.curvature + m_filter_config.deltaInvHelixDiameter;
0121
0122       for (std::size_t j = 0; j < result.size(); ++j) {
0123         if (i == j) {
0124           continue;
0125         }
0126
0127         const triplet& other_triplet = result[j];
0128         const sp_location& other_spT_idx = other_triplet.sp3;
0129         const edm::spacepoint_collection::const_device::const_proxy_type
0130             other_spT = spacepoints.at(
0131                 sp_grid.bin(other_spT_idx.bin_idx)[other_spT_idx.sp_idx]);
0132
0133         // compared top SP should have at least deltaRMin distance
0134         const scalar otherTop_r = other_spT.radius();
0135         const scalar deltaR = currentTop_r - otherTop_r;
0136         if (std::abs(deltaR) < m_filter_config.deltaRMin) {
0137           continue;
0138         }
0139
0140         // curvature difference within limits?
0141         // TODO: how much slower than sorting all vectors by curvature
0142         // and breaking out of loop? i.e. is vector size large (e.g. in
0143         // jets?)
0144         if (other_triplet.curvature < lowerLimitCurv) {
0145           continue;
0146         }
0147         if (other_triplet.curvature > upperLimitCurv) {
0148           continue;
0149         }
0150
0151         bool newCompSeed = true;
0152         for (scalar previousDiameter : compatibleSeedR) {
0153           // original ATLAS code uses higher min distance for 2nd
0154           // found compatible seed (20mm instead of 5mm) add new
0155           // compatible seed only if distance larger than rmin to all
0156           // other compatible seeds
0157           if (std::abs(previousDiameter - otherTop_r) <
0158               m_filter_config.deltaRMin) {
0159             newCompSeed = false;
0160             break;
0161           }
0162         }
0163
0164         if (newCompSeed) {
0165           compatibleSeedR.push_back(otherTop_r);
0166           current_triplet.weight += m_filter_config.compatSeedWeight;
0167         }
0168
0169         if (compatibleSeedR.size() >= m_filter_config.compatSeedLimit) {
0170           break;
0171         }
0172       }
0173     }
0174
0175     // Return the reconstructed triplets.
0176     return result;
0177   }
0178
0179  private:
0180   /// @name Triplet Finding Configuration
0181   /// @{
0182   seedfinder_config m_finding_config;
0183   seedfilter_config m_filter_config;
0184   /// @}
0185
0186   /// Memory resource to use
0187   std::reference_wrapper<vecmem::memory_resource> m_mr;
0188
0189 };  // struct triplet_finding
0190
0191 }  // namespace traccc::host::details