Warning, file /acts/Core/src/Definitions/ParticleData.cpp was not indexed
or was modified since last indexation (in which case cross-reference links may be missing, inaccurate or erroneous).
0001
0002
0003
0004
0005
0006
0007
0008
0009 #include "Acts/Definitions/ParticleData.hpp"
0010
0011 #include "Acts/Definitions/PdgParticle.hpp"
0012 #include "Acts/Definitions/Units.hpp"
0013
0014 #include <algorithm>
0015 #include <array>
0016 #include <cassert>
0017 #include <cmath>
0018 #include <cstdint>
0019 #include <iterator>
0020 #include <limits>
0021 #include <optional>
0022 #include <ostream>
0023 #include <type_traits>
0024 #include <utility>
0025
0026 #include "ParticleDataTable.hpp"
0027
0028 namespace {
0029
0030 enum class Type {
0031 Charge,
0032 Mass,
0033 Name,
0034 };
0035
0036
0037
0038 template <typename T, Acts::PdgParticle pdg, Type type>
0039 std::optional<T> findCachedImpl(const std::map<std::int32_t, T>& map) {
0040 const static std::optional<T> value = [&map]() -> std::optional<T> {
0041 const auto it = map.find(pdg);
0042 if (it == map.end()) {
0043 return std::nullopt;
0044 }
0045 return it->second;
0046 }();
0047
0048 return value;
0049 }
0050
0051
0052
0053 template <typename T, Type type>
0054 std::optional<T> findCached(Acts::PdgParticle pdg,
0055 const std::map<std::int32_t, T>& map) {
0056 using enum Acts::PdgParticle;
0057 switch (pdg) {
0058 case eElectron:
0059 return findCachedImpl<T, eElectron, type>(map);
0060 case ePositron:
0061 return findCachedImpl<T, ePositron, type>(map);
0062 case eMuon:
0063 return findCachedImpl<T, eMuon, type>(map);
0064 case eAntiMuon:
0065 return findCachedImpl<T, eAntiMuon, type>(map);
0066 case eTau:
0067 return findCachedImpl<T, eTau, type>(map);
0068 case eAntiTau:
0069 return findCachedImpl<T, eAntiTau, type>(map);
0070 case eGamma:
0071 return findCachedImpl<T, eGamma, type>(map);
0072 case ePionZero:
0073 return findCachedImpl<T, ePionZero, type>(map);
0074 case ePionPlus:
0075 return findCachedImpl<T, ePionPlus, type>(map);
0076 case ePionMinus:
0077 return findCachedImpl<T, ePionMinus, type>(map);
0078 case eKaonPlus:
0079 return findCachedImpl<T, eKaonPlus, type>(map);
0080 case eKaonMinus:
0081 return findCachedImpl<T, eKaonMinus, type>(map);
0082 case eNeutron:
0083 return findCachedImpl<T, eNeutron, type>(map);
0084 case eAntiNeutron:
0085 return findCachedImpl<T, eAntiNeutron, type>(map);
0086 case eProton:
0087 return findCachedImpl<T, eProton, type>(map);
0088 case eAntiProton:
0089 return findCachedImpl<T, eAntiProton, type>(map);
0090 case eLead:
0091 return findCachedImpl<T, eLead, type>(map);
0092 default:
0093 return std::nullopt;
0094 }
0095 }
0096
0097 }
0098
0099 std::optional<float> Acts::findCharge(Acts::PdgParticle pdg) {
0100 if (auto cached = findCached<float, Type::Charge>(pdg, kParticlesMapCharge);
0101 cached) {
0102 return cached;
0103 }
0104
0105 const auto it = kParticlesMapCharge.find(pdg);
0106 if (it == kParticlesMapCharge.end()) {
0107 if (isNucleus(pdg)) {
0108 return Acts::findChargeOfNucleus(pdg);
0109 } else {
0110 return std::nullopt;
0111 }
0112 }
0113 return it->second;
0114 }
0115
0116 float Acts::findChargeOfNucleus(Acts::PdgParticle pdg) {
0117 if (!isNucleus(pdg)) {
0118 throw std::invalid_argument("PDG must represent a nucleus");
0119 }
0120 auto pdgGround = makeNucleusGroundState(pdg);
0121 const auto it = kParticlesMapCharge.find(pdgGround);
0122 if (it == kParticlesMapCharge.end()) {
0123
0124 auto ZandA = extractNucleusZandA(pdg);
0125 return static_cast<float>(ZandA.first);
0126 }
0127 return it->second;
0128 }
0129
0130 std::optional<float> Acts::findMass(Acts::PdgParticle pdg) {
0131 if (auto cached = findCached<float, Type::Mass>(pdg, kParticlesMapMass);
0132 cached) {
0133 return cached;
0134 }
0135
0136 const auto it = kParticlesMapMass.find(pdg);
0137 if (it == kParticlesMapMass.end()) {
0138 if (isNucleus(pdg)) {
0139 return Acts::findMassOfNucleus(pdg);
0140 } else {
0141 return std::nullopt;
0142 }
0143 }
0144 return it->second;
0145 }
0146
0147 float Acts::findMassOfNucleus(Acts::PdgParticle pdg) {
0148 if (!isNucleus(pdg)) {
0149 throw std::invalid_argument("PDG must represent a nucleus");
0150 }
0151 auto pdgGround = makeNucleusGroundState(pdg);
0152 const auto it = kParticlesMapMass.find(pdgGround);
0153 if (it == kParticlesMapMass.end()) {
0154
0155 return calculateNucleusMass(pdg);
0156 }
0157 return it->second;
0158 }
0159
0160 float Acts::calculateNucleusMass(Acts::PdgParticle pdg) {
0161 if (!isNucleus(pdg)) {
0162 throw std::invalid_argument("PDG must represent a nucleus");
0163 }
0164
0165 auto ZandA = extractNucleusZandA(pdg);
0166 int Z = std::abs(ZandA.first);
0167 int A = ZandA.second;
0168
0169 float a_Vol = 15.260f * Acts::UnitConstants::MeV;
0170 float a_Surf = 16.267f * Acts::UnitConstants::MeV;
0171 float a_Col = 0.689f * Acts::UnitConstants::MeV;
0172 float a_Sym = 22.209f * Acts::UnitConstants::MeV;
0173 float a_Pair =
0174 (1 - 2 * (Z % 2)) * (1 - A % 2) * 10.076f * Acts::UnitConstants::MeV;
0175
0176 float massP = 0.938272f * Acts::UnitConstants::GeV;
0177 float massN = 0.939565f * Acts::UnitConstants::GeV;
0178
0179 float bindEnergy = 0.f;
0180 bindEnergy += a_Vol * A;
0181 bindEnergy -= a_Surf * std::pow(A, 2. / 3.);
0182 bindEnergy -= a_Col * Z * (Z - 1) * std::pow(A, -1. / 3.);
0183 bindEnergy -= a_Sym * std::pow(A - 2 * Z, 2.) / A;
0184 bindEnergy -= a_Pair * std::pow(A, -1. / 2.);
0185
0186 return massP * Z + massN * (A - Z) - bindEnergy;
0187 }
0188
0189 std::optional<std::string_view> Acts::findName(Acts::PdgParticle pdg) {
0190 if (auto cached =
0191 findCached<const char* const, Type::Name>(pdg, kParticlesMapName);
0192 cached) {
0193 return cached;
0194 }
0195
0196 const auto it = kParticlesMapName.find(pdg);
0197 if (it == kParticlesMapName.end()) {
0198 if (isNucleus(pdg)) {
0199 return Acts::findNameOfNucleus(pdg);
0200 } else {
0201 return std::nullopt;
0202 }
0203 }
0204 return it->second;
0205 }
0206
0207 std::optional<std::string_view> Acts::findNameOfNucleus(Acts::PdgParticle pdg) {
0208 if (!isNucleus(pdg)) {
0209 throw std::invalid_argument("PDG must represent a nucleus");
0210 }
0211
0212 auto pdgGround = makeNucleusGroundState(pdg);
0213 const auto it = kParticlesMapName.find(pdgGround);
0214 if (it == kParticlesMapName.end()) {
0215 return std::nullopt;
0216 }
0217 return it->second;
0218 }
0219
0220 std::optional<Acts::ParticleData> Acts::findParticleData(PdgParticle pdg) {
0221 auto charge = findCharge(pdg);
0222 auto mass = findMass(pdg);
0223 auto name = findName(pdg);
0224
0225 if (charge.has_value() && mass.has_value() && name.has_value()) {
0226 return ParticleData{charge.value(), mass.value(), name.value()};
0227 }
0228
0229 return std::nullopt;
0230 }
0231
0232 std::ostream& Acts::operator<<(std::ostream& os, Acts::PdgParticle pdg) {
0233 const auto name = Acts::findName(pdg);
0234 if (name) {
0235 os << *name;
0236 } else {
0237 os << static_cast<std::int32_t>(pdg);
0238 }
0239 return os;
0240 }
0241
0242 std::optional<std::string_view> Acts::pdgToShortAbsString(PdgParticle pdg) {
0243 pdg = makeAbsolutePdgParticle(pdg);
0244 if (pdg == eElectron) {
0245 return "e";
0246 }
0247 if (pdg == eMuon) {
0248 return "mu";
0249 }
0250 if (pdg == eTau) {
0251 return "t";
0252 }
0253 if (pdg == eGamma) {
0254 return "g";
0255 }
0256 if (pdg == ePionZero) {
0257 return "pi0";
0258 }
0259 if (pdg == ePionPlus) {
0260 return "pi";
0261 }
0262 if (pdg == eKaonPlus) {
0263 return "K";
0264 }
0265 if (pdg == eNeutron) {
0266 return "n";
0267 }
0268 if (pdg == eProton) {
0269 return "p";
0270 }
0271 if (pdg == eLead) {
0272 return "lead";
0273 }
0274 return std::nullopt;
0275 }