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