diff --git a/applications/CMakeLists.txt b/applications/CMakeLists.txt
index b9d47c0c96c056daf3199388d16d91c8fd7e1600..a7231f257d98c911ff6094b9d26755258ff75d0e 100644
--- a/applications/CMakeLists.txt
+++ b/applications/CMakeLists.txt
@@ -13,3 +13,10 @@ endif()
install (
TARGETS c8_air_shower DESTINATION bin
)
+
+add_executable (cross_section cross_section.cpp)
+target_link_libraries (cross_section CORSIKA8)
+
+install (
+ TARGETS cross_section DESTINATION bin
+)
diff --git a/applications/cross_section.cpp b/applications/cross_section.cpp
new file mode 100644
index 0000000000000000000000000000000000000000..1eb34c0388b45752c75e6d55d30a22fc5dc90b2f
--- /dev/null
+++ b/applications/cross_section.cpp
@@ -0,0 +1,153 @@
+/*
+ * (c) Copyright 2019 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.
+ */
+
+#include <corsika/modules/Epos.hpp>
+#include <corsika/modules/Sibyll.hpp>
+#include <corsika/modules/QGSJetII.hpp>
+#include <corsika/modules/Pythia8.hpp>
+
+#include <corsika/framework/core/ParticleProperties.hpp>
+#include <corsika/framework/core/PhysicalUnits.hpp>
+
+//#include <corsika/framework/geometry/Point.hpp>
+#include <corsika/framework/geometry/RootCoordinateSystem.hpp>
+#include <corsika/framework/geometry/Vector.hpp>
+
+#include <corsika/framework/random/RNGManager.hpp>
+
+// #include <corsika/media/Environment.hpp>
+// #include <corsika/media/HomogeneousMedium.hpp>
+// #include <corsika/media/MediumPropertyModel.hpp>
+// #include <corsika/media/NuclearComposition.hpp>
+// #include <corsika/media/ShowerAxis.hpp>
+// #include <corsika/media/UniformMagneticField.hpp>
+
+// #include <corsika/setup/SetupStack.hpp>
+
+#include <fstream>
+#include <string>
+#include <tuple>
+
+/*
+ NOTE, WARNING, ATTENTION
+
+ The Random.hpp implements the hook of epos to the C8 random
+ number generator. It has to occur excatly ONCE per linked
+ executable. If you include the header below in multiple "tests" and
+ link them togehter, it will fail.
+ */
+#include <corsika/modules/Random.hpp>
+
+using namespace corsika;
+// using namespace corsika::epos;
+
+template <typename TModel>
+void calculate_cross_sections(TModel & model, std::string const &model_name,
+ bool lab_frame = true) {
+
+ std::stringstream fname;
+ fname << "cross-sections-" << model_name << ".txt";
+ std::ofstream out(fname.str());
+ out << "# cross section table for " << model_name << "\n";
+ out << "# i, elab, sqsnn, sig_pp, sig_pO \n";
+
+ CoordinateSystemPtr const &cs = get_root_CoordinateSystem();
+
+ float const amin = 2;
+ float const amax = 12;
+ int const nebins = 15;
+ float const da = (amax - amin) / (float)nebins;
+
+ for (int i = 0; i < nebins; ++i) {
+
+ corsika::units::si::HEPMomentumType const setP =
+ pow(10, da * i + amin) * 1_GeV;
+
+ auto setE = calculate_total_energy(setP, Proton::mass);
+ corsika::units::si::HEPEnergyType const comEnn =
+ corsika::calculate_com_energy(setE, Proton::mass,
+ corsika::constants::nucleonMass);
+
+ //corsika::units::si::HEPMomentumType setP =
+ // corsika::calculate_momentum(setE, Proton::mass);
+
+ FourMomentum p4protonProj(0_eV, {cs, 0_GeV, 0_GeV, 0_GeV});
+ FourMomentum p4protonTarg(0_eV, {cs, 0_GeV, 0_GeV, 0_GeV});
+ FourMomentum p4oxygenTarg(0_eV, {cs, 0_GeV, 0_GeV, 0_GeV});
+
+ if (lab_frame) {
+ // four momenta in lab frame
+ p4protonProj =
+ corsika::FourMomentum(setE, MomentumVector(cs, {0_eV, 0_eV, setP}));
+ p4protonTarg = corsika::FourMomentum(
+ Proton::mass, MomentumVector(cs, {0_eV, 0_eV, 0_eV}));
+ p4oxygenTarg = corsika::FourMomentum(
+ Oxygen::mass, MomentumVector(cs, {0_eV, 0_eV, 0_eV}));
+ } else {
+ // four momenta in com frame
+ p4protonProj = corsika::FourMomentum(
+ corsika::calculate_total_energy(setP / 2, Proton::mass),
+ MomentumVector(cs, {0_eV, 0_eV, setP / 2}));
+ p4protonTarg = corsika::FourMomentum(
+ corsika::calculate_total_energy(setP / 2, Proton::mass),
+ MomentumVector(cs, {0_eV, 0_eV, -setP / 2}));
+ p4oxygenTarg = corsika::FourMomentum(
+ corsika::calculate_total_energy(setP / 2, Oxygen::mass),
+ MomentumVector(cs, {0_eV, 0_eV, -setP / 2}));
+ }
+
+ // p-p
+ CrossSectionType const xs_prod_pp = model->getCrossSection(
+ Code::Proton, Code::Proton, p4protonProj, p4protonTarg);
+ // p-O
+ CrossSectionType const xs_prod_pO = model->getCrossSection(
+ Code::Proton, Code::Oxygen, p4protonProj, p4oxygenTarg);
+
+ CORSIKA_LOG_INFO("Elab={:15.2f} GeV,Ecom={:15.2f} GeV, sig-p-p={:6.2f} mb, sig-p-O={:6.2f} mb",
+ setE / 1_GeV, comEnn / 1_GeV, xs_prod_pp / 1_mb,
+ xs_prod_pO / 1_mb);
+
+ out << i << " " << setE / 1_GeV << " " << comEnn / 1_GeV << " "
+ << xs_prod_pp / 1_mb << " " << xs_prod_pO / 1_mb << "\n";
+ }
+ out.close();
+}
+
+int main(int argc, char **argv) {
+
+ logging::set_level(logging::level::info);
+
+ if (argc != 2) {
+ std::cout
+ << "usage: check <interaction model> \n valid models are: sibyll, "
+ "epos, qgsII, pythia8"
+ << std::endl;
+ return 1;
+ }
+ std::string int_model_name = std::string(argv[1]);
+
+ if (int_model_name == "sibyll") {
+ RNGManager<>::getInstance().registerRandomStream("sibyll");
+ auto model = std::make_shared<corsika::sibyll::InteractionModel>(std::set<Code>{Code::Proton, Code::Oxygen});
+ calculate_cross_sections(model, int_model_name);
+ } else if (int_model_name == "pythia8") {
+ RNGManager<>::getInstance().registerRandomStream("pythia");
+ auto model = std::make_shared<corsika::pythia8::InteractionModel>();
+ calculate_cross_sections(model, int_model_name);
+
+ } else if (int_model_name == "epos") {
+ RNGManager<>::getInstance().registerRandomStream("epos");
+ auto model = std::make_shared<corsika::epos::InteractionModel>();
+ calculate_cross_sections(model, int_model_name);
+ }
+ else if (int_model_name == "qgsjet") {
+ RNGManager<>::getInstance().registerRandomStream("qgsjet");
+ auto model = std::make_shared<corsika::qgsjetII::InteractionModel>();
+ calculate_cross_sections(model, int_model_name);
+ }
+}