IAP GITLAB

Skip to content
Snippets Groups Projects
Commit 680e3e6d authored by Maximilian Sackel's avatar Maximilian Sackel
Browse files

add type checking before processing particle with proposal

parent bfc401d5
No related branches found
No related tags found
No related merge requests found
......@@ -188,7 +188,7 @@ int main(int argc, char** argv) {
PROPOSAL::InterpolationDef::path_to_tables = "~/.local/share/PROPOSAL/tables/";
PROPOSAL::InterpolationDef::path_to_tables_readonly = "~/.local/share/PROPOSAL/tables/";
process::particle_cut::ParticleCut cut{60_GeV, true, true};
process::particle_cut::ParticleCut cut{60_GeV, false, true};
process::proposal::Interaction proposal(env, cut);
process::proposal::ContinuousProcess em_continuous(env, cut);
process::interaction_counter::InteractionCounter proposalCounted(proposal);
......@@ -205,9 +205,9 @@ int main(int argc, char** argv) {
process::UrQMD::UrQMD urqmd;
process::interaction_counter::InteractionCounter urqmdCounted{urqmd};
process::conex_source_cut::CONEXSourceCut conexSource(
center, showerAxis, t, injectionHeight, E0,
particles::GetPDG(particles::Code::Proton));
/* process::conex_source_cut::CONEXSourceCut conexSource( */
/* center, showerAxis, t, injectionHeight, E0, */
/* particles::GetPDG(particles::Code::Proton)); */
// assemble all processes into an ordered process list
......@@ -218,9 +218,9 @@ int main(int argc, char** argv) {
//auto sequence = switchProcess << reset_particle_mass << decaySequence << conexSource
// << longprof << eLoss << cut << observationLevel;
auto sequence = switchProcess << reset_particle_mass << decaySequence << proposalCounted << cut << em_continuous
auto sequence = switchProcess << reset_particle_mass << decaySequence << proposalCounted << cut << em_continuous
<< longprof << observationLevel;
// define air shower object, run simulation
tracking_line::TrackingLine tracking;
cascade::Cascade EAS(env, tracking, sequence, stack);
......@@ -232,7 +232,7 @@ int main(int argc, char** argv) {
EAS.Run();
eLoss.PrintProfile(); // print longitudinal profile
conexSource.SolveCE();
/* conexSource.SolveCE(); */
cut.ShowResults();
const HEPEnergyType Efinal =
......
......@@ -23,7 +23,10 @@ namespace corsika::process::proposal {
void ContinuousProcess::BuildCalculator(particles::Code code,
environment::NuclearComposition const& comp) {
auto c = cross[code](media.at(&comp), cut);
auto p_cross = cross.find(code);
if (p_cross == cross.end())
throw std::runtime_error("PROPOSAL could not find corresponding builder");
auto c = p_cross->second(media.at(&comp), cut);
auto disp = PROPOSAL::make_displacement(c, true);
auto scatter = PROPOSAL::make_scattering("highland", particle[code], media.at(&comp));
calc[std::make_pair(&comp, code)] =
......@@ -83,13 +86,16 @@ namespace corsika::process::proposal {
template <>
units::si::LengthType ContinuousProcess::MaxStepLength(
setup::Stack::ParticleType const& vP, setup::Trajectory const& vT) {
if (!CanInteract(vP.GetPID()))
return units::si::meter * std::numeric_limits<double>::infinity();
auto energy_lim = 0.9 * vP.GetEnergy();
if (cut.GetECut() > energy_lim) energy_lim = cut.GetECut();
auto c = GetCalculator(vP, calc);
auto grammage = get<DISPLACEMENT>(c->second)->SolveTrackIntegral(
vP.GetEnergy() / 1_MeV, energy_lim / 1_MeV) *
1_g / square(1_cm);
return vP.GetNode()->GetModelProperties().ArclengthFromGrammage(vT, grammage) * 1.0001;
return vP.GetNode()->GetModelProperties().ArclengthFromGrammage(vT, grammage) *
1.0001;
}
} // namespace corsika::process::proposal
......@@ -31,7 +31,10 @@ namespace corsika::process::proposal {
void Interaction::BuildCalculator(particles::Code code,
environment::NuclearComposition const& comp) {
auto c = cross[code](media.at(&comp), cut);
auto p_cross = cross.find(code);
if (p_cross == cross.end())
throw std::runtime_error("PROPOSAL could not find corresponding builder");
auto c = p_cross->second(media.at(&comp), cut);
auto inter_types = PROPOSAL::CrossSectionVector::GetInteractionTypes(c);
calc[std::make_pair(&comp, code)] = std::make_tuple(
PROPOSAL::make_secondaries(inter_types, particle[code], media.at(&comp)),
......
......@@ -49,7 +49,7 @@ namespace corsika::process::proposal {
}
PROPOSAL::InterpolationDef::order_of_interpolation = 2;
PROPOSAL::InterpolationDef::nodes_cross_section = 100;
PROPOSAL::InterpolationDef::nodes_propagate = 100;
PROPOSAL::InterpolationDef::nodes_propagate = 1000;
}
size_t ProposalProcessBase::hash::operator()(const calc_key_t& p) const noexcept {
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment