|
|
|||
File indexing completed on 2025-11-07 09:15:08
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/Definitions/Algebra.hpp" 0012 #include "Acts/Material/IVolumeMaterial.hpp" 0013 #include "Acts/Material/Material.hpp" 0014 #include "Acts/Utilities/BinUtility.hpp" 0015 #include "Acts/Utilities/Interpolation.hpp" 0016 0017 #include <functional> 0018 #include <optional> 0019 0020 namespace Acts { 0021 0022 /// @brief Struct for mapping global 3D positions to material values 0023 /// 0024 /// Global 3D positions are transformed into a @c DIM_POS Dimensional vector 0025 /// which is used to look up the material classification value in the 0026 /// underlying material map. 0027 template <typename G> 0028 struct MaterialMapper { 0029 public: 0030 /// Type alias for material grid 0031 using Grid_t = G; 0032 /// Dimensionality of the position space for material interpolation 0033 static constexpr std::size_t DIM_POS = Grid_t::DIM; 0034 0035 /// @brief Struct representing smallest grid unit in material grid 0036 /// 0037 /// This type encapsulate all required information to perform linear 0038 /// interpolation of material classification values within a 3D volume. 0039 struct MaterialCell { 0040 /// Number of corner points defining the confining hyper-box 0041 static constexpr unsigned int N = 1 << DIM_POS; 0042 0043 /// @brief Default constructor 0044 /// 0045 /// @param [in] transformPos Mapping of global 3D coordinates onto grid 0046 /// space 0047 /// @param [in] lowerLeft Generalized lower-left corner of hyper box 0048 /// (containing the minima of the hyper box along 0049 /// each Dimension) 0050 /// @param [in] upperRight Generalized upper-right corner of hyper box 0051 /// (containing the maxima of the hyper box along 0052 /// each Dimension) 0053 /// @param [in] materialValues Material classification values at the hyper 0054 /// box corners sorted in the canonical order defined in Acts::interpolate 0055 MaterialCell( 0056 std::function<ActsVector<DIM_POS>(const Vector3&)> transformPos, 0057 std::array<double, DIM_POS> lowerLeft, 0058 std::array<double, DIM_POS> upperRight, 0059 std::array<Material::ParametersVector, N> materialValues) 0060 : m_transformPos(std::move(transformPos)), 0061 m_lowerLeft(std::move(lowerLeft)), 0062 m_upperRight(std::move(upperRight)), 0063 m_materialValues(std::move(materialValues)) {} 0064 0065 /// @brief Retrieve material at given position 0066 /// 0067 /// @param [in] position Global 3D position 0068 /// @return Material at the given position 0069 /// 0070 /// @pre The given @c position must lie within the current cell. 0071 Material getMaterial(const Vector3& position) const { 0072 // defined in Interpolation.hpp 0073 return Material(interpolate(m_transformPos(position), m_lowerLeft, 0074 m_upperRight, m_materialValues)); 0075 } 0076 0077 /// @brief Check whether given 3D position is inside this cell 0078 /// 0079 /// @param [in] position Global 3D position 0080 /// @return @c true if position is inside the current cell, 0081 /// otherwise @c false 0082 bool isInside(const Vector3& position) const { 0083 const auto& gridCoordinates = m_transformPos(position); 0084 for (unsigned int i = 0; i < DIM_POS; ++i) { 0085 if (gridCoordinates[i] < m_lowerLeft.at(i) || 0086 gridCoordinates[i] >= m_upperRight.at(i)) { 0087 return false; 0088 } 0089 } 0090 return true; 0091 } 0092 0093 private: 0094 /// Geometric transformation applied to global 3D positions 0095 std::function<ActsVector<DIM_POS>(const Vector3&)> m_transformPos; 0096 0097 /// Generalized lower-left corner of the confining hyper-box 0098 std::array<double, DIM_POS> m_lowerLeft; 0099 0100 /// Generalized upper-right corner of the confining hyper-box 0101 std::array<double, DIM_POS> m_upperRight; 0102 0103 /// @brief Material component vectors at the hyper-box corners 0104 /// 0105 /// @note These values must be order according to the prescription detailed 0106 /// in Acts::interpolate. 0107 std::array<Material::ParametersVector, N> m_materialValues; 0108 }; 0109 0110 /// @brief Default constructor 0111 /// 0112 /// @param [in] transformPos Mapping of global 3D coordinates (cartesian) 0113 /// onto grid space 0114 /// @param [in] grid Grid storing material classification values 0115 MaterialMapper( 0116 std::function<ActsVector<DIM_POS>(const Vector3&)> transformPos, 0117 Grid_t grid) 0118 : m_transformPos(std::move(transformPos)), m_grid(std::move(grid)) {} 0119 0120 /// @brief Retrieve binned material at given position 0121 /// 0122 /// @param [in] position Global 3D position 0123 /// @return Material at the given position 0124 /// 0125 /// @pre The given @c position must lie within the range of the underlying 0126 /// map. 0127 Material material(const Vector3& position) const { 0128 return Material(m_grid.atLocalBins( 0129 m_grid.localBinsFromLowerLeftEdge(m_transformPos(position)))); 0130 } 0131 0132 /// @brief Retrieve interpolated material at given position 0133 /// 0134 /// @param [in] position Global 3D position 0135 /// @return Material at the given position 0136 /// 0137 /// @pre The given @c position must lie within the range of the underlying 0138 /// map. 0139 Material getMaterial(const Vector3& position) const { 0140 return Material(m_grid.interpolate(m_transformPos(position))); 0141 } 0142 0143 /// @brief Retrieve material cell for given position 0144 /// 0145 /// @param [in] position Global 3D position 0146 /// @return material cell containing the given global position 0147 /// 0148 /// @pre The given @c position must lie within the range of the underlying 0149 /// map. 0150 MaterialCell getMaterialCell(const Vector3& position) const { 0151 const auto& gridPosition = m_transformPos(position); 0152 std::size_t bin = m_grid.globalBinFromPosition(gridPosition); 0153 const auto& indices = m_grid.localBinsFromPosition(bin); 0154 const auto& lowerLeft = m_grid.lowerLeftBinEdge(indices); 0155 const auto& upperRight = m_grid.upperRightBinEdge(indices); 0156 0157 // Loop through all corner points 0158 constexpr std::size_t nCorners = 1 << DIM_POS; 0159 std::array<Material::ParametersVector, nCorners> neighbors; 0160 const auto& cornerIndices = m_grid.closestPointsIndices(gridPosition); 0161 0162 std::size_t i = 0; 0163 for (std::size_t index : cornerIndices) { 0164 neighbors.at(i++) = m_grid.at(index); 0165 } 0166 0167 return MaterialCell(m_transformPos, lowerLeft, upperRight, 0168 std::move(neighbors)); 0169 } 0170 0171 /// @brief Get the number of bins for all axes of the map 0172 /// 0173 /// @return Vector returning number of bins for all map axes 0174 std::vector<std::size_t> getNBins() const { 0175 auto nBinsArray = m_grid.numLocalBins(); 0176 return std::vector<std::size_t>(nBinsArray.begin(), nBinsArray.end()); 0177 } 0178 0179 /// @brief Get the minimum value of all axes of the map 0180 /// 0181 /// @return Vector returning the minima of all map axes 0182 std::vector<double> getMin() const { 0183 auto minArray = m_grid.minPosition(); 0184 return std::vector<double>(minArray.begin(), minArray.end()); 0185 } 0186 0187 /// @brief Get the maximum value of all axes of the map 0188 /// 0189 /// @return Vector returning the maxima of all map axes 0190 std::vector<double> getMax() const { 0191 auto maxArray = m_grid.maxPosition(); 0192 return std::vector<double>(maxArray.begin(), maxArray.end()); 0193 } 0194 0195 /// @brief Check whether given 3D position is inside look-up domain 0196 /// 0197 /// @param [in] position Global 3D position 0198 /// @return @c true if position is inside the defined look-up grid, 0199 /// otherwise @c false 0200 bool isInside(const Vector3& position) const { 0201 return m_grid.isInside(m_transformPos(position)); 0202 } 0203 0204 /// @brief Get a const reference on the underlying grid structure 0205 /// 0206 /// @return Grid reference 0207 const Grid_t& getGrid() const { return m_grid; } 0208 0209 private: 0210 /// Geometric transformation applied to global 3D positions 0211 std::function<ActsVector<DIM_POS>(const Vector3&)> m_transformPos; 0212 /// Grid storing material values 0213 Grid_t m_grid; 0214 }; 0215 0216 /// @ingroup Material 0217 /// @brief Interpolate material classification values from material values on a 0218 /// given grid 0219 /// 0220 /// This class implements a material service which is initialized by a 0221 /// material map defined by: 0222 /// - a list of material values on a regular grid in some n-Dimensional space, 0223 /// - a transformation of global 3D coordinates onto this n-Dimensional 0224 /// space. 0225 /// - a transformation of local n-Dimensional material coordinates into 0226 /// global (cartesian) 3D coordinates 0227 /// 0228 /// The material value for a given global position is then determined by: 0229 /// - mapping the position onto the grid, 0230 /// - looking up the material classification values on the closest grid points, 0231 /// - doing a linear interpolation of these values. 0232 /// @warning Each classification number of the material is interpolated 0233 /// independently and thus does not consider any correlations that exists 0234 /// between these values. This might work out since the used material is already 0235 /// a mean of the materials in a certain bin and can therewith be treated as a 0236 /// collection of numbers. 0237 /// @tparam G Type of the grid 0238 template <typename Mapper_t> 0239 class InterpolatedMaterialMap : public IVolumeMaterial { 0240 public: 0241 /// @brief Temporary storage of a certain cell to improve material access 0242 struct Cache { 0243 /// Stored material cell 0244 std::optional<typename Mapper_t::MaterialCell> matCell; 0245 /// Boolean statement if the cell is initialized 0246 bool initialized = false; 0247 }; 0248 0249 /// @brief Create interpolated map 0250 /// 0251 /// @param [in] mapper Material map 0252 explicit InterpolatedMaterialMap(Mapper_t&& mapper) 0253 : m_mapper(std::move(mapper)) {} 0254 0255 /// @brief Create interpolated map 0256 /// 0257 /// @param [in] mapper Material map 0258 /// @param [in] bu @c BinUtility for build from 0259 InterpolatedMaterialMap(Mapper_t&& mapper, BinUtility bu) 0260 : m_mapper(std::move(mapper)), m_binUtility(std::move(bu)) {} 0261 0262 /// @brief Retrieve the binned material 0263 /// 0264 /// @param [in] position Global 3D position 0265 /// 0266 /// @return Material at given position 0267 const Material material(const Vector3& position) const override { 0268 return m_mapper.material(position); 0269 } 0270 0271 /// @brief Retrieve the interpolated material 0272 /// 0273 /// @param [in] position Global 3D position 0274 /// 0275 /// @return material at given position 0276 Material getMaterial(const Vector3& position) const { 0277 return m_mapper.getMaterial(position); 0278 } 0279 0280 /// @brief Retrieve material 0281 /// 0282 /// @param [in] position Global 3D position 0283 /// @param [in,out] cache Cache object. Contains material cell used for 0284 /// interpolation 0285 /// 0286 /// @return material at given position 0287 Material getMaterial(const Vector3& position, Cache& cache) const { 0288 if (!cache.initialized || !(*cache.matCell).isInside(position)) { 0289 cache.matCell = getMaterialCell(position); 0290 cache.initialized = true; 0291 } 0292 return (*cache.matCell).getMaterial(position); 0293 } 0294 0295 /// @brief Retrieve material value & its "gradient" 0296 /// 0297 /// @param [in] position Global 3D position 0298 /// @return Material 0299 /// 0300 /// @note Currently the derivative is not calculated 0301 /// @todo return derivative 0302 Material getMaterialGradient(const Vector3& position, 0303 ActsMatrix<5, 5>& /*derivative*/) const { 0304 return m_mapper.getMaterial(position); 0305 } 0306 0307 /// @brief Retrieve material value & its "gradient" 0308 /// 0309 /// @param [in] position Global 3D position 0310 /// @return Material 0311 /// 0312 /// @note Currently the derivative is not calculated 0313 /// @note Cache is not used currently 0314 /// @todo return derivative 0315 Material getMaterialGradient(const Vector3& position, 0316 ActsMatrix<5, 5>& /*derivative*/, 0317 Cache& /*cache*/) const { 0318 return m_mapper.getMaterial(position); 0319 } 0320 0321 /// @brief Convenience method to access underlying material mapper 0322 /// 0323 /// @return The material mapper 0324 const Mapper_t& getMapper() const { return m_mapper; } 0325 0326 /// @brief Check whether given 3D position is inside look-up domain 0327 /// 0328 /// @param [in] position Global 3D position 0329 /// @return @c true if position is inside the defined map, otherwise @c false 0330 bool isInside(const Vector3& position) const { 0331 return m_mapper.isInside(position); 0332 } 0333 0334 /// Return the BinUtility 0335 /// @return Const reference to the bin utility for the material map 0336 const BinUtility& binUtility() const { return m_binUtility; } 0337 0338 /// Output Method for std::ostream 0339 /// 0340 /// @param sl The outoput stream 0341 /// @return Reference to the output stream for method chaining 0342 std::ostream& toStream(std::ostream& sl) const override { 0343 sl << "Acts::InterpolatedMaterialMap : " << std::endl; 0344 sl << " - Number of Material bins [0,1] : " << m_binUtility.max(0) + 1 0345 << " / " << m_binUtility.max(1) + 1 << std::endl; 0346 sl << " - Parse full update material : " << std::endl; 0347 return sl; 0348 } 0349 0350 private: 0351 /// @brief Retrieve cell for given position 0352 /// 0353 /// @param [in] position Global 3D position 0354 /// @return Material cell containing the given global position 0355 /// 0356 /// @pre The given @c position must lie within the range of the underlying 0357 /// map. 0358 typename Mapper_t::MaterialCell getMaterialCell( 0359 const Vector3& position) const { 0360 return m_mapper.getMaterialCell(position); 0361 } 0362 0363 /// @brief object for global coordinate transformation and interpolation 0364 /// 0365 /// This object performs the mapping of the global 3D coordinates onto the 0366 /// material grid and the interpolation of the material component values on 0367 /// close-by grid points. 0368 Mapper_t m_mapper; 0369 0370 BinUtility m_binUtility{}; 0371 }; 0372 } // namespace Acts
| [ Source navigation ] | [ Diff markup ] | [ Identifier search ] | [ general search ] |
|
This page was automatically generated by the 2.3.7 LXR engine. The LXR team |
|