diff --git a/corsika/detail/modules/fluka/InteractionModel.inl b/corsika/detail/modules/fluka/InteractionModel.inl
index f7ae6fd321d7eeb355280bccffb5c40ebe5b7fd2..197adb54bb2856cc2deed53e7ad2a654a3f705a7 100644
--- a/corsika/detail/modules/fluka/InteractionModel.inl
+++ b/corsika/detail/modules/fluka/InteractionModel.inl
@@ -23,6 +23,7 @@
#include <corsika/framework/core/PhysicalUnits.hpp>
#include <FLUKA.hpp>
+#include <ParticleConversion.hpp>
namespace corsika::fluka {
template <typename TEnvironment>
@@ -32,11 +33,8 @@ namespace corsika::fluka {
}
}
- inline bool InteractionModel::isValid(Code projectileID, Code targetID /*, HEPEnergyType sqrtS*/) {
- static std::array<Code> constexpr validProjectiles{Code::Proton, Code::AntiProton, Code::Neutron, Code::AntiNeutron,
- Code::PiPlus, Code::PiMinus, Code::KPlus, Code::KMinus, Code::K0Long, Code::K0Short}
-
- if (std::find(validProjectiles.cbegin(), validProjectiles.cend(), projectileID) == validProjectiles.cend()) {
+ inline bool InteractionModel::isValid(Code projectileID, Code targetID, HEPEnergyType /*sqrtS*/) const {
+ if (fluka::canInteract(projectileID)) {
// invalid projectile
return false;
}
@@ -47,10 +45,11 @@ namespace corsika::fluka {
}
// TODO: check validity range of sqrtS
+ return true;
}
inline int InteractionModel::getMaterialIndex(Code targetID) const {
- if (auto it = std::find(materials_.cbegin(), materials_.cend(), [=targetID](std::pair<Code, int>& p) {return p.first == targetID;});
+ if (auto it = std::find(materials_.cbegin(), materials_.cend(), [=](std::pair<Code, int> const& p) {return p.first == targetID;});
it == materials_.cend()) {
return -1;
} else {
@@ -61,12 +60,12 @@ namespace corsika::fluka {
inline CrossSectionType InteractionModel::getCrossSection(
Code const projectileId, Code const targetId, FourMomentum const& projectileP4,
FourMomentum const& targetP4) const {
-
- if (!isValid(projectileId, targetId)) { return CrossSectionType::zero(); }
+ HEPEnergyType const sqrtS = (projectileP4 + targetP4).getNorm();
+ if (!isValid(projectileId, targetId, sqrtS)) { return CrossSectionType::zero(); }
- COMBoost const targetRestBoost{projecileP4.getSpaceLikeComponents(), get_mass(targetId)};
- FourMomentum const projectileLab4mom = targetRestBoost.toCOM(projectileP4);
- projec
+ COMBoost const targetRestBoost{projectileP4.getSpaceLikeComponents(), get_mass(targetId)};
+ FourMomentum const projectileLab4mom = targetRestBoost.toCoM(projectileP4);
+ //~ projec
int const iflxyz = 1;
diff --git a/modules/fluka/ParticleConversion.hpp b/modules/fluka/ParticleConversion.hpp
new file mode 100644
index 0000000000000000000000000000000000000000..64e97dea0b75d93a75c89060923617ac4674b176
--- /dev/null
+++ b/modules/fluka/ParticleConversion.hpp
@@ -0,0 +1,48 @@
+/*
+ * (c) Copyright 2023 CORSIKA Project, corsika-project@lists.kit.edu
+ *
+ * This software is distributed under the terms of the GNU General Public
+ * Licence version 3 (GPL Version 3). See file LICENSE for a full version of
+ * the license.
+ */
+
+#pragma once
+
+#include <array>
+#include <cstdint>
+#include <type_traits>
+#include <vector>
+
+#include <corsika/framework/core/ParticleProperties.hpp>
+#include <fmt/core.h>
+
+#include <string>
+
+namespace corsika::fluka {
+
+ enum class FLUKACode : int;
+ using FLUKACodeIntType = std::underlying_type<FLUKACode>::type;
+
+#include <corsika/modules/fluka/Generated.inc>
+
+ FLUKACode constexpr convertToFluka(Code const c8id) {
+ if (is_nucleus(c8id)) {
+ throw std::runtime_error{"nucleus conversion to FLUKA not implemented"};
+ }
+
+ FLUKACode const flukaID = corsika2fluka[static_cast<corsika::CodeIntType>(c8id)];
+ if (flukaID == FLUKACode::Unknown) {
+ throw std::runtime_error{fmt::format("no correspondig FLUKA id for {}", get_name(c8id)).c_str()};
+ }
+
+ return flukaID;
+ }
+
+ int constexpr convertToFlukaRaw(Code const code) {
+ return static_cast<FLUKACodeIntType>(convertToFluka(code));
+ }
+
+ bool const canInteract(Code const code) {
+ return flukaCanInteract[static_cast<CodeIntType>(code)];
+ }
+} // namespace corsika::fluka
diff --git a/src/CMakeLists.txt b/src/CMakeLists.txt
index f8d5df55a95bd8e4547d58322fecb6cc23e5f262..34fe8bdd3040ce796ec5a6c435c4fa5680f2c0e0 100644
--- a/src/CMakeLists.txt
+++ b/src/CMakeLists.txt
@@ -1,6 +1,7 @@
-add_subdirectory (framework/core)
+add_subdirectory (framework/core)
add_subdirectory (media)
-add_subdirectory (modules/sibyll)
-add_subdirectory (modules/sophia)
-add_subdirectory (modules/qgsjetII)
-add_subdirectory (modules/epos)
+add_subdirectory (modules/sibyll)
+add_subdirectory (modules/sophia)
+add_subdirectory (modules/qgsjetII)
+add_subdirectory (modules/epos)
+add_subdirectory (modules/fluka)
diff --git a/src/modules/fluka/code_generator.py b/src/modules/fluka/code_generator.py
index 690d1ca32f4a709cf0e2048703360e66e8325578..3574fa25756193f42f21166088e70506a457ede2 100755
--- a/src/modules/fluka/code_generator.py
+++ b/src/modules/fluka/code_generator.py
@@ -21,22 +21,23 @@ def load_particledb(filename):
return particle_db
-def read_epos_codes(filename, particle_db):
+def read_fluka_codes(filename, particle_db):
'''
- reads to epos codes data file
+ reads to FLUKA codes data file
- For particls known to epos, add 'epos_code' and 'epos_xsType' to particle_db
+ For particles that can interact in FLUKA, add can-interact flag
'''
with open(filename) as f:
for line in f:
line = line.strip()
if len(line) == 0 or line[0].startswith('#'):
continue
+ line = line.split()
if len(line) == 2:
canInteractFlag = "no"
- identifier, fluka_code = line.split()
+ identifier, fluka_code = line
else:
- identifier, fluka_code, canInteractFlag = line.split()
+ identifier, fluka_code, canInteractFlag = line
try:
particle_db[identifier]["fluka_code"] = int(fluka_code)
@@ -49,7 +50,8 @@ def generate_fluka_enum(particle_db):
'''
generates the enum to access epos particles by readable names
'''
- output = "enum class FLUKACode : int {\n"
+ output = ("enum class FLUKACode : int {\n"
+ " Unknown = 0,\n")
for identifier, pData in particle_db.items():
if 'fluka_code' in pData:
output += " {:s} = {:d},\n".format(identifier, pData['fluka_code'])
@@ -61,7 +63,7 @@ def generate_corsika2fluka(particle_db):
'''
generates the look-up table to convert corsika codes to epos codes
'''
- string = "std::array<FLUKACode, {:d}> constexpr corsika2epos = {{\n".format(len(particle_db))
+ string = "static inline std::array<FLUKACode, {:d}> constexpr corsika2fluka = {{\n".format(len(particle_db))
for identifier, pData in particle_db.items():
if pData['isNucleus']: continue
if 'fluka_code' in pData:
@@ -76,10 +78,10 @@ def generate_fluka_canInteract(particle_db):
'''
generates the look-up table to convert corsika codes to epos codes
'''
- string = "std::vector<bool> const flukaCanInteract = {\n"
+ string = "static inline std::vector<bool> const flukaCanInteract = {\n"
for identifier, pData in particle_db.items():
canInteract = pData.get("fluka_canInteract", False)
- string += " {:s}, // {:s}\n".format(str(canInteract).lower(), "identifier")
+ string += " {:s}, // {:s}\n".format(str(canInteract).lower(), identifier)
string += "};\n"
return string
@@ -89,11 +91,10 @@ if __name__ == "__main__":
print("usage: {:s} <particle_db.pkl> <fluka_codes.dat>".format(sys.argv[0]), file=sys.stderr)
sys.exit(1)
- print("code_generator.py for EPOS")
+ print("code_generator.py for FLUKA")
particle_db = load_particledb(sys.argv[1])
read_fluka_codes(sys.argv[2], particle_db)
- # ~ set_default_epos_definition(particle_db)
with open("Generated.inc", "w") as f:
print("// this file is automatically generated\n// edit at your own risk!\n", file=f)
diff --git a/src/modules/fluka/fluka_codes.dat b/src/modules/fluka/fluka_codes.dat
index e2c1c319ccd06f194f7a54f66ff10a2b3e8c5f59..71bad2741b1846f9765ba38847a3a508d054e835 100644
--- a/src/modules/fluka/fluka_codes.dat
+++ b/src/modules/fluka/fluka_codes.dat
@@ -14,7 +14,7 @@ Neutron 8 yes
Proton 1 yes
AntiProton 2 yes
K0Short 19 yes
-Eta 0
+# Eta 0
Lambda0 17 yes
SigmaPlus 21
Sigma0 22
@@ -31,48 +31,48 @@ Xi0Bar 34
XiPlusBar 37
OmegaPlusBar 39
-Omega 0
-Rho0 0
-RhoPlus 0
-RhoRhoMinus 0
-DeltaPlusPlus 0
-DeltaPlus 0
-Delta0 0
-DeltaMinus 0
-DeltaMinusMinusBar 0
-DeltaMinusBar 0
-Delta0Bar 0
-DeltaPlusBar 0
-KStar0 0
-KStarPlus 0
-KStarMinusBar 0
-KStar0Bar 0
+# Omega 0
+# Rho0 0
+# RhoPlus 0
+# RhoMinus 0
+# DeltaPlusPlus 0
+# DeltaPlus 0
+# Delta0 0
+# DeltaMinus 0
+# DeltaMinusMinusBar 0
+# DeltaMinusBar 0
+# Delta0Bar 0
+# DeltaPlusBar 0
+# KStar0 0
+# KStarPlus 0
+# KStarMinusBar 0
+# KStar0Bar 0
NuE 5
NuEBar 6
NuMu 27
-NuMubar 28
+NuMuBar 28
D0 16
DPlus 25
-DMinusBar 24
+DMinus 24
D0Bar 15
-DSPlus 0
-DSMinusBar 0
-EtaC 0
-DStar0 0
-DStarPlus 0
-DStarMinusBar 0
-DStar0Bar 0
-DStarSPlus 0
-DStarsMinusBar 0
- 0
-JPsi 0
+# DSPlus 0
+# DSMinusBar 0
+# EtaC 0
+# DStar0 0
+# DStarPlus 0
+# DStarMinusBar 0
+# DStar0Bar 0
+# DStarSPlus 0
+# DStarsMinusBar 0
+# 0
+# JPsi 0
TauPlus 41
TauMinus 42
NuTau 43
NuTauBar 44
- 0
- 0
+# 0
+# 0
LambdaCPlus 17
XiCPlus 34
XiC0 36
@@ -93,31 +93,31 @@ XiPrimeCMinusBar 35
XiPrimeC0Bar 37
OmegaC0Bar 39
-SigmaStarCPlusPlus 0
-SigmaStarCPlus 0
-SigmaStarC0 0
+# SigmaStarCPlusPlus 0
+# SigmaStarCPlus 0
+# SigmaStarC0 0
-SigmaStarCBar-- 0
-SigmaStarCMinusBar 0
-SigmaStarC0Bar 0
+# SigmaStarCBar-- 0
+# SigmaStarCMinusBar 0
+# SigmaStarC0Bar 0
-B0 0
-BPlus 0
-BMinusBar 0
-B0Bar 0
-BS0 0
-BS0Bar 0
-BCPlus 0
-BCMinusBar 0
-LambdaB0 0
-SigmaBMinus 0
-SigmaBPlus 0
-XiB0 0
-XiBMinus 0
-OmegaBMinus 0
-LambdaB0Bar 0
-SigmaBPlusBar 0
-SigmaBMinusBar 0
-XiB0Bar 0
-XiBPlusBar 0
-OmegaBPlusBar 0
+# B0 0
+# BPlus 0
+# BMinusBar 0
+# B0Bar 0
+# BS0 0
+# BS0Bar 0
+# BCPlus 0
+# BCMinusBar 0
+# LambdaB0 0
+# SigmaBMinus 0
+# SigmaBPlus 0
+# XiB0 0
+# XiBMinus 0
+# OmegaBMinus 0
+# LambdaB0Bar 0
+# SigmaBPlusBar 0
+# SigmaBMinusBar 0
+# XiB0Bar 0
+# XiBPlusBar 0
+# OmegaBPlusBar 0