CONTRIBUTING.md

 # Guidelines for code development, structure, formating etc.
 
-The CORSIKA Project very much welcomes contributions. Here we outlined
+The CORSIKA Project very much welcomes contributions. Here we outline
 how you can find the right place to contribute, and how to do that.
-Connect to https://gitlab.ikp.kit.edu and corsika-devel@lists.kit.edu (self-register at https://www.lists.kit.edu/sympa/subscribe/corsika-devel) to get in touch with the project.
-The CORSIKA Project decides on the [GUIDELINES](CONTRIBUTING.md) and can decide to
+Connect to https://gitlab.iap.kit.edu and corsika-devel@lists.kit.edu (self-register at https://www.lists.kit.edu/sympa/subscribe/corsika-devel) to get in touch with the project.
+The CORSIKA Project decides on the [contributing guidelines](CONTRIBUTING.md) and can decide to
 change/improve them.
 
 # How to contribute
 
-  - We organize all development via `Issues` that may be feature requests,
+  - We organize all development via gitlab `Issues` that may be feature requests,
     ideas, discussions, or bugs fix requests. 
   - New issues can be created, or existing issues
     picked up or contributed to. 
@@ -19,7 +19,8 @@ change/improve them.
     created directly via the gitlab web interface. 
   - Proposed code to close one issue (located in a specific git
     branch) is reviewed, discussed, and eventually merged
