Acts/Material/MaterialComposition.hpp

0001 // This file is part of the Acts project.
0002 //
0003 // Copyright (C) 2018-2020 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 http://mozilla.org/MPL/2.0/.
0008
0009 #pragma once
0010
0011 #include <algorithm>
0012 #include <cassert>
0013 #include <cstdint>
0014 #include <vector>
0015
0016 namespace Acts {
0017
0018 /// Memory-efficient storage of the relative fraction of an element.
0019 ///
0020 /// This can be used to define materials that are compounds of multiple elements
0021 /// with varying fractions. The element is identified by its atomic number
0022 /// stored as a single byte (allows up to 256 elements; more than we need).
0023 /// Its fraction is also stored as a single byte with values between 0 and
0024 /// 255. This gives an accuracy of 1/256 ~ 0.5 %.
0025 ///
0026 /// The element fraction allows you to store element composition in merged
0027 /// materials with a large number of bins. Depending on the
0028 /// detector and the description granularity this can be a lot of information
0029 /// and thus requires the reduced memory footprint. This is really only needed
0030 /// for nuclear interaction in the fast simulation where the reduced fractional
0031 /// accuracy is not a problem. The fractional accuracy should be much better
0032 /// than the parametrization uncertainty for hadronic interactions.
0033 class ElementFraction {
0034  public:
0035   /// Construct from atomic number and relative fraction.
0036   ///
0037   /// @param e is the atomic number of the element
0038   /// @param f is the relative fraction and must be a value in [0,1]
0039   constexpr ElementFraction(unsigned int e, float f)
0040       : m_element(static_cast<uint8_t>(e)),
0041         m_fraction(static_cast<uint8_t>(f * UINT8_MAX)) {
0042     assert((0u < e) && ("The atomic number must be positive"));
0043     assert((0.0f <= f) && (f <= 1.0f) && "Relative fraction must be in [0,1]");
0044   }
0045   /// Construct from atomic number and integer weight.
0046   ///
0047   /// @param e is the atomic number of the element
0048   /// @param w is the integer weight and must be a value in [0,256)
0049   constexpr explicit ElementFraction(unsigned int e, unsigned int w)
0050       : m_element(static_cast<uint8_t>(e)),
0051         m_fraction(static_cast<uint8_t>(w)) {
0052     assert((0u < e) && ("The atomic number must be positive"));
0053     assert((w < 256u) && "Integer weight must be in [0,256)");
0054   }
0055
0056   /// Must always be created with valid data.
0057   ElementFraction() = delete;
0058   ElementFraction(ElementFraction&&) = default;
0059   ElementFraction(const ElementFraction&) = default;
0060   ~ElementFraction() = default;
0061   ElementFraction& operator=(ElementFraction&&) = default;
0062   ElementFraction& operator=(const ElementFraction&) = default;
0063
0064   /// The element atomic number.
0065   constexpr uint8_t element() const { return m_element; }
0066   /// The relative fraction of this element.
0067   constexpr float fraction() const {
0068     return static_cast<float>(m_fraction) / UINT8_MAX;
0069   }
0070
0071  private:
0072   // element atomic number
0073   uint8_t m_element;
0074   // element fraction in the compound scaled to the [0,256) range.
0075   uint8_t m_fraction;
0076
0077   friend constexpr bool operator==(ElementFraction lhs, ElementFraction rhs) {
0078     return (lhs.m_fraction == rhs.m_fraction) &&
0079            (lhs.m_element == rhs.m_element);
0080   }
0081   /// Sort by fraction for fastest access to the most probable element.
0082   friend constexpr bool operator<(ElementFraction lhs, ElementFraction rhs) {
0083     return lhs.m_fraction < rhs.m_fraction;
0084   }
0085   friend class MaterialComposition;
0086 };
0087
0088 /// Material composed from multiple elements with varying factions.
0089 ///
0090 /// @see ElementFraction for details.
0091 class MaterialComposition {
0092  public:
0093   /// Construct an empty composition corresponding to vacuum.
0094   MaterialComposition() = default;
0095   /// Constructor from element fractions.
0096   ///
0097   /// Rescales the fractions so they all add up to unity within the accuracy.
0098   MaterialComposition(std::vector<ElementFraction> elements)
0099       : m_elements(std::move(elements)) {
0100     std::sort(m_elements.begin(), m_elements.end());
0101     // compute the total weight first
0102     unsigned total = 0u;
0103     for (auto element : m_elements) {
0104       total += element.m_fraction;
0105     }
0106     // compute scale factor into the [0, 256) range
0107     float scale = float(UINT8_MAX) / float(total);
0108     for (auto& element : m_elements) {
0109       element.m_fraction = static_cast<uint8_t>(element.m_fraction * scale);
0110     }
0111   }
0112
0113   MaterialComposition(MaterialComposition&&) = default;
0114   MaterialComposition(const MaterialComposition&) = default;
0115   ~MaterialComposition() = default;
0116   MaterialComposition& operator=(MaterialComposition&&) = default;
0117   MaterialComposition& operator=(const MaterialComposition&) = default;
0118
0119   // Support range-based iteration over contained elements.
0120   auto begin() const { return m_elements.begin(); }
0121   auto end() const { return m_elements.end(); }
0122
0123   /// Check if the composed material is valid, i.e. it is not vacuum.
0124   operator bool() const { return !m_elements.empty(); }
0125   /// Return the number of elements.
0126   std::size_t size() const { return m_elements.size(); }
0127
0128  private:
0129   std::vector<ElementFraction> m_elements;
0130
0131   friend inline bool operator==(const MaterialComposition& lhs,
0132                                 const MaterialComposition& rhs) {
0133     return (lhs.m_elements == rhs.m_elements);
0134   }
0135 };
0136
0137 }  // namespace Acts