From d34cdcc8eb6e4a98533a6505cd59fb552f48f1ad Mon Sep 17 00:00:00 2001
From: Felix Riehn <felix@matilda>
Date: Fri, 30 Apr 2021 18:53:49 +0100
Subject: [PATCH] change threshold to kinetic energy for proposal and energy
 loss modules

---
 .../modules/energy_loss/BetheBlochPDG.inl       | 17 +++++++++--------
 .../modules/proposal/ContinuousProcess.inl      |  5 +++--
 corsika/detail/modules/proposal/Interaction.inl |  5 +++--
 3 files changed, 15 insertions(+), 12 deletions(-)

diff --git a/corsika/detail/modules/energy_loss/BetheBlochPDG.inl b/corsika/detail/modules/energy_loss/BetheBlochPDG.inl
index 310613d24..693814c4d 100644
--- a/corsika/detail/modules/energy_loss/BetheBlochPDG.inl
+++ b/corsika/detail/modules/energy_loss/BetheBlochPDG.inl
@@ -180,14 +180,15 @@ namespace corsika {
     auto constexpr dX = 1_g / square(1_cm);
     auto const dEdX = -getTotalEnergyLoss(vParticle, dX) / dX;
     auto const energy = vParticle.getEnergy();
-    auto const energy_lim = std::max(
-        energy * 0.9,                            // either 10% relative loss max., or
-        get_energy_threshold(vParticle.getPID()) // energy thresholds globally defined for
-                                                 // individual particles
-            * 0.99999 // need to go slightly below global e-cut to assure removal in
-                      // ParticleCut. This does not matter since at cut-time the entire
-                      // energy is removed.
-    );
+    auto const energy_lim =
+        std::max(energy * 0.9, // either 10% relative loss max., or
+                 get_kinetic_energy_threshold(vParticle.getPID()) +
+                     get_mass(vParticle.getPID()) // energy thresholds globally defined
+                                                  // for individual particles
+                         * 0.99999 // need to go slightly below global e-cut to assure removal in
+                                // ParticleCut. The 1% does not matter since at cut-time
+                                // the entire energy is removed.
+        );
     auto const maxGrammage = (vParticle.getEnergy() - energy_lim) / dEdX;
 
     return vParticle.getNode()->getModelProperties().getArclengthFromGrammage(
diff --git a/corsika/detail/modules/proposal/ContinuousProcess.inl b/corsika/detail/modules/proposal/ContinuousProcess.inl
index 0e230a762..c22226ac5 100644
--- a/corsika/detail/modules/proposal/ContinuousProcess.inl
+++ b/corsika/detail/modules/proposal/ContinuousProcess.inl
@@ -33,8 +33,9 @@ namespace corsika::proposal {
     // interpolate the crosssection for given media and energy cut. These may
     // take some minutes if you have to build the tables and cannot read the
     // from disk
-    auto const emCut = get_energy_threshold(
-        code); //! energy thresholds globally defined for individual particles
+    auto const emCut =
+        get_kinetic_energy_threshold(code) +
+        get_mass(code); //! energy thresholds globally defined for individual particles
     auto c = p_cross->second(media.at(comp.getHash()), emCut);
 
     // Build displacement integral and scattering object and interpolate them too and
diff --git a/corsika/detail/modules/proposal/Interaction.inl b/corsika/detail/modules/proposal/Interaction.inl
index fb6c157da..868b00810 100644
--- a/corsika/detail/modules/proposal/Interaction.inl
+++ b/corsika/detail/modules/proposal/Interaction.inl
@@ -36,8 +36,9 @@ namespace corsika::proposal {
     // interpolate the crosssection for given media and energy cut. These may
     // take some minutes if you have to build the tables and cannot read the
     // from disk
-    auto const emCut = get_energy_threshold(
-        code); //! energy thresholds globally defined for individual particles
+    auto const emCut =
+        get_kinetic_energy_threshold(code) +
+        get_mass(code); //! energy thresholds globally defined for individual particles
 
     auto c = p_cross->second(media.at(comp.getHash()), emCut);
 
-- 
GitLab