-    into the master branch to close the issue.
+    into the master branch via a merge-request (MR) to close the issue.
+  - all merge request will undergo a code review, and must be approved before merge, in order to ensure high code qualtiy: [Code Approval Procedure](https://gitlab.iap.kit.edu/AirShowerPhysics/corsika/-/wikis/Code-Approval-Procedure)
 
 
 ## Code formatting
@@ -31,9 +32,9 @@ formatting. We provide a script `do-clang-format.sh`, which can be
 very useful. But we urge everybody to integrate `clang-format` already
 on the level of your source code editor. See [the official
 page](https://clang.llvm.org/docs/ClangFormat.html) for information
-about `clang-format` and its editor integration. At least: run
-`do-clang-format.sh` prior to any `git add` command. Code with
-improper formatting will not be accepted for merging. 
+about `clang-format` and its editor integrations. At least: run
+`do-clang-format.sh` prior to any `git add/commit` command. Code with
+improper formatting will not be accepted for merging. It will trigger automatic warnings by the continuous integration (CI) system.
 
 The definition of source code format is written down in the file
 [.clang-format](.clang-format) and can be changed, if the CORSIKA
@@ -42,87 +43,35 @@ e.g. [link1](https://clangformat.com/) or
 [link2](https://zed0.co.uk/clang-format-configurator/).
 
 
-## Naming conventions
-
-While `clang-format` does the structural formatting, we still need to agree on naming conventions:
-
-  - Classes and structs start with capital letters
-  - Class member variables start with "f"
-  - Any static variable has a "g" prefix. A static member variable starts with "fg"
-  - Class member functions start with capital letters
-  - Any class getter begins with "Get", and setter with "Set". Logical getters start with "Is" or "Has".
-  - enums should be "enum class" and start with a capital "E"
-  - Function parameter names start with "v"
-
-  - We use namespaces to avoid clashes and to structure code
-    - *Everything* is part of the corsika namespace
-    - All classes and objects are encapsulated into suited sub-namespaces,
-      thus corsika::geometry, corsika::processes, corsika::units, etc.
-    - Namespace names do not use capital letters.  
-    - Every header file is copied during build and install into
-      "include/corsika/[namespace]" which also means, each header file
-      can only provide definitions for a _single_ namespace. It is one
-      main purpose of namespaces to structure the location of header
-      files.
-    - Each header file uses an include protection that includes at
-      least the namespace name, and header file name, thus, `#ifndef
-      __include_geometry_Point_h__` or `#ifndef __geometry_Point_h__`,
-      or similar are acceptable.
-    - Header files should always be included with `<..>`, thus,
-      `#include <corsika/geometry/Point.h>` since the build system
-      will always provide the correct include directives (and files
-      anyway cannot be found in file-system paths, which generally do
-      not follow the namespace naming conventions outlined
-      here).
-
-  - Header files are named after the main class (or object) they
-    define. This also means each header file name starts with a
-    capital letter.
-  
-
-## Coding rules
-
-  - Code may not introduce any compiler errors, or warnings
-  - All unit tests must succeed at all times. If tests fail, code is not merged.
-  - We use C++17 concepts wherever useful and helpful
-  - On any major error or malfunction we throw an exception. This is needed and required for complex physics and shower debugging.
-  - We never catch exceptions for error handling, there might be very few special exceptions from this. We need to discuss such cases.
-  - Everything that should not change should be `const`
-  - Everything that does not need to be visible to the outside of a class/struct should be `private` or `protected`
-  - We prefer the use of references, wherever useful
-  - There cannot be any pointer in the interface of any class or object
-    exposed to outside users, there might be pointers for very special cases
-    inside of classes.
-  - When you contribute new code, or extend existing code, at the same time provide unit-tests for all functionality.
-  - When you contribute new physics algorithms, in addition you also need to provide a validation module
-  - Code must be documented with `doxygen` commands
-  
+## Coding rules, conventions and guidelines
+
+Please read the [Coding wiki page](https://gitlab.iap.kit.edu/AirShowerPhysics/corsika/-/wikis/Coding-Conventions-and-Guidelines).
+
+
+
 
 ## Release versioning scheme
 
-Releases of CORSIKA are thought to be the baseline for larger scale
-testing, and full production.  The releases are numbered as x.y.z,
-starting with x=8 form the gitlab c++ version. X will only be
-incremented for major design or physics changes. The y index is
-updated for new official releases that normally contain improved or
-enhanced physics performance, and may also contain normal interface
-changes to accomodate improvements. The z index can be updated more
+Releases of CORSIKA 8 are thought to be the baseline for larger scale
+validation, and full production.  The releases are numbered as CORSIKA 8 vx.y.z,
+starting with x=1. The x index is updated for new releases that normally contain improved or
+enhanced physics performance, and also interface
+changes to accomodate improvements. The y and z indices can be updated more
 frequently for bug fixes or new features. Changes in z will not
-contain (major) interface changes, thus, production code will remain
-fully compatible within changes of z.  Special releases of CORSIKA can
-also have a tag name from git, e.g. as in the "milestone1" release.
+contain interface changes, thus, production code will remain
+fully compatible within changes of z.  Special releases of CORSIKA might
+also have a release name.
 
 
 # How to become scientific author of the CORSIKA Project
 
 The CORSIKA Project decides on who becomes scientific author. The
-following conditions are clearly sufficient, but not all of them are
+following conditions are sufficient, but not all of them are
 required all the time:
   - responsibility for a particular functionality or software/management part 
-  - have read and follow these [GUIDELINES](CONTRIBUTING.md)
+  - have read and follow our [contributing guidelines](CONTRIBUTING.md)
   - active in the CORSIKA Project, that means responsive to
-    discussions and problems in corsika-devel@list.kit.edu or on https//gitlab.ikp.kit.edu, 
-    of relevant *Issues*,
-    or in (phone) meetings
-  - agreement to the [COLLABORATION_AGREEMENT](COLLABORATION_AGREEMENT.md) is strictly required
-  - the members of the CORSIKA Project panel agree
+    discussions and problems in corsika-devel@list.kit.edu or on https//gitlab.iap.kit.edu, 
+    of relevant *Issues*, in (phone) meetings or our Mattermost workspace
+  - agreement to the [using and collaborating agreement](USING_COLLABORATING.md) 
+  - the members of the CORSIKA Project must agree

Documentation/CMakeLists.txt

-
-add_subdirectory (Doxygen)
-add_subdirectory (Examples)
-

Documentation/Examples/CMakeLists.txt

-add_executable (helix_example helix_example.cc)
-target_compile_options(helix_example PRIVATE -g) # do not skip asserts
-target_link_libraries(helix_example CORSIKAgeometry CORSIKAunits)
-install (TARGETS helix_example DESTINATION share/examples)
-CORSIKA_ADD_TEST (helix_example)
-
-add_executable (geometry_example geometry_example.cc)
-target_compile_options(geometry_example PRIVATE -g) # do not skip asserts
-target_link_libraries (geometry_example CORSIKAgeometry CORSIKAunits)
-install (TARGETS geometry_example DESTINATION share/examples)
-CORSIKA_ADD_TEST (geometry_example)
-
-add_executable (logger_example logger_example.cc)
-target_compile_options(logger_example PRIVATE -g) # do not skip asserts
-target_link_libraries (logger_example CORSIKAunits CORSIKAlogging)
-install (TARGETS logger_example DESTINATION share/examples)
-CORSIKA_ADD_TEST (logger_example)
-
-add_executable (stack_example stack_example.cc)
-target_compile_options(stack_example PRIVATE -g) # do not skip asserts
-target_link_libraries (stack_example SuperStupidStack CORSIKAunits
-  CORSIKAlogging)
-CORSIKA_ADD_TEST (stack_example)
-
-add_executable (cascade_example cascade_example.cc)
-target_compile_options(cascade_example PRIVATE -g) # do not skip asserts
-target_link_libraries (cascade_example SuperStupidStack CORSIKAunits
-  CORSIKAlogging
-  CORSIKArandom
-  ProcessSibyll
-  ProcessPythia
-  CORSIKAcascade
-  ProcessEnergyLoss
-  ProcessStackInspector
-  ProcessTrackWriter
-  ProcessTrackingLine
-  ProcessHadronicElasticModel
-  CORSIKAprocesses
-  CORSIKAcascade
-  CORSIKAparticles
-  CORSIKAgeometry
-  CORSIKAenvironment
-  CORSIKAprocesssequence
-  )
-install (TARGETS cascade_example DESTINATION share/examples)
-CORSIKA_ADD_TEST (cascade_example)
-
-add_executable (staticsequence_example staticsequence_example.cc)
-target_compile_options(staticsequence_example PRIVATE -g) # do not skip asserts
-target_link_libraries (staticsequence_example
-  CORSIKAprocesssequence
-  CORSIKAunits
-  CORSIKAgeometry
-  CORSIKAlogging)
-install (TARGETS staticsequence_example DESTINATION share/examples)
-CORSIKA_ADD_TEST (staticsequence_example)

Documentation/Examples/cascade_example.cc

-
-/*
- * (c) Copyright 2018 CORSIKA Project, corsika-project@lists.kit.edu
- *
- * See file AUTHORS for a list of contributors.
- *
- * 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/cascade/Cascade.h>
-#include <corsika/process/ProcessSequence.h>
-#include <corsika/process/energy_loss/EnergyLoss.h>
-#include <corsika/process/hadronic_elastic_model/HadronicElasticModel.h>
-#include <corsika/process/stack_inspector/StackInspector.h>
-#include <corsika/process/tracking_line/TrackingLine.h>
-
-#include <corsika/setup/SetupStack.h>
-#include <corsika/setup/SetupTrajectory.h>
-
-#include <corsika/environment/Environment.h>
-#include <corsika/environment/HomogeneousMedium.h>
-#include <corsika/environment/NuclearComposition.h>
-
-#include <corsika/geometry/Sphere.h>
-
-#include <corsika/process/sibyll/Decay.h>
-#include <corsika/process/sibyll/Interaction.h>
-#include <corsika/process/sibyll/NuclearInteraction.h>
-
-#include <corsika/process/pythia/Decay.h>
-
-#include <corsika/process/track_writer/TrackWriter.h>
-
-#include <corsika/units/PhysicalUnits.h>
-
-#include <corsika/random/RNGManager.h>
-
-#include <corsika/utl/CorsikaFenv.h>
-
-#include <boost/type_index.hpp>
-using boost::typeindex::type_id_with_cvr;
-
-#include <iostream>
-#include <limits>
-#include <typeinfo>
-
-using namespace corsika;
-using namespace corsika::process;
-using namespace corsika::units;
-using namespace corsika::particles;
-using namespace corsika::random;
-using namespace corsika::setup;
-using namespace corsika::geometry;
-using namespace corsika::environment;
-
-using namespace std;
-using namespace corsika::units::si;
-
-class ProcessCut : public process::ContinuousProcess<ProcessCut> {
-
-  HEPEnergyType fECut;
-
-  HEPEnergyType fEnergy = 0_GeV;
-  HEPEnergyType fEmEnergy = 0_GeV;
-  int fEmCount = 0;
-  HEPEnergyType fInvEnergy = 0_GeV;
-  int fInvCount = 0;
-
-public:
-  ProcessCut(const HEPEnergyType cut)
-      : fECut(cut) {}
-
-  template <typename Particle>
-  bool isBelowEnergyCut(Particle& p) const {
-    // nuclei
-    if (p.GetPID() == particles::Code::Nucleus) {
-      auto const ElabNuc = p.GetEnergy() / p.GetNuclearA();
-      auto const EcmNN = sqrt(2. * ElabNuc * 0.93827_GeV);
-      if (ElabNuc < fECut || EcmNN < 10_GeV)
-        return true;
-      else
-        return false;
-    } else {
-      // TODO: center-of-mass energy hard coded
-      const HEPEnergyType Ecm = sqrt(2. * p.GetEnergy() * 0.93827_GeV);
-      if (p.GetEnergy() < fECut || Ecm < 10_GeV)
-        return true;
-      else
-        return false;
-    }
-  }
-
-  bool isEmParticle(Code pCode) const {
-    bool is_em = false;
-    // FOR NOW: switch
-    switch (pCode) {
-      case Code::Electron:
-        is_em = true;
-        break;
-      case Code::Positron:
-        is_em = true;
-        break;
-      case Code::Gamma:
-        is_em = true;
-        break;
-      default:
-        break;
-    }
-    return is_em;
-  }
-
-  void defineEmParticles() const {
-    // create bool array identifying em particles
-  }
-
-  bool isInvisible(Code pCode) const {
-    bool is_inv = false;
-    // FOR NOW: switch
-    switch (pCode) {
-      case Code::NuE:
-        is_inv = true;
-        break;
-      case Code::NuEBar:
-        is_inv = true;
-        break;
-      case Code::NuMu:
-        is_inv = true;
-        break;
-      case Code::NuMuBar:
-        is_inv = true;
-        break;
-      case Code::MuPlus:
-        is_inv = true;
-        break;
-      case Code::MuMinus:
-        is_inv = true;
-        break;
-
-      case Code::Neutron:
-        is_inv = true;
-        break;
-
-      case Code::AntiNeutron:
-        is_inv = true;
-        break;
-
-      default:
-        break;
-    }
-    return is_inv;
-  }
-
-  template <typename Particle>
-  LengthType MaxStepLength(Particle& p, setup::Trajectory&) const {
-    cout << "ProcessCut: MinStep: pid: " << p.GetPID() << endl;
-    cout << "ProcessCut: MinStep: energy (GeV): " << p.GetEnergy() / 1_GeV << endl;
-    const Code pid = p.GetPID();
-    if (isEmParticle(pid) || isInvisible(pid) || isBelowEnergyCut(p)) {
-      cout << "ProcessCut: MinStep: next cut: " << 0. << endl;
-      return 0_m;
-    } else {
-      LengthType next_step = 1_m * std::numeric_limits<double>::infinity();
-      cout << "ProcessCut: MinStep: next cut: " << next_step << endl;
-      return next_step;
-    }
-  }
-
-  template <typename Particle, typename Stack>
-  EProcessReturn DoContinuous(Particle& p, setup::Trajectory&, Stack&) {
-    const Code pid = p.GetPID();
-    HEPEnergyType energy = p.GetEnergy();
-    cout << "ProcessCut: DoContinuous: " << pid << " E= " << energy
-         << ", EcutTot=" << (fEmEnergy + fInvEnergy + fEnergy) / 1_GeV << " GeV" << endl;
-    EProcessReturn ret = EProcessReturn::eOk;
-    if (isEmParticle(pid)) {
-      cout << "removing em. particle..." << endl;
-      fEmEnergy += energy;
-      fEmCount += 1;
-      // p.Delete();
-      ret = EProcessReturn::eParticleAbsorbed;
-    } else if (isInvisible(pid)) {
-      cout << "removing inv. particle..." << endl;
-      fInvEnergy += energy;
-      fInvCount += 1;
-      // p.Delete();
-      ret = EProcessReturn::eParticleAbsorbed;
-    } else if (isBelowEnergyCut(p)) {
-      cout << "removing low en. particle..." << endl;
-      fEnergy += energy;
-      // p.Delete();
-      ret = EProcessReturn::eParticleAbsorbed;
-    }
-    return ret;
-  }
-
-  void Init() {
-    fEmEnergy = 0. * 1_GeV;
-    fEmCount = 0;
-    fInvEnergy = 0. * 1_GeV;
-    fInvCount = 0;
-    fEnergy = 0. * 1_GeV;
-    // defineEmParticles();
-  }
-
-  void ShowResults() {
-    cout << " ******************************" << endl
-         << " ParticleCut: " << endl
-         << " energy in em.  component (GeV):  " << fEmEnergy / 1_GeV << endl
-         << " no. of em.  particles injected:  " << fEmCount << endl
-         << " energy in inv. component (GeV):  " << fInvEnergy / 1_GeV << endl
-         << " no. of inv. particles injected:  " << fInvCount << endl
-         << " energy below particle cut (GeV): " << fEnergy / 1_GeV << endl
-         << " ******************************" << endl;
-  }
-
-  HEPEnergyType GetInvEnergy() const { return fInvEnergy; }
-  HEPEnergyType GetCutEnergy() const { return fEnergy; }
-  HEPEnergyType GetEmEnergy() const { return fEmEnergy; }
-};
-
-//
-// The example main program for a particle cascade
-//
-int main() {
-  feenableexcept(FE_INVALID);
-  // initialize random number sequence(s)
-  random::RNGManager::GetInstance().RegisterRandomStream("cascade");
-
-  // setup environment, geometry
-  environment::Environment env;
-  auto& universe = *(env.GetUniverse());
-
-  auto theMedium = environment::Environment::CreateNode<Sphere>(
-      Point{env.GetCoordinateSystem(), 0_m, 0_m, 0_m},
-      1_km * std::numeric_limits<double>::infinity());
-
-  // fraction of oxygen
-  const float fox = 0.20946;
-  using MyHomogeneousModel = environment::HomogeneousMedium<environment::IMediumModel>;
-  theMedium->SetModelProperties<MyHomogeneousModel>(
-      1_kg / (1_m * 1_m * 1_m),
-      environment::NuclearComposition(
-          std::vector<particles::Code>{particles::Code::Nitrogen,
-                                       particles::Code::Oxygen},
-          std::vector<float>{(float)1. - fox, fox}));
-
-  universe.AddChild(std::move(theMedium));
-
-  const CoordinateSystem& rootCS = env.GetCoordinateSystem();
-
-  // setup processes, decays and interactions
-  tracking_line::TrackingLine<setup::Stack, setup::Trajectory> tracking(env);
-  stack_inspector::StackInspector<setup::Stack> stackInspect(true);
-
-  const std::vector<particles::Code> trackedHadrons = {
-      particles::Code::PiPlus, particles::Code::PiMinus, particles::Code::KPlus,
-      particles::Code::KMinus, particles::Code::K0Long,  particles::Code::K0Short};
-
-  random::RNGManager::GetInstance().RegisterRandomStream("s_rndm");
-  process::sibyll::Interaction sibyll(env);
-  process::sibyll::NuclearInteraction sibyllNuc(env, sibyll);
-  process::sibyll::Decay decay(trackedHadrons);
-  // random::RNGManager::GetInstance().RegisterRandomStream("pythia");
-  // process::pythia::Decay decay(trackedHadrons);
-  ProcessCut cut(20_GeV);
-
-  // random::RNGManager::GetInstance().RegisterRandomStream("HadronicElasticModel");
-  // process::HadronicElasticModel::HadronicElasticInteraction
-  // hadronicElastic(env);
-
-  process::TrackWriter::TrackWriter trackWriter("tracks.dat");
-  process::EnergyLoss::EnergyLoss eLoss;
-
-  // assemble all processes into an ordered process list
-  // auto sequence = stackInspect << sibyll << decay << hadronicElastic << cut <<
-  // trackWriter; auto sequence = stackInspect << sibyll << sibyllNuc << decay << eLoss <<
-  // cut << trackWriter;
-  // auto sequence = sibyll << sibyllNuc << decay << eLoss << cut;
-  auto sequence = stackInspect << sibyll << sibyllNuc << decay << eLoss << cut;
-
-  // cout << "decltype(sequence)=" << type_id_with_cvr<decltype(sequence)>().pretty_name()
-  // << "\n";
-
-  // setup particle stack, and add primary particle
-  setup::Stack stack;
-  stack.Clear();
-  const Code beamCode = Code::Nucleus;
-  const int nuclA = 4;
-  const int nuclZ = int(nuclA / 2.15 + 0.7);
-  const HEPMassType mass = GetNucleusMass(nuclA, nuclZ);
-  const HEPEnergyType E0 = nuclA * 10_TeV;
-  double theta = 0.;
-  double phi = 0.;
-
-  {
-    auto elab2plab = [](HEPEnergyType Elab, HEPMassType m) {
-      return sqrt((Elab - m) * (Elab + m));
-    };
-    HEPMomentumType P0 = elab2plab(E0, mass);
-    auto momentumComponents = [](double theta, double phi, HEPMomentumType ptot) {
-      return std::make_tuple(ptot * sin(theta) * cos(phi), ptot * sin(theta) * sin(phi),
-                             -ptot * cos(theta));
-    };
-    auto const [px, py, pz] =
-        momentumComponents(theta / 180. * M_PI, phi / 180. * M_PI, P0);
-    auto plab = corsika::stack::MomentumVector(rootCS, {px, py, pz});
-    cout << "input particle: " << beamCode << endl;
-    cout << "input angles: theta=" << theta << " phi=" << phi << endl;
-    cout << "input momentum: " << plab.GetComponents() / 1_GeV << endl;
-    Point pos(rootCS, 0_m, 0_m, 0_m);
-    stack.AddParticle(std::tuple<particles::Code, units::si::HEPEnergyType,
-                                 corsika::stack::MomentumVector, geometry::Point,
-                                 units::si::TimeType, unsigned short, unsigned short>{
-        beamCode, E0, plab, pos, 0_ns, nuclA, nuclZ});
-  }
-
-  // define air shower object, run simulation
-  cascade::Cascade EAS(env, tracking, sequence, stack);
-  EAS.Init();
-  EAS.Run();
-
-  eLoss.PrintProfile(); // print longitudinal profile
-
-  cut.ShowResults();
-  const HEPEnergyType Efinal =
-      cut.GetCutEnergy() + cut.GetInvEnergy() + cut.GetEmEnergy();
-  cout << "total cut energy (GeV): " << Efinal / 1_GeV << endl
-       << "relative difference (%): " << (Efinal / E0 - 1) * 100 << endl;
-  cout << "total dEdX energy (GeV): " << eLoss.GetTotal() / 1_GeV << endl
-       << "relative difference (%): " << eLoss.GetTotal() / E0 * 100 << endl;
-}

Documentation/Examples/logger_example.cc

-
-/*
- * (c) Copyright 2018 CORSIKA Project, corsika-project@lists.kit.edu
- *
- * See file AUTHORS for a list of contributors.
- *
- * 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/logging/Logger.h>
-#include <boost/format.hpp>
-#include <fstream>
-#include <iostream>
-#include <string>
-
-using namespace std;
-using namespace corsika::logging;
-
-int main() {
-  {
-    cout << "writing to \"another.log\"" << endl;
-    ofstream logfile("another.log");
-    sink::SinkStream unbuffered_sink(logfile);
-    sink::BufferedSinkStream sink(logfile, sink::StdBuffer(10000));
-    Logger<MessageOn, sink::BufferedSinkStream> info("\033[32m", "info", sink);
-    Logger<MessageOn, sink::BufferedSinkStream> err("\033[31m", "error", sink);
-    // logger<ostream,messageconst,StdBuffer> info(std::cout, StdBuffer(10000));
-
-    /*
-      Logging& logs = Logging::GetInstance();
-      logs.AddLogger<>("info", info);
-      auto& log_1 = logs.GetLogger("info"); // no so useful, since type of log_1 is
-      std::any
-    */
-
-    for (int i = 0; i < 10000; ++i) {
-      LOG(info, "irgendwas", " ", string("and more"), " ",
-          boost::format("error: %i message: %s. done."), i, "stupido");
-      LOG(err, "Fehler");
-    }
-  }
-
-  {
-    sink::NoSink off;
-    Logger<MessageOff> info("", "", off);
-
-    for (int i = 0; i < 10000; ++i) {
-      LOG(info, "irgendwas", string("and more"),
-          boost::format("error: %i message: %s. done."), i, "stupido", "a-number:", 8.99,
-          "ENDE");
-    }
-  }
-
-  return 0;
-}

Documentation/Examples/stack_example.cc

-
-/*
- * (c) Copyright 2018 CORSIKA Project, corsika-project@lists.kit.edu
- *
- * See file AUTHORS for a list of contributors.
- *
- * 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/particles/ParticleProperties.h>
-#include <corsika/stack/super_stupid/SuperStupidStack.h>
-
-#include <corsika/geometry/Point.h>
-#include <corsika/geometry/RootCoordinateSystem.h>
-
-#include <cassert>
-#include <iomanip>
-#include <iostream>
-
-using namespace corsika;
-using namespace corsika::units::si;
-using namespace corsika::stack;
-using namespace corsika::geometry;
-using namespace std;
-
-void fill(corsika::stack::super_stupid::SuperStupidStack& s) {
-  const geometry::CoordinateSystem& rootCS =
-      geometry::RootCoordinateSystem::GetInstance().GetRootCoordinateSystem();
-  for (int i = 0; i < 11; ++i) {
-    s.AddParticle(
-        std::tuple<particles::Code, units::si::HEPEnergyType,
-                   corsika::stack::MomentumVector, geometry::Point, units::si::TimeType>{
-            particles::Code::Electron, 1.5_GeV * i,
-            corsika::stack::MomentumVector(rootCS, {0_GeV, 0_GeV, 1_GeV}),
-            geometry::Point(rootCS, 0_m, 0_m, 0_m), 0_ns});
-  }
-}
-
-void read(corsika::stack::super_stupid::SuperStupidStack& s) {
-  assert(s.GetSize() == 11); // stack has 11 particles
-
-  HEPEnergyType total_energy;
-  int i = 0;
-  for (auto& p : s) {
-    total_energy += p.GetEnergy();
-    // particles are electrons with 1.5 GeV energy times i
-    assert(p.GetPID() == particles::Code::Electron);
-    assert(p.GetEnergy() == 1.5_GeV * (i++));
-  }
-}
-
-int main() {
-  corsika::stack::super_stupid::SuperStupidStack s;
-  fill(s);
-  read(s);
-  return 0;
-}

Documentation/Examples/staticsequence_example.cc

-
-/*
- * (c) Copyright 2018 CORSIKA Project, corsika-project@lists.kit.edu
- *
- * See file AUTHORS for a list of contributors.
- *
- * 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 <array>
-#include <iomanip>
-#include <iostream>
-
-#include <corsika/process/ProcessSequence.h>
-
-#include <corsika/geometry/Point.h>
-#include <corsika/geometry/RootCoordinateSystem.h>
-#include <corsika/geometry/Vector.h>
-
-using namespace corsika;
-using namespace corsika::units::si;
-using namespace corsika::process;
-using namespace std;
-
-const int nData = 10;
-
-class Process1 : public BaseProcess<Process1> {
-public:
-  Process1() {}
-  template <typename D, typename T, typename S>
-  EProcessReturn DoContinuous(D& d, T&, S&) const {
-    for (int i = 0; i < nData; ++i) d.p[i] += 1;
-    return EProcessReturn::eOk;
-  }
-};
-
-class Process2 : public BaseProcess<Process2> {
-public:
-  Process2() {}
-
-  template <typename D, typename T, typename S>
-  inline EProcessReturn DoContinuous(D& d, T&, S&) const {
-    for (int i = 0; i < nData; ++i) d.p[i] -= 0.1 * i;
-    return EProcessReturn::eOk;
-  }
-};
-
-class Process3 : public BaseProcess<Process3> {
-public:
-  Process3() {}
-
-  template <typename D, typename T, typename S>
-  inline EProcessReturn DoContinuous(D&, T&, S&) const {
-    return EProcessReturn::eOk;
-  }
-};
-
-class Process4 : public BaseProcess<Process4> {
-public:
-  Process4(const double v)
-      : fV(v) {}
-  template <typename D, typename T, typename S>
-  inline EProcessReturn DoContinuous(D& d, T&, S&) const {
-    for (int i = 0; i < nData; ++i) d.p[i] *= fV;
-    return EProcessReturn::eOk;
-  }
-
-private:
-  double fV;
-};
-
-struct DummyData {
-  double p[nData] = {0, 0, 0, 0, 0, 0, 0, 0, 0, 0};
-};
-struct DummyStack {
-  void clear() {}
-};
-struct DummyTrajectory {};
-
-void modular() {
-
-  Process1 m1;
-  Process2 m2;
-  Process3 m3;
-  Process4 m4(0.9);
-
-  auto sequence = m1 << m2 << m3 << m4;
-
-  DummyData p;
-  DummyStack s;
-  DummyTrajectory t;
-
-  const int n = 1000;
-  for (int i = 0; i < n; ++i) { sequence.DoContinuous(p, t, s); }
-
-  for (int i = 0; i < nData; ++i) {
-    // cout << p.p[i] << endl;
-    // assert(p.p[i] == n-i*100);
-  }
-
-  cout << " done (nothing...) " << endl;
-}
-
-int main() {
-  modular();
-  return 0;
-}

Environment/CMakeLists.txt

-set (
-  ENVIRONMENT_HEADERS
-  VolumeTreeNode.h
-  IMediumModel.h
-  NuclearComposition.h
-  HomogeneousMedium.h
-  InhomogeneousMedium.h
-  HomogeneousMedium.h
-  LinearApproximationIntegrator.h
-  DensityFunction.h
-  Environment.h
-  FlatExponential.h
-  )
-
-set (
-  ENVIRONMENT_NAMESPACE
-  corsika/environment
-  )
-
-add_library (CORSIKAenvironment INTERFACE)
-CORSIKA_COPY_HEADERS_TO_NAMESPACE (CORSIKAenvironment ${ENVIRONMENT_NAMESPACE} ${ENVIRONMENT_HEADERS})
-
-# target dependencies on other libraries (also the header onlys)
-target_link_libraries (
-  CORSIKAenvironment
-  INTERFACE
-  CORSIKAgeometry
-  CORSIKAparticles
-  CORSIKAunits
-  CORSIKArandom
-  )
-
-target_include_directories (
-  CORSIKAenvironment
-  INTERFACE
-  $<BUILD_INTERFACE:${PROJECT_BINARY_DIR}/include>
-  $<INSTALL_INTERFACE:include>
-  )
-
-install (
-  TARGETS CORSIKAenvironment
-  LIBRARY DESTINATION lib
-  ARCHIVE DESTINATION lib
-  PUBLIC_HEADER DESTINATION include/${ENVIRONMENT_NAMESPACE}
-  )
-
-# --------------------
-# code unit testing
-add_executable (testEnvironment testEnvironment.cc)
-
-target_link_libraries (
-  testEnvironment
-  CORSIKAsetup
-  CORSIKAenvironment
-  CORSIKAthirdparty # for catch2
-  )
-CORSIKA_ADD_TEST(testGeometry)

Environment/Environment.h

-
-/*
- * (c) Copyright 2018 CORSIKA Project, corsika-project@lists.kit.edu
- *
- * See file AUTHORS for a list of contributors.
- *
- * 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.
- */
-
-#ifndef _include_Environment_h
-#define _include_Environment_h
-
-#include <corsika/environment/IMediumModel.h>
-#include <corsika/environment/VolumeTreeNode.h>
-#include <corsika/geometry/Point.h>
-#include <corsika/geometry/RootCoordinateSystem.h>
-#include <corsika/geometry/Sphere.h>
-#include <corsika/setup/SetupEnvironment.h>
-#include <limits>
-
-namespace corsika::environment {
-  using BaseNodeType = VolumeTreeNode<corsika::setup::IEnvironmentModel>;
-
-  struct Universe : public corsika::geometry::Sphere {
-    Universe(corsika::geometry::CoordinateSystem const& pCS)
-        : corsika::geometry::Sphere(
-              corsika::geometry::Point{pCS, 0 * corsika::units::si::meter,
-                                       0 * corsika::units::si::meter,
-                                       0 * corsika::units::si::meter},
-              corsika::units::si::meter * std::numeric_limits<double>::infinity()) {}
-
-    bool Contains(corsika::geometry::Point const&) const override { return true; }
-  };
-
-  // template <typename IEnvironmentModel>
-  class Environment {
-  public:
-    Environment()
-        : fCoordinateSystem{corsika::geometry::RootCoordinateSystem::GetInstance()
-                                .GetRootCoordinateSystem()}
-        , fUniverse(std::make_unique<BaseNodeType>(
-              std::make_unique<Universe>(fCoordinateSystem))) {}
-
-    using IEnvironmentModel = corsika::setup::IEnvironmentModel;
-
-    auto& GetUniverse() { return fUniverse; }
-    auto const& GetUniverse() const { return fUniverse; }
-
-    auto const& GetCoordinateSystem() const { return fCoordinateSystem; }
-
-    // factory method for creation of VolumeTreeNodes
-    template <class TVolumeType, typename... TVolumeArgs>
-    static auto CreateNode(TVolumeArgs&&... args) {
-      static_assert(std::is_base_of_v<corsika::geometry::Volume, TVolumeType>,
-                    "unusable type provided, needs to be derived from "
-                    "\"corsika::geometry::Volume\"");
-
-      return std::make_unique<BaseNodeType>(
-          std::make_unique<TVolumeType>(std::forward<TVolumeArgs>(args)...));
-    }
-
-  private:
-    corsika::geometry::CoordinateSystem const& fCoordinateSystem;
-    BaseNodeType::VTNUPtr fUniverse;
-  };
-
-} // namespace corsika::environment
-
-#endif

Environment/FlatAtmosphere/FlatAtmosphere.h

-
-/*
- * (c) Copyright 2018 CORSIKA Project, corsika-project@lists.kit.edu
- *
- * See file AUTHORS for a list of contributors.
- *
- * 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.
- */

Environment/FlatExponential.h

-/*
- * (c) Copyright 2019 CORSIKA Project, corsika-project@lists.kit.edu
- *
- * See file AUTHORS for a list of contributors.
- *
- * 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.
- */
-
-#ifndef _include_Environment_FlatExponential_h_
-#define _include_Environment_FlatExponential_h_
-
-#include <corsika/environment/NuclearComposition.h>
-#include <corsika/geometry/Line.h>
-#include <corsika/geometry/Point.h>
-#include <corsika/geometry/Trajectory.h>
-#include <corsika/particles/ParticleProperties.h>
-#include <corsika/random/RNGManager.h>
-#include <corsika/units/PhysicalUnits.h>
-
-#include <cassert>
-#include <limits>
-
-/**
- *
- */
-
-namespace corsika::environment {
-
-  template <class T>
-  class FlatExponential : public T {
-    corsika::units::si::MassDensityType const fRho0;
-    units::si::LengthType const fLambda;
-    units::si::InverseLengthType const fInvLambda;
-    NuclearComposition const fNuclComp;
-    geometry::Vector<units::si::dimensionless_d> const fAxis;
-    geometry::Point const fP0;
-
-  public:
-    FlatExponential(geometry::Point const& p0,
-                    geometry::Vector<units::si::dimensionless_d> const& axis,
-                    units::si::MassDensityType rho, units::si::LengthType lambda,
-                    NuclearComposition pNuclComp)
-        : fRho0(rho)
-        , fLambda(lambda)
-        , fInvLambda(1 / lambda)
-        , fNuclComp(pNuclComp)
-        , fAxis(axis)
-        , fP0(p0) {}
-
-    corsika::units::si::MassDensityType GetMassDensity(
-        corsika::geometry::Point const& p) const override {
-      return fRho0 * exp(fInvLambda * (p - fP0).dot(fAxis));
-    }
-    NuclearComposition const& GetNuclearComposition() const override { return fNuclComp; }
-
-    corsika::units::si::GrammageType IntegratedGrammage(
-        corsika::geometry::Trajectory<corsika::geometry::Line> const& line,
-        corsika::units::si::LengthType pTo) const override {
-      auto const vDotA = line.NormalizedDirection().dot(fAxis).magnitude();
-
-      if (vDotA == 0) {
-        return pTo * GetMassDensity(line.GetR0());
-      } else {
-        return GetMassDensity(line.GetR0()) * (fLambda / vDotA) *
-               (exp(vDotA * pTo / fLambda) - 1);
-      }
-    }
-
-    corsika::units::si::LengthType ArclengthFromGrammage(
-        corsika::geometry::Trajectory<corsika::geometry::Line> const& line,
-        corsika::units::si::GrammageType pGrammage) const override {
-      auto const vDotA = line.NormalizedDirection().dot(fAxis).magnitude();
-
-      if (vDotA == 0) {
-        return pGrammage / GetMassDensity(line.GetR0());
-      } else {
-        auto const logArg = pGrammage * vDotA / (fRho0 * fLambda) + 1;
-        if (logArg > 0) {
-          return fLambda / vDotA * log(pGrammage * vDotA / (fRho0 * fLambda) + 1);
-        } else {
-          return std::numeric_limits<typename decltype(
-                     pGrammage)::value_type>::infinity() *
-                 corsika::units::si::meter;
-        }
-      }
-    }
-  };
-
-} // namespace corsika::environment
-#endif

Environment/HomogeneousMedium.h

-
-/*
- * (c) Copyright 2018 CORSIKA Project, corsika-project@lists.kit.edu
- *
- * See file AUTHORS for a list of contributors.
- *
- * 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.
- */
-
-#ifndef _include_HomogeneousMedium_h_
-#define _include_HomogeneousMedium_h_
-
-#include <corsika/environment/NuclearComposition.h>
-#include <corsika/geometry/Line.h>
-#include <corsika/geometry/Point.h>
-#include <corsika/geometry/Trajectory.h>
-#include <corsika/particles/ParticleProperties.h>
-#include <corsika/random/RNGManager.h>
-#include <corsika/units/PhysicalUnits.h>
-
-#include <cassert>
-
-/**
- * a homogeneous medium
- */
-
-namespace corsika::environment {
-
-  template <class T>
-  class HomogeneousMedium : public T {
-    corsika::units::si::MassDensityType const fDensity;
-    NuclearComposition const fNuclComp;
-
-  public:
-    HomogeneousMedium(corsika::units::si::MassDensityType pDensity,
-                      NuclearComposition pNuclComp)
-        : fDensity(pDensity)
-        , fNuclComp(pNuclComp) {}
-
-    corsika::units::si::MassDensityType GetMassDensity(
-        corsika::geometry::Point const&) const override {
-      return fDensity;
-    }
-    NuclearComposition const& GetNuclearComposition() const override { return fNuclComp; }
-
-    corsika::units::si::GrammageType IntegratedGrammage(
-        corsika::geometry::Trajectory<corsika::geometry::Line> const&,
-        corsika::units::si::LengthType pTo) const override {
-      using namespace corsika::units::si;
-      return pTo * fDensity;
-    }
-
-    corsika::units::si::LengthType ArclengthFromGrammage(
-        corsika::geometry::Trajectory<corsika::geometry::Line> const&,
-        corsika::units::si::GrammageType pGrammage) const override {
-      return pGrammage / fDensity;
-    }
-  };
-
-} // namespace corsika::environment
-#endif

Environment/IMediumModel.h

-
-/*
- * (c) Copyright 2018 CORSIKA Project, corsika-project@lists.kit.edu
- *
- * See file AUTHORS for a list of contributors.
- *
- * 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.
- */
-
-#ifndef _include_IMediumModel_h
-#define _include_IMediumModel_h
-
-#include <corsika/environment/NuclearComposition.h>
-#include <corsika/geometry/Line.h>
-#include <corsika/geometry/Point.h>
-#include <corsika/geometry/Trajectory.h>
-#include <corsika/units/PhysicalUnits.h>
-
-namespace corsika::environment {
-
-  class IMediumModel {
-  public:
-    virtual ~IMediumModel() = default;
-
-    virtual corsika::units::si::MassDensityType GetMassDensity(
-        corsika::geometry::Point const&) const = 0;
-
-    // todo: think about the mixin inheritance of the trajectory vs the BaseTrajectory
-    // approach for now, only lines are supported
-    virtual corsika::units::si::GrammageType IntegratedGrammage(
-        corsika::geometry::Trajectory<corsika::geometry::Line> const&,
-        corsika::units::si::LengthType) const = 0;
-
-    virtual corsika::units::si::LengthType ArclengthFromGrammage(
-        corsika::geometry::Trajectory<corsika::geometry::Line> const&,
-        corsika::units::si::GrammageType) const = 0;
-
-    virtual NuclearComposition const& GetNuclearComposition() const = 0;
-  };
-
-} // namespace corsika::environment
-
-#endif

Environment/InhomogeneousMedium.h

-
-/*
- * (c) Copyright 2018 CORSIKA Project, corsika-project@lists.kit.edu
- *
- * See file AUTHORS for a list of contributors.
- *
- * 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.
- */
-
-#ifndef _include_environment_InhomogeneousMedium_h_
-#define _include_environment_InhomogeneousMedium_h_
-
-#include <corsika/environment/NuclearComposition.h>
-#include <corsika/geometry/Line.h>
-#include <corsika/geometry/Point.h>
-#include <corsika/geometry/Trajectory.h>
-#include <corsika/particles/ParticleProperties.h>
-#include <corsika/random/RNGManager.h>
-#include <corsika/units/PhysicalUnits.h>
-
-/**
- * A general inhomogeneous medium. The mass density distribution TDensityFunction must be
- * a \f$C^2\f$-function.
- */
-
-namespace corsika::environment {
-
-  template <class T, class TDensityFunction>
-  class InhomogeneousMedium : public T {
-    NuclearComposition const fNuclComp;
-    TDensityFunction const fDensityFunction;
-
-  public:
-    template <typename... Args>
-    InhomogeneousMedium(NuclearComposition pNuclComp, Args&&... rhoArgs)
-        : fNuclComp(pNuclComp)
-        , fDensityFunction(rhoArgs...) {}
-
-    corsika::units::si::MassDensityType GetMassDensity(
-        corsika::geometry::Point const& p) const override {
-      return fDensityFunction.EvaluateAt(p);
-    }
-    NuclearComposition const& GetNuclearComposition() const override { return fNuclComp; }
-
-    corsika::units::si::GrammageType IntegratedGrammage(
-        corsika::geometry::Trajectory<corsika::geometry::Line> const& pLine,
-        corsika::units::si::LengthType pTo) const override {
-      return fDensityFunction.IntegrateGrammage(pLine, pTo);
-    }
-
-    corsika::units::si::LengthType ArclengthFromGrammage(
-        corsika::geometry::Trajectory<corsika::geometry::Line> const& pLine,
-        corsika::units::si::GrammageType pGrammage) const override {
-      return fDensityFunction.ArclengthFromGrammage(pLine, pGrammage);
-    }
-  };
-
-} // namespace corsika::environment
-#endif

Environment/LinearApproximationIntegrator.h

-
-/*
- * (c) Copyright 2018 CORSIKA Project, corsika-project@lists.kit.edu
- *
- * See file AUTHORS for a list of contributors.
- *
- * 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.
- */
-
-#ifndef _include_environment_LinearApproximationIntegrator_h_
-#define _include_environment_LinearApproximationIntegrator_h_
-
-#include <limits>
-
-#include <corsika/geometry/Line.h>
-#include <corsika/geometry/Point.h>
-#include <corsika/geometry/Trajectory.h>
-
-namespace corsika::environment {
-  template <class TDerived>
-  class LinearApproximationIntegrator {
-    auto const& GetImplementation() const { return *static_cast<TDerived const*>(this); }
-
-  public:
-    auto IntegrateGrammage(
-        corsika::geometry::Trajectory<corsika::geometry::Line> const& line,
-        corsika::units::si::LengthType length) const {
-      auto const c0 = GetImplementation().EvaluateAt(line.GetPosition(0));
-      auto const c1 = GetImplementation().fRho.FirstDerivative(
-          line.GetPosition(0), line.NormalizedDirection());
-      return (c0 + 0.5 * c1 * length) * length;
-    }
-
-    auto ArclengthFromGrammage(
-        corsika::geometry::Trajectory<corsika::geometry::Line> const& line,
-        corsika::units::si::GrammageType grammage) const {
-      auto const c0 = GetImplementation().fRho(line.GetPosition(0));
-      auto const c1 = GetImplementation().fRho.FirstDerivative(
-          line.GetPosition(0), line.NormalizedDirection());
-
-      return (1 - 0.5 * grammage * c1 / (c0 * c0)) * grammage / c0;
-    }
-
-    auto MaximumLength(corsika::geometry::Trajectory<corsika::geometry::Line> const& line,
-                       [[maybe_unused]] double relError) const {
-      using namespace corsika::units::si;
-      [[maybe_unused]] auto const c1 = GetImplementation().fRho.SecondDerivative(
-          line.GetPosition(0), line.NormalizedDirection());
-
-      // todo: provide a real, working implementation
-      return 1_m * std::numeric_limits<double>::infinity();
-    }
-  };
-} // namespace corsika::environment
-
-#endif

Environment/NuclearComposition.h

-
-/*
- * (c) Copyright 2018 CORSIKA Project, corsika-project@lists.kit.edu
- *
- * See file AUTHORS for a list of contributors.
- *
- * 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.
- */
-
-#ifndef _include_NuclearComposition_h
-#define _include_NuclearComposition_h
-
-#include <corsika/particles/ParticleProperties.h>
-#include <cassert>
-#include <numeric>
-#include <random>
-#include <stdexcept>
-#include <vector>
-
-namespace corsika::environment {
-  class NuclearComposition {
-    std::vector<float> const fNumberFractions; //!< relative fractions of number density
-    std::vector<corsika::particles::Code> const
-        fComponents; //!< particle codes of consitutents
-
-    double const fAvgMassNumber;
-
-    template <class AConstIterator, class BConstIterator>
-    class WeightProviderIterator {
-      AConstIterator fAIter;
-      BConstIterator fBIter;
-
-    public:
-      using value_type = double;
-      using iterator_category = std::input_iterator_tag;
-      using pointer = value_type*;
-      using reference = value_type&;
-      using difference_type = ptrdiff_t;
-
-      WeightProviderIterator(AConstIterator a, BConstIterator b)
-          : fAIter(a)
-          , fBIter(b) {}
-
-      value_type operator*() const { return ((*fAIter) * (*fBIter)).magnitude(); }
-
-      WeightProviderIterator& operator++() { // prefix ++
-        ++fAIter;
-        ++fBIter;
-        return *this;
-      }
-
-      auto operator==(WeightProviderIterator other) { return fAIter == other.fAIter; }
-
-      auto operator!=(WeightProviderIterator other) { return !(*this == other); }
-    };
-
-  public:
-    NuclearComposition(std::vector<corsika::particles::Code> pComponents,
-                       std::vector<float> pFractions)
-        : fNumberFractions(pFractions)
-        , fComponents(pComponents)
-        , fAvgMassNumber(std::inner_product(
-              pComponents.cbegin(), pComponents.cend(), pFractions.cbegin(), 0.,
-              [](double x, double y) { return x + y; },
-              [](auto const& compID, auto const& fraction) {
-                return corsika::particles::GetNucleusA(compID) * fraction;
-              })) {
-      assert(pComponents.size() == pFractions.size());
-      auto const sumFractions =
-          std::accumulate(pFractions.cbegin(), pFractions.cend(), 0.f);
-
-      if (!(0.999f < sumFractions && sumFractions < 1.001f)) {
-        throw std::runtime_error("element fractions do not add up to 1");
-      }
-    }
-
-    auto size() const { return fNumberFractions.size(); }
-
-    auto const& GetFractions() const { return fNumberFractions; }
-    auto const& GetComponents() const { return fComponents; }
-    auto const GetAverageMassNumber() const { return fAvgMassNumber; }
-
-    template <class TRNG>
-    corsika::particles::Code SampleTarget(
-        std::vector<corsika::units::si::CrossSectionType> const& sigma,
-        TRNG& randomStream) const {
-      using namespace corsika::units::si;
-
-      assert(sigma.size() == fNumberFractions.size());
-
-      std::discrete_distribution channelDist(
-          WeightProviderIterator<decltype(fNumberFractions.begin()),
-                                 decltype(sigma.begin())>(fNumberFractions.begin(),
-                                                          sigma.begin()),
-          WeightProviderIterator<decltype(fNumberFractions.begin()),
-                                 decltype(sigma.end())>(fNumberFractions.end(),
-                                                        sigma.end()));
-
-      auto const iChannel = channelDist(randomStream);
-      return fComponents[iChannel];
-    }
-  };
-
-} // namespace corsika::environment
-
-#endif

Environment/VolumeTreeNode.h

-
-/*
- * (c) Copyright 2018 CORSIKA Project, corsika-project@lists.kit.edu
- *
- * See file AUTHORS for a list of contributors.
- *
- * 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.
- */
-
-#ifndef _include_VolumeTreeNode_H
-#define _include_VolumeTreeNode_H
-
-#include <corsika/geometry/Volume.h>
-#include <memory>
-#include <vector>
-
-namespace corsika::environment {
-
-  class Empty {}; //<! intended for usage as default template argument
-
-  template <typename IModelProperties = Empty>
-  class VolumeTreeNode {
-  public:
-    using VTNUPtr = std::unique_ptr<VolumeTreeNode<IModelProperties>>;
-    using IMPSharedPtr = std::shared_ptr<IModelProperties>;
-    using VolUPtr = std::unique_ptr<corsika::geometry::Volume>;
-
-    VolumeTreeNode(VolUPtr pVolume = nullptr)
-        : fGeoVolume(std::move(pVolume)) {}
-
-    bool Contains(corsika::geometry::Point const& p) const {
-      return fGeoVolume->Contains(p);
-    }
-
-    VolumeTreeNode<IModelProperties> const* Excludes(
-        corsika::geometry::Point const& p) const {
-      auto exclContainsIter =
-          std::find_if(fExcludedNodes.cbegin(), fExcludedNodes.cend(),
-                       [&](auto const& s) { return bool(s->Contains(p)); });
-
-      return exclContainsIter != fExcludedNodes.cend() ? *exclContainsIter : nullptr;
-    }
-
-    /** returns a pointer to the sub-VolumeTreeNode which is "responsible" for the given
-     * \class Point \arg p, or nullptr iff \arg p is not contained in this volume.
-     */
-    VolumeTreeNode<IModelProperties> const* GetContainingNode(
-        corsika::geometry::Point const& p) const {
-      if (!Contains(p)) { return nullptr; }
-
-      if (auto const childContainsIter =
-              std::find_if(fChildNodes.cbegin(), fChildNodes.cend(),
-                           [&](auto const& s) { return bool(s->Contains(p)); });
-          childContainsIter == fChildNodes.cend()) // not contained in any of the children
-      {
-        if (auto const exclContainsIter = Excludes(p)) // contained in any excluded nodes
-        {
-          return exclContainsIter->GetContainingNode(p);
-        } else {
-          return this;
-        }
-      } else {
-        return (*childContainsIter)->GetContainingNode(p);
-      }
-    }
-
-    /**
-     * Traverses the VolumeTree pre- or post-order and calls the functor  \p func for each
-     * node. \p func takes a reference to VolumeTreeNode as argument. The return value \p
-     * func is ignored.
-     */
-    template <typename TCallable, bool preorder = true>
-    void walk(TCallable func) {
-      if constexpr (preorder) {
-        func(*this);
-        std::for_each(fChildNodes.begin(), fChildNodes.end(),
-                      [&](auto& v) { v->walk(func); });
-      } else {
-        std::for_each(fChildNodes.begin(), fChildNodes.end(),
-                      [&](auto& v) { v->walk(func); });
-        func(*this);
-      }
-    }
-
-    void AddChild(VTNUPtr pChild) {
-      pChild->fParentNode = this;
-      fChildNodes.push_back(std::move(pChild));
-      // It is a bad idea to return an iterator to the inserted element
-      // because it might get invalidated when the vector needs to grow
-      // later and the caller won't notice.
-    }
-
-    void ExcludeOverlapWith(VTNUPtr const& pNode) {
-      fExcludedNodes.push_back(pNode.get());
-    }
-
-    auto* GetParent() const { return fParentNode; };
-
-    auto const& GetChildNodes() const { return fChildNodes; }
-
-    auto const& GetExcludedNodes() const { return fExcludedNodes; }
-
-    auto const& GetVolume() const { return *fGeoVolume; }
-
-    auto const& GetModelProperties() const { return *fModelProperties; }
-
-    IMPSharedPtr GetModelPropertiesPtr() const { return fModelProperties; }
-
-    template <typename TModelProperties, typename... Args>
-    auto SetModelProperties(Args&&... args) {
-      static_assert(std::is_base_of_v<IModelProperties, TModelProperties>,
-                    "unusable type provided");
-
-      fModelProperties = std::make_shared<TModelProperties>(std::forward<Args>(args)...);
-      return fModelProperties;
-    }
-
-    void SetModelProperties(IMPSharedPtr ptr) { fModelProperties = ptr; }
-
-    template <class MediumType, typename... Args>
-    static auto CreateMedium(Args&&... args) {
-      static_assert(std::is_base_of_v<IMediumModel, MediumType>,
-                    "unusable type provided, needs to be derived from \"IMediumModel\"");
-
-      return std::make_shared<MediumType>(std::forward<Args>(args)...);
-    }
-
-  private:
-    std::vector<VTNUPtr> fChildNodes;
-    std::vector<VolumeTreeNode<IModelProperties> const*> fExcludedNodes;
-    VolumeTreeNode<IModelProperties> const* fParentNode = nullptr;
-    VolUPtr fGeoVolume;
-    IMPSharedPtr fModelProperties;
-  };
-
-} // namespace corsika::environment
-
-#endif

Environment/testEnvironment.cc

-
-/*
- * (c) Copyright 2018 CORSIKA Project, corsika-project@lists.kit.edu
- *
- * See file AUTHORS for a list of contributors.
- *
- * 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.
- */
-
-#define CATCH_CONFIG_MAIN // This tells Catch to provide a main() - only do this in one
-                          // cpp file
-
-#include <corsika/environment/DensityFunction.h>
-#include <corsika/environment/FlatExponential.h>
-#include <corsika/environment/HomogeneousMedium.h>
-#include <corsika/environment/IMediumModel.h>
-#include <corsika/environment/InhomogeneousMedium.h>
-#include <corsika/environment/LinearApproximationIntegrator.h>
-#include <corsika/environment/NuclearComposition.h>
-#include <corsika/environment/VolumeTreeNode.h>
-#include <corsika/geometry/Line.h>
-#include <corsika/geometry/RootCoordinateSystem.h>
-#include <corsika/geometry/Vector.h>
-#include <corsika/particles/ParticleProperties.h>
-#include <corsika/units/PhysicalUnits.h>
-#include <catch2/catch.hpp>
-
-using namespace corsika::geometry;
-using namespace corsika::environment;
-using namespace corsika::particles;
-using namespace corsika::units::si;
-using namespace corsika;
-
-CoordinateSystem const& gCS =
-    RootCoordinateSystem::GetInstance().GetRootCoordinateSystem();
-
-Point const gOrigin(gCS, {0_m, 0_m, 0_m});
-
-TEST_CASE("HomogeneousMedium") {
-  NuclearComposition const protonComposition(std::vector<Code>{Code::Proton},
-                                             std::vector<float>{1.f});
-  HomogeneousMedium<IMediumModel> const medium(19.2_g / cube(1_cm), protonComposition);
-}
-
-TEST_CASE("FlatExponential") {
-  NuclearComposition const protonComposition(std::vector<Code>{Code::Proton},
-                                             std::vector<float>{1.f});
-
-  Vector const axis(gCS, QuantityVector<dimensionless_d>(0, 0, 1));
-  LengthType const lambda = 3_m;
-  auto const rho0 = 1_g / units::si::detail::static_pow<3>(1_cm);
-  FlatExponential<IMediumModel> const medium(gOrigin, axis, rho0, lambda,
-                                             protonComposition);
-  auto const tEnd = 5_s;
-
-  SECTION("horizontal") {
-    Line const line(gOrigin, Vector<SpeedType::dimension_type>(
-                                 gCS, {20_cm / second, 0_m / second, 0_m / second}));
-    Trajectory<Line> const trajectory(line, tEnd);
-
-    REQUIRE((medium.IntegratedGrammage(trajectory, 2_m) / (rho0 * 2_m)) == Approx(1));
-    REQUIRE((medium.ArclengthFromGrammage(trajectory, rho0 * 5_m) / 5_m) == Approx(1));
-  }
-
-  SECTION("vertical") {
-    Line const line(gOrigin, Vector<SpeedType::dimension_type>(
-                                 gCS, {0_m / second, 0_m / second, 5_m / second}));
-    Trajectory<Line> const trajectory(line, tEnd);
-    LengthType const length = 2 * lambda;
-    GrammageType const exact = rho0 * lambda * (exp(length / lambda) - 1);
-
-    REQUIRE((medium.IntegratedGrammage(trajectory, length) / exact) == Approx(1));
-    REQUIRE((medium.ArclengthFromGrammage(trajectory, exact) / length) == Approx(1));
-  }
-
-  SECTION("escape grammage") {
-    Line const line(gOrigin, Vector<SpeedType::dimension_type>(
-                                 gCS, {0_m / second, 0_m / second, -5_m / second}));
-    Trajectory<Line> const trajectory(line, tEnd);
-
-    GrammageType const escapeGrammage = rho0 * lambda;
-
-    REQUIRE(trajectory.NormalizedDirection().dot(axis).magnitude() < 0);
-    REQUIRE(medium.ArclengthFromGrammage(trajectory, 1.2 * escapeGrammage) ==
-            std::numeric_limits<typename GrammageType::value_type>::infinity() * 1_m);
-  }
-
-  SECTION("inclined") {
-    Line const line(gOrigin, Vector<SpeedType::dimension_type>(
-                                 gCS, {0_m / second, 5_m / second, 5_m / second}));
-    Trajectory<Line> const trajectory(line, tEnd);
-    double const cosTheta = M_SQRT1_2;
-    LengthType const length = 2 * lambda;
-    GrammageType const exact =
-        rho0 * lambda * (exp(cosTheta * length / lambda) - 1) / cosTheta;
-    REQUIRE((medium.IntegratedGrammage(trajectory, length) / exact) == Approx(1));
-    REQUIRE((medium.ArclengthFromGrammage(trajectory, exact) / length) == Approx(1));
-  }
-}
-
-auto constexpr rho0 = 1_kg / 1_m / 1_m / 1_m;
-
-struct Exponential {
-  auto operator()(corsika::geometry::Point const& p) const {
-    return exp(p.GetCoordinates()[0] / 1_m) * rho0;
-  }
-
-  template <int N>
-  auto Derivative(Point const& p, Vector<dimensionless_d> const& v) const {
-    return v.GetComponents()[0] * (*this)(p) /
-           corsika::units::si::detail::static_pow<N>(1_m);
-  }
-
-  auto FirstDerivative(Point const& p, Vector<dimensionless_d> const& v) const {
-    return Derivative<1>(p, v);
-  }
-
-  auto SecondDerivative(Point const& p, Vector<dimensionless_d> const& v) const {
-    return Derivative<2>(p, v);
-  }
-};
-
-TEST_CASE("InhomogeneousMedium") {
-  Vector direction(gCS, QuantityVector<dimensionless_d>(1, 0, 0));
-
-  Line line(gOrigin, Vector<SpeedType::dimension_type>(
-                         gCS, {20_m / second, 0_m / second, 0_m / second}));
-
-  auto const tEnd = 5_s;
-  Trajectory<Line> const trajectory(line, tEnd);
-
-  Exponential const e;
-  DensityFunction<decltype(e), LinearApproximationIntegrator> const rho(e);
-
-  SECTION("DensityFunction") {
-    REQUIRE(e.Derivative<1>(gOrigin, direction) / (1_kg / 1_m / 1_m / 1_m / 1_m) ==
-            Approx(1));
-    REQUIRE(rho.EvaluateAt(gOrigin) == e(gOrigin));
-  }
-
-  auto const exactGrammage = [](auto l) { return 1_m * rho0 * (exp(l / 1_m) - 1); };
-  auto const exactLength = [](auto X) { return 1_m * log(1 + X / (rho0 * 1_m)); };
-
-  auto constexpr l = 15_cm;
-
-  NuclearComposition const composition{{Code::Proton}, {1.f}};
-  InhomogeneousMedium<IMediumModel, decltype(rho)> const inhMedium(composition, rho);
-
-  SECTION("Integration") {
-    REQUIRE(rho.IntegrateGrammage(trajectory, l) / exactGrammage(l) ==
-            Approx(1).epsilon(1e-2));
-    REQUIRE(rho.ArclengthFromGrammage(trajectory, exactGrammage(l)) /
-                exactLength(exactGrammage(l)) ==
-            Approx(1).epsilon(1e-2));
-    REQUIRE(rho.MaximumLength(trajectory, 1e-2) >
-            l); // todo: write reasonable test when implementation is working
-
-    REQUIRE(rho.IntegrateGrammage(trajectory, l) ==
-            inhMedium.IntegratedGrammage(trajectory, l));
-    REQUIRE(rho.ArclengthFromGrammage(trajectory, 20_g / (1_cm * 1_cm)) ==
-            inhMedium.ArclengthFromGrammage(trajectory, 20_g / (1_cm * 1_cm)));
-  }
-}

Framework/CMakeLists.txt

-
-add_subdirectory (Utilities)
-add_subdirectory (Units)
-add_subdirectory (Geometry)
-add_subdirectory (Particles)
-add_subdirectory (Logging)
-add_subdirectory (StackInterface)
-add_subdirectory (ProcessSequence)
-add_subdirectory (Cascade)
-add_subdirectory (Random)

Framework/Cascade/CMakeLists.txt

-
-# namespace of library -> location of header files
-set (
-  CORSIKAcascade_NAMESPACE
-  corsika/cascade
-  )
-
-# header files of this library
-set (
-  CORSIKAcascade_HEADERS
-  Cascade.h
-  )
-
-add_library (CORSIKAcascade INTERFACE)
-
-CORSIKA_COPY_HEADERS_TO_NAMESPACE (CORSIKAcascade ${CORSIKAcascade_NAMESPACE} ${CORSIKAcascade_HEADERS})
-
-# include directive for upstream code
-target_include_directories (
-  CORSIKAcascade
-  INTERFACE
-  $<BUILD_INTERFACE:${PROJECT_BINARY_DIR}/include>
-  $<INSTALL_INTERFACE:include/>
-  )
-
-# install library
-install (
-  FILES ${CORSIKAcascade_HEADERS}
-  DESTINATION include/${CORSIKAcascade_NAMESPACE}
-  )
-
-# ----------------
-# code unit testing
-add_executable (
-  testCascade
-  testCascade.cc
-  )
-
-target_link_libraries (
-  testCascade
-  #  CORSIKAutls
-  CORSIKArandom
-  ProcessSibyll
-  CORSIKAcascade
-  ProcessStackInspector
-  ProcessTrackingLine
-  CORSIKAstackinterface
-  CORSIKAprocesses
-  CORSIKAparticles
-  CORSIKAgeometry
-  CORSIKAenvironment
-  CORSIKAprocesssequence
-  CORSIKAunits
-  CORSIKAthirdparty # for catch2
-  )
-CORSIKA_ADD_TEST(testCascade)