Back to home page

EIC code displayed by LXR

 
 

    


File indexing completed on 2025-06-30 07:51:44

0001 // This file is part of the ACTS project.
0002 //
0003 // Copyright (C) 2016 CERN for the benefit of the ACTS project
0004 //
0005 // This Source Code Form is subject to the terms of the Mozilla Public
0006 // License, v. 2.0. If a copy of the MPL was not distributed with this
0007 // file, You can obtain one at https://mozilla.org/MPL/2.0/.
0008 
0009 #pragma once
0010 
0011 #include "Acts/Geometry/GeometryIdentifier.hpp"
0012 
0013 #include <algorithm>
0014 #include <cassert>
0015 #include <initializer_list>
0016 #include <iterator>
0017 #include <stdexcept>
0018 #include <utility>
0019 #include <vector>
0020 
0021 namespace Acts {
0022 
0023 /// Store values mapped into the geometry hierarchy.
0024 ///
0025 /// @tparam value_t stored value type
0026 ///
0027 /// The core functionality is to find an equivalent element, i.e. an
0028 /// identifier-value pair, for a given geometry identifier via
0029 ///
0030 ///     auto it = container.find(GeometryIdentifier(...));
0031 ///     if (it != container.end()) {
0032 ///         ...
0033 ///     }
0034 ///
0035 /// Trailing zero levels of stored geometry identifiers are used as
0036 /// broadcast values to refer to higher-level objects within the geometry, e.g.
0037 /// a geometry identifier with vanishing approach and sensitive index identifies
0038 /// a layer. An entry will all geometry identifier levels set to zero acts as
0039 /// the global default value.
0040 ///
0041 /// The container also supports range-based iteration over all stored elements
0042 ///
0043 ///     for (const auto& element : container) {
0044 ///         ...
0045 ///     }
0046 ///
0047 /// and index-based access to stored elements and associated geometry
0048 /// identifiers
0049 ///
0050 ///     GeometryIdentifier id3 = container.idAt(3);
0051 ///     const auto& element4 = container.valueAt(4);
0052 ///
0053 /// @note No guarantees are given for the element order when using range-based
0054 ///   or index-based access. Any apparent ordering must be considered an
0055 ///   implementation detail and might change.
0056 ///
0057 /// Adding elements is potentially expensive as the internal lookup structure
0058 /// must be updated. In addition, modifying an element in-place could change its
0059 /// identifier which would also break the lookup. Thus, the container can not be
0060 /// modified after construction to prevent misuse.
0061 template <typename value_t>
0062 class GeometryHierarchyMap {
0063  public:
0064   /// Combined geometry identifier and value element. Only used for input.
0065   using InputElement = typename std::pair<GeometryIdentifier, value_t>;
0066   using Iterator = typename std::vector<value_t>::const_iterator;
0067   using Value = value_t;
0068 
0069   /// Construct the container from the given elements.
0070   ///
0071   /// @param elements input elements (must be unique with respect to identifier)
0072   explicit GeometryHierarchyMap(std::vector<InputElement> elements);
0073 
0074   /// Construct the container from an initializer list.
0075   ///
0076   /// @param elements input initializer list
0077   GeometryHierarchyMap(std::initializer_list<InputElement> elements);
0078 
0079   // defaulted constructors and assignment operators
0080   GeometryHierarchyMap() = default;
0081   GeometryHierarchyMap(const GeometryHierarchyMap&) = default;
0082   GeometryHierarchyMap(GeometryHierarchyMap&&) noexcept = default;
0083   ~GeometryHierarchyMap() = default;
0084   GeometryHierarchyMap& operator=(const GeometryHierarchyMap&) = default;
0085   GeometryHierarchyMap& operator=(GeometryHierarchyMap&&) noexcept = default;
0086 
0087   /// Return an iterator pointing to the beginning of the stored values.
0088   Iterator begin() const { return m_values.begin(); }
0089 
0090   /// Return an iterator pointing to the end of the stored values.
0091   Iterator end() const { return m_values.end(); }
0092 
0093   /// Check if any elements are stored.
0094   bool empty() const { return m_values.empty(); }
0095 
0096   /// Return the number of stored elements.
0097   std::size_t size() const { return m_values.size(); }
0098 
0099   /// Access the geometry identifier for the i-th element with bounds check.
0100   ///
0101   /// @throws std::out_of_range for invalid indices
0102   GeometryIdentifier idAt(std::size_t index) const {
0103     return GeometryIdentifier(m_ids.at(index));
0104   }
0105 
0106   /// Access the value of the i-th element in the container with bounds check.
0107   ///
0108   /// @throws std::out_of_range for invalid indices
0109   const Value& valueAt(std::size_t index) const { return m_values.at(index); }
0110 
0111   /// Find the most specific value for a given geometry identifier.
0112   ///
0113   /// This can be either from the element matching exactly to the given geometry
0114   /// id, if it exists, or from the element for the next available higher level
0115   /// within the geometry hierarchy.
0116   ///
0117   /// @param id geometry identifier for which information is requested
0118   /// @retval iterator to an existing value
0119   /// @retval `.end()` iterator if no matching element exists
0120   Iterator find(const GeometryIdentifier& id) const;
0121 
0122   /// Check if the most specific value exists for a given geometry identifier.
0123   ///
0124   /// This function checks if there is an element matching exactly the given
0125   /// geometry id, or from the element for the next available higher level
0126   /// within the geometry hierarchy.
0127   ///
0128   /// @param id geometry identifier for which existence is being checked
0129   /// @retval `true` if a matching element exists
0130   /// @retval `false` if no matching element exists
0131   bool contains(const GeometryIdentifier& id) const;
0132 
0133  private:
0134   // NOTE this class assumes that it knows the ordering of the levels within
0135   //      the geometry id. if the geometry id changes, this code has to be
0136   //      adapted too. the asserts ensure that such a change is caught.
0137   static_assert(GeometryIdentifier().withVolume(1).value() <
0138                     GeometryIdentifier().withVolume(1).withBoundary(1).value(),
0139                 "Incompatible GeometryIdentifier hierarchy");
0140   static_assert(GeometryIdentifier().withBoundary(1).value() <
0141                     GeometryIdentifier().withBoundary(1).withLayer(1).value(),
0142                 "Incompatible GeometryIdentifier hierarchy");
0143   static_assert(GeometryIdentifier().withLayer(1).value() <
0144                     GeometryIdentifier().withLayer(1).withApproach(1).value(),
0145                 "Incompatible GeometryIdentifier hierarchy");
0146   static_assert(
0147       GeometryIdentifier().withApproach(1).value() <
0148           GeometryIdentifier().withApproach(1).withSensitive(1).value(),
0149       "Incompatible GeometryIdentifier hierarchy");
0150 
0151   using Identifier = GeometryIdentifier::Value;
0152 
0153   // encoded ids for all elements for faster lookup.
0154   std::vector<Identifier> m_ids;
0155   // validity bit masks for the ids: which parts to use for comparison
0156   std::vector<Identifier> m_masks;
0157   std::vector<Value> m_values;
0158 
0159   /// Construct a mask where all leading non-zero levels are set.
0160   static constexpr Identifier makeLeadingLevelsMask(GeometryIdentifier id) {
0161     // construct id from encoded value with all bits set
0162     const auto allSet = GeometryIdentifier(~GeometryIdentifier::Value{0u});
0163     // manually iterate over identifier levels starting from the lowest
0164     if (id.sensitive() != 0u) {
0165       // all levels are valid; keep all bits set.
0166       return allSet.withExtra(0u).value();
0167     }
0168     if (id.approach() != 0u) {
0169       return allSet.withExtra(0u).withSensitive(0u).value();
0170     }
0171     if (id.layer() != 0u) {
0172       return allSet.withExtra(0u).withSensitive(0u).withApproach(0u).value();
0173     }
0174     if (id.boundary() != 0u) {
0175       return allSet.withExtra(0u)
0176           .withSensitive(0u)
0177           .withApproach(0u)
0178           .withLayer(0u)
0179           .value();
0180     }
0181     if (id.volume() != 0u) {
0182       return allSet.withExtra(0u)
0183           .withSensitive(0u)
0184           .withApproach(0u)
0185           .withLayer(0u)
0186           .withBoundary(0u)
0187           .value();
0188     }
0189     // no valid levels; all bits are zero.
0190     return Identifier{0u};
0191   }
0192 
0193   /// Construct a mask where only the highest level is set.
0194   static constexpr Identifier makeHighestLevelMask() {
0195     return makeLeadingLevelsMask(GeometryIdentifier(0u).withVolume(1u));
0196   }
0197 
0198   /// Compare the two identifiers only within the masked bits.
0199   static constexpr bool equalWithinMask(Identifier lhs, Identifier rhs,
0200                                         Identifier mask) {
0201     return (lhs & mask) == (rhs & mask);
0202   }
0203 
0204   /// Ensure identifier ordering and uniqueness.
0205   static void sortAndCheckDuplicates(std::vector<InputElement>& elements);
0206 
0207   /// Fill the container from the input elements.
0208   ///
0209   /// This assumes that the elements are ordered and unique with respect to
0210   /// their identifiers.
0211   void fill(const std::vector<InputElement>& elements);
0212 };
0213 
0214 // implementations
0215 
0216 template <typename value_t>
0217 inline GeometryHierarchyMap<value_t>::GeometryHierarchyMap(
0218     std::vector<InputElement> elements) {
0219   sortAndCheckDuplicates(elements);
0220   fill(elements);
0221 }
0222 
0223 template <typename value_t>
0224 inline GeometryHierarchyMap<value_t>::GeometryHierarchyMap(
0225     std::initializer_list<InputElement> elements)
0226     : GeometryHierarchyMap(
0227           std::vector<InputElement>(elements.begin(), elements.end())) {}
0228 
0229 template <typename value_t>
0230 inline void GeometryHierarchyMap<value_t>::sortAndCheckDuplicates(
0231     std::vector<InputElement>& elements) {
0232   // ensure elements are sorted by identifier
0233   std::ranges::sort(elements, [=](const auto& lhs, const auto& rhs) {
0234     return lhs.first < rhs.first;
0235   });
0236 
0237   // Check that all elements have unique identifier
0238   auto dup = std::ranges::adjacent_find(
0239       elements,
0240       [](const auto& lhs, const auto& rhs) { return lhs.first == rhs.first; });
0241 
0242   if (dup != elements.end()) {
0243     throw std::invalid_argument("Input elements contain duplicates");
0244   }
0245 }
0246 
0247 template <typename value_t>
0248 inline void GeometryHierarchyMap<value_t>::fill(
0249     const std::vector<InputElement>& elements) {
0250   m_ids.clear();
0251   m_masks.clear();
0252   m_values.clear();
0253 
0254   m_ids.reserve(elements.size());
0255   m_masks.reserve(elements.size());
0256   m_values.reserve(elements.size());
0257 
0258   for (const auto& element : elements) {
0259     m_ids.push_back(element.first.value());
0260     m_masks.push_back(
0261         makeLeadingLevelsMask(GeometryIdentifier(element.first.value())));
0262     m_values.push_back(std::move(element.second));
0263   }
0264 }
0265 
0266 template <typename value_t>
0267 inline auto GeometryHierarchyMap<value_t>::find(
0268     const GeometryIdentifier& id) const -> Iterator {
0269   assert((m_ids.size() == m_values.size()) &&
0270          "Inconsistent container state: #ids != # values");
0271   assert((m_masks.size() == m_values.size()) &&
0272          "Inconsistent container state: #masks != #values");
0273 
0274   // we can not search for the element directly since the relevant one
0275   // might be stored at a higher level. ids for higher levels would always
0276   // be sorted before the requested id. searching for the first element
0277   // after the requested ensures that we include the full hierarchy.
0278   const auto it = std::upper_bound(m_ids.begin(), m_ids.end(), id.value());
0279   auto i = std::distance(m_ids.begin(), it);
0280 
0281   // now go up the hierarchy to find the first matching element.
0282   // example: the container stores four identifiers
0283   //
0284   //     2|x|x (volume-only)
0285   //     2|2|1 (volume, layer, and sensitive)
0286   //     2|3|x (volume and layer)
0287   //     2|3|4 (volume, layer, and sensitive)
0288   //
0289   // where | marks level boundaries. searching for either 2|3|4, 2|3|7, or
0290   // 2|4|x would first point to 2|3|4 and thus needs to go up the hierarchy.
0291   while (0 < i) {
0292     // index always starts after item of interest due to upper bound search
0293     --i;
0294 
0295     // if the input id does not even match at the highest hierarchy level
0296     // with the current comparison id, then have reached the end of this
0297     // hierarchy. having a special check for the highest level avoids an
0298     // unbounded search window all the way to the beginning of the container for
0299     // the global default entry.
0300     if (!equalWithinMask(id.value(), m_ids[i], makeHighestLevelMask())) {
0301       // check if a global default entry exists
0302       if (m_ids.front() == Identifier{0u}) {
0303         return begin();
0304       } else {
0305         return end();
0306       }
0307     }
0308 
0309     // since the search is going backwards in the sorted container, it
0310     // progresses from more specific to less specific elements. the first
0311     // match is automatically the appropriate one.
0312     if (equalWithinMask(id.value(), m_ids[i], m_masks[i])) {
0313       return std::next(begin(), i);
0314     }
0315   }
0316 
0317   // all options are exhausted and no matching element was found.
0318   return end();
0319 }
0320 
0321 template <typename value_t>
0322 inline auto GeometryHierarchyMap<value_t>::contains(
0323     const GeometryIdentifier& id) const -> bool {
0324   return this->find(id) != this->end();
0325 }
0326 
0327 }  // namespace Acts