src/Material/Material.cpp

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 #include "Acts/Material/Material.hpp"
0010
0011 #include "Acts/Definitions/Units.hpp"
0012
0013 #include <cmath>
0014 #include <ostream>
0015
0016 namespace Acts {
0017
0018 namespace {
0019 enum MaterialClassificationNumberIndices {
0020   eRadiationLength = 0,
0021   eInteractionLength = 1,
0022   eRelativeAtomicMass = 2,
0023   eNuclearCharge = 3,
0024   eMolarDensity = 4,
0025 };
0026
0027 // Avogadro constant
0028 constexpr double kAvogadro = 6.02214076e23 / UnitConstants::mol;
0029
0030 constexpr float calculateMolarElectronDensity(float z, float molarRho) {
0031   return z * molarRho;
0032 }
0033
0034 constexpr float approximateMeanExcitationEnergy(float z) {
0035   using namespace UnitLiterals;
0036
0037   // use approximative computation as defined in ATL-SOFT-PUB-2008-003
0038   return 16_eV * std::pow(z, 0.9f);
0039 }
0040 }  // namespace
0041
0042 Material Material::fromMassDensity(float x0, float l0, float ar, float z,
0043                                    float massRho) {
0044   using namespace UnitLiterals;
0045
0046   // mass density is defined as
0047   //
0048   //     mass-density = atomic-mass * number-of-atoms / volume
0049   //                  = atomic-mass * molar-density * avogadro-constant
0050   // -> molar-density = mass-density / (atomic-mass * avogadro-constant)
0051   //
0052   // with the atomic mass given by
0053   //
0054   //      atomic-mass = relative-atomic-mass * atomic-mass-unit
0055   //
0056   // perform computations in double precision to avoid loss of precision
0057   const double atomicMass = static_cast<double>(ar) * 1_u;
0058   float molarRho = static_cast<float>(massRho / (atomicMass * kAvogadro));
0059
0060   return Material::fromMolarDensity(x0, l0, ar, z, molarRho);
0061 }
0062
0063 Material Material::fromMolarDensity(float x0, float l0, float ar, float z,
0064                                     float molarRho) {
0065   return Material::fromMolarDensity(x0, l0, ar, z, molarRho,
0066                                     calculateMolarElectronDensity(z, molarRho),
0067                                     std::nullopt);
0068 }
0069
0070 Material Material::fromMolarDensity(float x0, float l0, float ar, float z,
0071                                     float molarRho, float molarElectronRho,
0072                                     std::optional<float> meanExcitationEnergy) {
0073   Material mat;
0074   mat.m_x0 = x0;
0075   mat.m_l0 = l0;
0076   mat.m_ar = ar;
0077   mat.m_z = z;
0078   mat.m_molarRho = molarRho;
0079   mat.m_molarElectronRho = molarElectronRho;
0080   mat.m_meanExcitationEnergy =
0081       meanExcitationEnergy.value_or(approximateMeanExcitationEnergy(z));
0082   return mat;
0083 }
0084
0085 Material::Material(const ParametersVector& parameters)
0086     : m_x0(parameters[eRadiationLength]),
0087       m_l0(parameters[eInteractionLength]),
0088       m_ar(parameters[eRelativeAtomicMass]),
0089       m_z(parameters[eNuclearCharge]),
0090       m_molarRho(parameters[eMolarDensity]),
0091       m_molarElectronRho(calculateMolarElectronDensity(m_z, m_molarRho)),
0092       m_meanExcitationEnergy(approximateMeanExcitationEnergy(m_z)) {}
0093
0094 float Material::massDensity() const {
0095   using namespace UnitLiterals;
0096
0097   // perform computations in double precision to avoid loss of precision
0098   const double atomicMass = static_cast<double>(m_ar) * 1_u;
0099   const double numberDensity = static_cast<double>(m_molarRho) * kAvogadro;
0100   return atomicMass * numberDensity;
0101 }
0102
0103 Material::ParametersVector Material::parameters() const {
0104   ParametersVector parameters;
0105   parameters[eRadiationLength] = m_x0;
0106   parameters[eInteractionLength] = m_l0;
0107   parameters[eRelativeAtomicMass] = m_ar;
0108   parameters[eNuclearCharge] = m_z;
0109   parameters[eMolarDensity] = m_molarRho;
0110   return parameters;
0111 }
0112
0113 std::ostream& operator<<(std::ostream& os, const Material& material) {
0114   if (material.isVacuum()) {
0115     os << "vacuum";
0116   } else {
0117     os << "x0=" << material.X0();
0118     os << "|l0=" << material.L0();
0119     os << "|ar=" << material.Ar();
0120     os << "|z=" << material.Z();
0121     os << "|molar_rho=" << material.molarDensity();
0122     os << "|molar_e_rho=" << material.molarElectronDensity();
0123     os << "|mean_excitation_energy=" << material.meanExcitationEnergy();
0124   }
0125   return os;
0126 }
0127
0128 }  // namespace Acts