From ff5931d10f8e80832af0fdd12b2789e350815372 Mon Sep 17 00:00:00 2001
From: ralfulrich <ralf.ulrich@kit.edu>
Date: Fri, 2 Oct 2020 10:41:56 +0200
Subject: [PATCH] pass em cut energy in HEPEnergyUnits to PROPOSAL, not the
entire ParticleCut process.
---
Documentation/Examples/em_shower.cc | 4 ++--
Documentation/Examples/vertical_EAS.cc | 4 ++--
Processes/Proposal/ContinuousProcess.cc | 12 ++++++------
Processes/Proposal/ContinuousProcess.h | 3 +--
Processes/Proposal/Interaction.cc | 7 +++----
Processes/Proposal/Interaction.h | 3 +--
Processes/Proposal/ProposalProcessBase.cc | 4 ++--
Processes/Proposal/ProposalProcessBase.h | 17 ++++++++---------
8 files changed, 25 insertions(+), 29 deletions(-)
diff --git a/Documentation/Examples/em_shower.cc b/Documentation/Examples/em_shower.cc
index e84931014..9bb1daf3f 100644
--- a/Documentation/Examples/em_shower.cc
+++ b/Documentation/Examples/em_shower.cc
@@ -131,8 +131,8 @@ int main(int argc, char** argv) {
// PROPOSAL processs proposal{...};
process::particle_cut::ParticleCut cut(10_GeV, false, true);
- process::proposal::Interaction proposal(env, cut);
- process::proposal::ContinuousProcess em_continuous(env, cut);
+ process::proposal::Interaction proposal(env, cut.GetECut());
+ process::proposal::ContinuousProcess em_continuous(env, cut.GetECut());
process::interaction_counter::InteractionCounter proposalCounted(proposal);
process::track_writer::TrackWriter trackWriter("tracks.dat");
diff --git a/Documentation/Examples/vertical_EAS.cc b/Documentation/Examples/vertical_EAS.cc
index 5f7f7d343..b639f6374 100644
--- a/Documentation/Examples/vertical_EAS.cc
+++ b/Documentation/Examples/vertical_EAS.cc
@@ -185,8 +185,8 @@ int main(int argc, char** argv) {
// PROPOSAL processs proposal{...};
process::particle_cut::ParticleCut cut{60_GeV, false, true};
- process::proposal::Interaction proposal(env, cut);
- process::proposal::ContinuousProcess em_continuous(env, cut);
+ process::proposal::Interaction proposal(env, cut.GetECut());
+ process::proposal::ContinuousProcess em_continuous(env, cut.GetECut());
process::interaction_counter::InteractionCounter proposalCounted(proposal);
process::on_shell_check::OnShellCheck reset_particle_mass(1.e-3, 1.e-1, false);
diff --git a/Processes/Proposal/ContinuousProcess.cc b/Processes/Proposal/ContinuousProcess.cc
index ebd0ba087..0337bc8b0 100644
--- a/Processes/Proposal/ContinuousProcess.cc
+++ b/Processes/Proposal/ContinuousProcess.cc
@@ -8,7 +8,6 @@
#include <PROPOSAL/PROPOSAL.h>
#include <corsika/environment/IMediumModel.h>
-#include <corsika/process/particle_cut/ParticleCut.h>
#include <corsika/process/proposal/ContinuousProcess.h>
#include <corsika/process/proposal/Interaction.h>
#include <corsika/setup/SetupEnvironment.h>
@@ -33,7 +32,7 @@ namespace corsika::process::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 c = p_cross->second(media.at(&comp), cut);
+ auto c = p_cross->second(media.at(&comp), emCut_);
// Build displacement integral and scattering object and interpolate them too and
// saved in the calc map by a key build out of a hash of composed of the component and
@@ -46,8 +45,8 @@ namespace corsika::process::proposal {
template <>
ContinuousProcess::ContinuousProcess(setup::SetupEnvironment const& _env,
- particle_cut::ParticleCut& _cut)
- : ProposalProcessBase(_env, _cut) {}
+ corsika::units::si::HEPEnergyType _emCut)
+ : ProposalProcessBase(_env, _emCut) {}
template <>
void ContinuousProcess::Scatter(setup::Stack::ParticleType& vP,
@@ -108,7 +107,7 @@ namespace corsika::process::proposal {
// Update the energy and absorbe the particle if it's below the energy
// threshold, because it will no longer propagated.
vP.SetEnergy(final_energy);
- if (vP.GetEnergy() <= cut.GetECut())
+ if (vP.GetEnergy() <= emCut_)
return process::EProcessReturn::eParticleAbsorbed;
vP.SetMomentum(vP.GetMomentum() * vP.GetEnergy() / vP.GetMomentum().GetNorm());
return process::EProcessReturn::eOk;
@@ -122,7 +121,8 @@ namespace corsika::process::proposal {
// Limit the step size of a conitnous loss. The maximal continuous loss seems to be a
// hyper parameter which must be adjusted.
- auto energy_lim = 0.9 * vP.GetEnergy();
+ // in any case: never go below emCut_
+ auto energy_lim = std::max(0.9 * vP.GetEnergy(), emCut_);;
// solving the track integral for giving energy lim
auto c = GetCalculator(vP, calc);
diff --git a/Processes/Proposal/ContinuousProcess.h b/Processes/Proposal/ContinuousProcess.h
index 6a3dd177d..e2624ed48 100644
--- a/Processes/Proposal/ContinuousProcess.h
+++ b/Processes/Proposal/ContinuousProcess.h
@@ -12,7 +12,6 @@
#include <corsika/environment/Environment.h>
#include <corsika/particles/ParticleProperties.h>
#include <corsika/process/ContinuousProcess.h>
-#include <corsika/process/particle_cut/ParticleCut.h>
#include <corsika/process/proposal/ProposalProcessBase.h>
#include <corsika/random/RNGManager.h>
#include <corsika/random/UniformRealDistribution.h>
@@ -49,7 +48,7 @@ namespace corsika::process::proposal {
//! compositions and stochastic description limited by the particle cut.
//!
template <typename TEnvironment>
- ContinuousProcess(TEnvironment const&, particle_cut::ParticleCut&);
+ ContinuousProcess(TEnvironment const&, corsika::units::si::HEPEnergyType _emCut);
//!
diff --git a/Processes/Proposal/Interaction.cc b/Processes/Proposal/Interaction.cc
index b88653360..aeeaae484 100644
--- a/Processes/Proposal/Interaction.cc
+++ b/Processes/Proposal/Interaction.cc
@@ -8,7 +8,6 @@
#include <corsika/environment/IMediumModel.h>
#include <corsika/environment/NuclearComposition.h>
-#include <corsika/process/particle_cut/ParticleCut.h>
#include <corsika/process/proposal/Interaction.h>
#include <corsika/setup/SetupEnvironment.h>
#include <corsika/setup/SetupStack.h>
@@ -26,8 +25,8 @@ namespace corsika::process::proposal {
template <>
Interaction::Interaction(setup::SetupEnvironment const& _env,
- particle_cut::ParticleCut& _cut)
- : ProposalProcessBase(_env, _cut) {}
+ corsika::units::si::HEPEnergyType _emCut)
+ : ProposalProcessBase(_env, _emCut) {}
void Interaction::BuildCalculator(particles::Code code,
environment::NuclearComposition const& comp) {
@@ -39,7 +38,7 @@ namespace corsika::process::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 c = p_cross->second(media.at(&comp), cut);
+ auto c = p_cross->second(media.at(&comp), emCut_);
// Look which interactions take place and build the corresponding
// interaction and secondarie builder. The interaction integral will
diff --git a/Processes/Proposal/Interaction.h b/Processes/Proposal/Interaction.h
index c96cf00df..0a4aed0a4 100644
--- a/Processes/Proposal/Interaction.h
+++ b/Processes/Proposal/Interaction.h
@@ -13,7 +13,6 @@
#include <corsika/environment/Environment.h>
#include <corsika/particles/ParticleProperties.h>
#include <corsika/process/InteractionProcess.h>
-#include <corsika/process/particle_cut/ParticleCut.h>
#include <corsika/process/proposal/ProposalProcessBase.h>
#include <corsika/random/RNGManager.h>
#include <corsika/random/UniformRealDistribution.h>
@@ -48,7 +47,7 @@ namespace corsika::process::proposal {
//! compositions and stochastic description limited by the particle cut.
//!
template <typename TEnvironment>
- Interaction(TEnvironment const& env, particle_cut::ParticleCut& cut);
+ Interaction(TEnvironment const& env, corsika::units::si::HEPEnergyType emCut);
//!
//! Calculate the rates for the different targets and interactions. Sample a
diff --git a/Processes/Proposal/ProposalProcessBase.cc b/Processes/Proposal/ProposalProcessBase.cc
index b24143570..0211b4788 100644
--- a/Processes/Proposal/ProposalProcessBase.cc
+++ b/Processes/Proposal/ProposalProcessBase.cc
@@ -27,8 +27,8 @@ namespace corsika::process::proposal {
}
ProposalProcessBase::ProposalProcessBase(setup::SetupEnvironment const& _env,
- particle_cut::ParticleCut& _cut)
- : cut(_cut)
+ corsika::units::si::HEPEnergyType _emCut)
+ : emCut_(_emCut)
, fRNG(corsika::random::RNGManager::GetInstance().GetRandomStream("proposal")) {
auto all_compositions = std::vector<const environment::NuclearComposition*>();
_env.GetUniverse()->walk([&](auto& vtn) {
diff --git a/Processes/Proposal/ProposalProcessBase.h b/Processes/Proposal/ProposalProcessBase.h
index b2f503e77..018fd1b4f 100644
--- a/Processes/Proposal/ProposalProcessBase.h
+++ b/Processes/Proposal/ProposalProcessBase.h
@@ -2,10 +2,9 @@
#include <PROPOSAL/PROPOSAL.h>
-#include <corsika/environment/Environment.h>
#include <corsika/particles/ParticleProperties.h>
-#include <corsika/process/particle_cut/ParticleCut.h>
#include <corsika/random/RNGManager.h>
+#include <corsika/setup/SetupEnvironment.h>
#include <array>
namespace corsika::process::proposal {
@@ -41,9 +40,9 @@ namespace corsika::process::proposal {
//! Crosssection factories for different particle types.
//!
template <typename T>
- static auto cross_builder = [](PROPOSAL::Medium& m, particle_cut::ParticleCut& cut) {
+ static auto cross_builder = [](PROPOSAL::Medium& m, corsika::units::si::HEPEnergyType emCut) {
auto p_cut = std::make_shared<const PROPOSAL::EnergyCutSettings>(
- cut.GetECut() / 1_MeV, 1, true);
+ emCut / 1_MeV, 1, true);
return PROPOSAL::DefaultCrossSections<T>::template Get<std::false_type>(T(), m, p_cut,
true);
};
@@ -54,7 +53,7 @@ namespace corsika::process::proposal {
//!
static std::map<particles::Code, std::function<PROPOSAL::crosssection_list_t<
PROPOSAL::ParticleDef, PROPOSAL::Medium>(
- PROPOSAL::Medium&, particle_cut::ParticleCut&)>>
+ PROPOSAL::Medium&, corsika::units::si::HEPEnergyType)>>
cross = {{particles::Code::Gamma, cross_builder<PROPOSAL::GammaDef>},
{particles::Code::Electron, cross_builder<PROPOSAL::EMinusDef>},
{particles::Code::Positron, cross_builder<PROPOSAL::EPlusDef>},
@@ -69,8 +68,8 @@ namespace corsika::process::proposal {
//!
class ProposalProcessBase {
protected:
- particle_cut::ParticleCut& cut; //!< Stochastic losses smaller than the given cut
- //!< will be handeled continuously.
+ corsika::units::si::HEPEnergyType emCut_; //!< Stochastic losses smaller than the given cut
+ //!< will be handeled continuously.
corsika::random::RNG& fRNG; //!< random number generator used by proposal
std::unordered_map<const environment::NuclearComposition*, PROPOSAL::Medium>
@@ -81,8 +80,8 @@ namespace corsika::process::proposal {
//! Store cut and nuclear composition of the whole universe in media which are
//! required for creating crosssections by proposal.
//!
- ProposalProcessBase(setup::SetupEnvironment const& _env,
- particle_cut::ParticleCut& _cut);
+ ProposalProcessBase(corsika::setup::SetupEnvironment const& _env,
+ corsika::units::si::HEPEnergyType _emCut);
//!
//! Checks if a particle can be processed by proposal
--
GitLab