Python/README.md

-# Python for CORSIKA8
-
-![License](https://img.shields.io/badge/license-GPL3-blue)
-![Python](https://img.shields.io/badge/python-3.6%20%7C%203.7%20%7C%203.8-blue)
-[![Code style: black](https://img.shields.io/badge/code%20style-black-000000.svg)](https://github.com/psf/black)
-
-## Installation
-
-To install the CORSIKA 8 Python library, you will need Python >= 3.6. The below
-instructions are assuming that `python` refers to Python 3.\*. If `python` still
-refers to Python 2.\*, please replace `python` with `python3` and `pip` with
-`pip3`.
-
-To install this package, change into the `corsika/Python` directory and run
-
-    pip install --user -e .
-
-You can verify that the installation was successful by trying to import `corsika`
-
-    python -c 'import corsika'
-
-## Running the Tests
-    
-If you wish to develop new features in `corsika`, you will also need to install
-some additional dependencies so you can run our unit tests. These can be
-installed with
-
-    pip install --user -e '.[test]'
-
-Once that is completed, you can run the unit tests directory from the `corsika` directory
-
-    python -m pytest tests
-
-## Development Guide
-
-`corsika` is a type Python package with type-checking using
-[Mypy](mypy-lang.org) and all code is formatted with the
-[Black](https://black.readthedocs.io/en/stable/) formatter.
-
-All code must pass these checks, as well as several others, before it can be
-merged into CORSIKA. A `Makefile` is provided to simplify running these in
-tasks.
-
-    make mypy   # this will run the mypy type checker
-    
-    make isort  # this will automatically sort `import` statements
-    
-    make black  # this will automatically format your code
-    
-    make flake  # this checks that all code passes our linting rules
-    
-    make tests  # this runs our test suite
-    
-Lastly, we provide an `all` target that runs *all* of the above in the sequence shown above:
-
-    make all  # run all of the above tests

Python/corsika/__init__.py

-"""
- (c) Copyright 2018 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.
-"""
-
-__version__: str = "8.0.0-alpha"

Python/tests/test_corsika.py

-"""
- (c) Copyright 2018 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.
-"""
-import corsika
-
-
-def test_corsika_version() -> None:
-    """
-    Check the current CORSIKA version.
-    """
-    assert corsika.__version__ == "8.0.0-alpha"

README.md

-# CORSIKA 8 Framework for Particle Cascades in Astroparticle Physics 
+# CORSIKA 8 Framework for Particle Cascades in Astroparticle Physics
 
-The purpose of CORSIKA is to simulate any particle cascades in
-astroparticle physics or astrophysical context. A lot of emphasis is
+The purpose of CORSIKA 8 is to simulate any particle cascades in
+astroparticle physics or astrophysical context. A lot of emphasis has been
 put on modularity, flexibility, completeness, validation and
-correctness. To boost computational efficiency different techniques
+correctness. To boost computational efficiency, different techniques
 are provided, like thinning or cascade equations. The aim is that
-CORSIKA remains the most comprehensive framework for simulating
+CORSIKA 8 remains the most comprehensive framework for simulating
 particle cascades with stochastic and continuous processes.
 
 The software makes extensive use of static design patterns and
 compiler optimization. Thus, the most fundamental configuration
-decision of the user must be performed at compile time. At run time
-only specific model parameters can still be changed.
+decisions of the user must be performed at compile time. At run time,
+model parameters can still be changed.
 
-CORSIKA 8 is by default released under the GPLv3 license. See [license
-file](https://gitlab.ikp.kit.edu/AirShowerPhysics/corsika/blob/master/LICENSE)
+CORSIKA 8 is by default released under the BSD 3-Clause License. See [license
+file](https://gitlab.iap.kit.edu/AirShowerPhysics/corsika/blob/master/LICENSE)
 which is part of every release and the source code.
 
 If you use, or want to refer to, CORSIKA 8 please cite ["Towards a Next
 Generation of CORSIKA: A Framework for the Simulation of Particle
-Cascades in Astroparticle Physics", Comput.Softw.Big Sci. 3 (2019)
-2](https://doi.org/10.1007/s41781-018-0013-0). We kindly ask (and
-require) any relevant improvement or addition to be offered or
+Cascades in Astroparticle Physics", Comput. Softw. Big Sci. 3 (2019)
+2](https://doi.org/10.1007/s41781-018-0013-0) as well as
+["Simulating radio emission from particle cascades with CORSIKA 8", Astropart. Phys. 166 (2025)
+103072](https://doi.org/10.1016/j.astropartphys.2024.103072).
+
+We kindly ask (and require) any relevant improvement or addition to be offered or
 contributed to the main CORSIKA 8 repository for the benefit of the
 whole community.
 
-When you plan to contribute to CORSIKA 8 check the guidelines outlined here:
+CORSIKA 8 makes use of various third-party code, in particular interaction
+models. Please check the [using and collaborating
+agreement](https://gitlab.iap.kit.edu/AirShowerPhysics/corsika/blob/master/USING_COLLABORATING.md)
+for further information on this topic.
+
+If you plan to contribute to CORSIKA 8, please check the guidelines outlined here:
 [coding
-guidelines](https://gitlab.ikp.kit.edu/AirShowerPhysics/corsika/blob/master/CONTRIBUTING.md). Code
-that fails the review by the CORSIKA author group must be improved
+guidelines](https://gitlab.iap.kit.edu/AirShowerPhysics/corsika/blob/master/CONTRIBUTING.md). Code
+that fails the review by the CORSIKA 8 author group must be improved
 before it can be merged in the official code base. After your code has
-been accepted and merged you become a contributor of the CORSIKA 8
-project (code author). 
-
-IMPORTANT: Before you contribute, you need to read and agree to the
-[collaboration
-agreement](https://gitlab.ikp.kit.edu/AirShowerPhysics/corsika/blob/master/COLLABORATION_AGREEMENT.md). The
-agreement can be discussed, and eventually improved.
+been accepted and merged, you become a contributor of the CORSIKA 8
+project (code author).
 
-We also want to point you to the [MCnet
-guidelines](https://gitlab.ikp.kit.edu/AirShowerPhysics/corsika/blob/master/MCNET_GUIDELINES),
-which are very useful also for us.
+IMPORTANT: Before you contribute, you need to read and agree to the conditions set out in the
+[using and collaborating
+agreement](https://gitlab.iap.kit.edu/AirShowerPhysics/corsika/blob/master/USING_COLLABORATING.md).
+The agreement can be discussed, and eventually improved if necessary.
 
 
 ## Get in contact
-  * Connect to https://gitlab.ikp.kit.edu register yourself and join the "Air Shower Physics" group
+  * Join our chat threads using Mattermost via this [invite link](https://mattermost.hzdr.de/signup_user_complete/?id=xtdd8jyt6trbiezt71gaz3z4ge&md=link&sbr=su). Click the `GitLab` button, then `Sign in with Helmholtz ID`. You will be able to make an account by either finding your institution, or using your e.g. ORCID, GitHub, or Google account.
+  * Connect to https://gitlab.iap.kit.edu, register yourself and join the "Air Shower Physics" group. Write to us on Mattermost (in the User Questions channel), or directly contact one of the [steering comittee members](https://gitlab.iap.kit.edu/AirShowerPhysics/corsika/-/wikis/Steering-Committee) in case there are problems with that.
   * Connect to corsika-devel@lists.kit.edu (self-register at
     https://www.lists.kit.edu/sympa/subscribe/corsika-devel) to get in
-    touch with the project
+    touch with the project.
 
 
-## Prerequisites
+## Installation
 
-CORSIKA 8 is tested regularly via gitlab-CI using recent gcc and clang
-versions.  Additional software prerequisites: cmake, g++, git.
-Furthermore, eigen3, boost, catch2, spdlog are shipped in the
-ThirdParty directory, so an installation on the system is optional.
-Also Pythia 8, CONEX and PROPOSAL are distributed in the ThirdParty
-folder. You may also install those packages on your system and use
-those; we test with Pythia version 8.235.
+CORSIKA 8 is tested regularly at least on `gcc11.0.0` and `clang-14.0.0`.
 
-On a bare Ubuntu 18.04, just add:
-```
-sudo apt install cmake g++ git
+### Prerequisites
+
+You will also need:
+
+- Python 3 (supported versions are Python >= 3.6), with pip
+- cmake > 3.4
+- git
+- g++, gfortran, binutils, make
+- optional: FLUKA (see below)
+
+
+On a bare Ubuntu machine, just add:
+``` shell
+sudo apt-get install python3 python3-pip cmake g++ gfortran git doxygen graphviz
 ```
-add ```libeigen3-dev``` if you want to use system version of eigen3.
 
-If you work with FreeBSD, run:
+### Creating a virtual environment and Conan
+
+It is recommended that you install CORSIKA 8 and its dependencies within a python3 virtual environment.
+To do so, you can run the following.
+``` shell
+# Create the environment using your native python3 binary
+python3 -m venv /path/to/new/virtual/environment/corsika-8
+# Load the environment (should be run each time you open a new terminal)
+source /path/to/new/virtual/environment/corsika-8/bin/activate
+
 ```
-pkg install git cmake python3 flang
+
+You will need to load the environment each time that you open a new terminal.
+
+CORSIKA 8 uses the [conan](https://conan.io/) package manager to
+manage our dependencies. Currently, version 2.50.0 or higher is required.
+**Note**: if you are NOT using a virtual environment, you may want to use the `pip install --user` flag.
+
+``` shell
+pip install conan particle==0.25.1 numpy
 ```
-or add ```boost-libs``` and ```eigen``` if you want to use the system versions. 
 
-## Installation
+### Enabling FLUKA support
+
+For legal reasons we do not distribute/bundle FLUKA together with CORSIKA 8.
+As FLUKA is the standard low-energy hadronic interaction model for CORSIKA 8, you have to download
+it separately from (http://www.fluka.org/), which requires registering there as FLUKA user.
+The following should be done *before* compiling CORSIKA 8:
 
+ 1. Note your system's version of gfortran (`gfortran --version`) and glibc (`ldd --version`)
+ 2. Download the FLUKA __binary__, ensuring that it matches the versions you found above
+ 3. Download the FLUKA __data file__ (will be named something similar to __fluka20xy.z-data.tar.gz__).
+ 4. Un-tar the files that you downloaded using `tar -xf <filename>`
+ 5. Set environmental variables `export FLUFOR=gfortran` and `export FLUPRO=<path to where you unzipped the files>`. Note that the `FLUPRO` directory should contain __libflukahp.a__. Both of these variables will have to be set every time you open a new terminal.
+ 6. Go to the `FLUPRO` directory and run `make`. This will compile an exe, __flukahp__, in your current directory.
+ 7. Follow the normal steps to compile CORSIKA 8 (see below).
 
-Follow these steps to download and install CORSIKA 8 milestone2
+When you later install CORSIKA 8, you should see a message during the __cmake__ step indicating the FLUKA was correctly found.
+
+``` shell
+libflukahp.a found in directory <some location here> via FLUPRO environment variable
+FLUKA support is enabled.
 ```
-git clone --recursive https://gitlab.ikp.kit.edu/AirShowerPhysics/corsika.git
-cd corsika
-mkdir ../corsika-build
-cd ../corsika-build
-cmake ../corsika -DCMAKE_INSTALL_PREFIX=../corsika-install
-make -j8
+
+### Compiling CORSIKA 8
+
+Once Conan is installed and FLUKA provided, follow these steps to download and install CORSIKA 8:
+
+``` shell
+cd ./top/directory/for/corsika/installation
+git clone --recursive git@gitlab.iap.kit.edu:AirShowerPhysics/corsika.git
+# Or for https: git clone --recursive https://gitlab.iap.kit.edu/AirShowerPhysics/corsika.git
+mkdir corsika-build
+cd corsika-build
+../corsika/conan-install.sh --source-directory ../corsika --release-with-debug
+# conan-install.sh takes required options from command line to install dependencies for 'Debug', 'Release' and 'RelWithDebInfo' builds.
+../corsika/corsika-cmake.sh -c "-DCMAKE_BUILD_TYPE="RelWithDebInfo" -DWITH_FLUKA=ON -DCMAKE_INSTALL_PREFIX=../corsika-install"
+make -j4  #The number should match the number of available cores on your machine
 make install
-make test
 ```
-and if you want to see how the first simple hadron cascade develops, 
-see `Documentation/Examples/cascade_example.cc` for a starting point. 
 
-Run the cascade_example with: 
+## Alternate installation using docker containers
+
+There are docker containers prepared that bring all the environment and packages you need to run CORSIKA. See [docker hub](https://hub.docker.com/repository/docker/corsika/devel) for a complete overview.
+
+### Prerequisites
+
+You only need docker, e.g. on Ubuntu: `sudo apt-get install docker` and of course root access.
+
+### Compiling
+
+Follow these steps to download and install CORSIKA 8, master development version
+```shell
+cd ./top/directory/for/corsika/installation
+git clone --recursive git@gitlab.iap.kit.edu:AirShowerPhysics/corsika.git
+sudo docker run -v $PWD:/corsika -it corsika/devel:clang-8 /bin/bash
+mkdir corsika-build
+cd corsika-build
+../corsika/conan-install.sh --source-directory ../corsika --release-with-debug
+# conan-install.sh takes required options from command line to install dependencies for 'Debug', 'Release' and 'RelWithDebInfo' builds.
+../corsika/corsika-cmake.sh -c "-DCMAKE_BUILD_TYPE="RelWithDebInfo" -DWITH_FLUKA=ON -DCMAKE_INSTALL_PREFIX=../corsika-install"
+make -j4  #The number should match the number of available cores on your machine
+make install
 ```
-cd ../corsika-install
-share/examples/cascade_example
+
+## Running Unit Tests
+
+To run the unit tests, do the following.
+
+```shell
+cd ./corsika-build
+ctest -j4  #The number should match the number of available cores on your machine
 ```
 
-Visualize output (needs gnuplot installed): 
+## Running applications and examples
+
+### Standard applications
+
+Applications for standard use-cases are located in the `applications` directory.
+These are example scripts that can be used directly or slightly modified for your use case.
+See [applications/README.md] for more.
+The applications are compiled automatically after running `make` and will appear your `corsika-build/bin` directory.
+After running `make install` the binaries will also be copied into your `corsika-install/bin` directory as well.
+
+
+For example, from inside your `corsika-install/bin` directory, run
+```shell
+c8_air_shower --pdg 2212 -E 1e5 -f my_shower
 ```
-bash share/tools/plot_tracks.sh tracks.dat 
-firefox tracks.dat.gif 
+This will run a vertical 100 TeV proton shower and will create and put the output into `./my_shower`.
+
+
+### Building the examples
+
+Unlike the applications, the examples must be compiled as a second step.
+From your top corsika directory, (the one that includes `corsika-build` and `corsika-install`) run
+```shell
+export CONAN_DEPENDENCIES=$PWD/corsika-install/lib/cmake/dependencies
+cmake -DCMAKE_TOOLCHAIN_FILE=${CONAN_DEPENDENCIES}/conan_toolchain.cmake -DCMAKE_PREFIX_PATH=${CONAN_DEPENDENCIES} -DCMAKE_POLICY_DEFAULT_CMP0091=NEW -DCMAKE_BUILD_TYPE=RelWithDebInfo -Dcorsika_DIR=$PWD/corsika-build -DWITH_FLUKA=ON -S $PWD/corsika/examples -B $PWD/corsika-build-examples
+cd corsika-build-examples
+make -j4 #The number should match the number of available cores on your machine
 ```
 
-Or also consider the `vertical_eas` example in the same directory, which can 
-be configured with command line options. 
+You can run the examples by using the binaries in `corsika-build-examples/bin/`.
+For example:
+```shell
+corsika-build-examples/bin/known_particles
+```
+This will print out all of the particles that are known by CORSIKA.
 
 
 ### Generating doxygen documentation
 
-To generate the documentation, you need doxygen and graphviz. On Ubuntu 18.04, do:
-```
+To generate the documentation, you need doxygen and graphviz. If you work with
+the docker corsika/devel containers this is already included.
+Otherwise, e.g. on Ubuntu machines, do:
+```shell
 sudo apt-get install doxygen graphviz
 ```
-Switch to the corsika build directory and do
-```
-make doxygen
+Switch to the `corsika-build` directory and do
+```shell
+make docs
 make install
 ```
-browse with firefox:
+open with firefox:
+```shell
+firefox ../corsika-install/share/corsika/doc/html/index.html
 ```
-firefox ../corsika-install/share/doc/html/index.html
-```
-

SCIENTIFIC_AUTHORS

-This is currently derived from
-https://gitlab.ikp.kit.edu/AirShowerPhysics/corsika/wikis/ICRC2019
-and will be revised whenever needed:
-
-- Luisa Arrabito (1)
-- Dominik Baack (2)
-- Johan Bregeon (1)
-- Hans Dembinski (3)
-- Ralph Engel (4)
-- Dieter Heck (4)
-- Tim Huege (4,5)
-- Lukas Nellen (6)
-- Tanguy Pierog (4)
-- Maximilian Reininghaus (4,7)
-- Felix Riehn (8)
-- Ralf Ulrich (4)
-- Darko Veberic (4)
-
-Affiliations:
-
-- University of Montpellier, France
-- TU Dortmund University, Germany
-- Max Planck Institute for Nuclear Physics, Heidelberg, Germany
-- Karlsruhe Institute of Technology (KIT), Germany
-- Vrije Universiteit Brussel, Belgium
-- National Autonomous University of Mexico (UNAM)
-- Instituto de Tecnologías en Detección y Astropartículas (CNEA, CONICET, UNSAM), Buenos Aires, Argentina
-- Laboratory of Instrumentation and Experimental Particles (LIP), Portugal

Setup/CMakeLists.txt

-set (
-  SETUP_HEADERS
-  SetupStack.h
-  SetupEnvironment.h
-  SetupTrajectory.h
-  )
-
-set (
-  SETUP_NAMESPACE
-  corsika/setup
-  )
-
-add_library (CORSIKAsetup INTERFACE)
-CORSIKA_COPY_HEADERS_TO_NAMESPACE (CORSIKAsetup ${SETUP_NAMESPACE} ${SETUP_HEADERS})
-
-target_link_libraries (
-  CORSIKAsetup
-  INTERFACE
-  CORSIKAgeometry
-  NuclearStackExtension
-  GeometryNodeStackExtension
-  CORSIKAhistory
-  )
-
-if (WITH_HISTORY)
-  target_compile_definitions (CORSIKAsetup INTERFACE "WITH_HISTORY")
-endif (WITH_HISTORY)
-
-target_include_directories (
-  CORSIKAsetup
-  INTERFACE 
-  $<BUILD_INTERFACE:${PROJECT_BINARY_DIR}/include>
-  $<INSTALL_INTERFACE:include/include>
-  )
-
-install (
-  FILES ${SETUP_HEADERS} 
-  DESTINATION include/${SETUP_NAMESPACE}
-  )

Setup/SetupEnvironment.h

-/*
- * (c) Copyright 2018 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.
- */
-
-#pragma once
-
-#include <corsika/environment/Environment.h>
-#include <corsika/environment/IMediumModel.h>
-#include <corsika/environment/NameModel.h>
-
-namespace corsika::setup {
-  using IEnvironmentModel = environment::IMediumModel;
-  using SetupEnvironment = environment::Environment<IEnvironmentModel>;
-} // namespace corsika::setup

Setup/SetupStack.h

-/*
- * (c) Copyright 2018 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.
- */
-
-#pragma once
-
-#include <corsika/stack/CombinedStack.h>
-#include <corsika/stack/node/GeometryNodeStackExtension.h>
-#include <corsika/stack/nuclear_extension/NuclearStackExtension.h>
-#include <corsika/stack/history/HistorySecondaryProducer.hpp>
-#include <corsika/stack/history/HistoryStackExtension.hpp>
-
-#include <corsika/setup/SetupEnvironment.h>
-
-namespace corsika::setup {
-
-  namespace detail {
-
-    // ------------------------------------------
-    // add geometry node tracking data to stack:
-
-    // the GeometryNode stack needs to know the type of geometry-nodes from the
-    // environment:
-    template <typename TStackIter>
-    using SetupGeometryDataInterface =
-        typename stack::node::MakeGeometryDataInterface<TStackIter,
-                                                        setup::SetupEnvironment>::type;
-
-    // combine particle data stack with geometry information for tracking
-    template <typename TStackIter>
-    using StackWithGeometryInterface = corsika::stack::CombinedParticleInterface<
-        stack::nuclear_extension::ParticleDataStack::MPIType, SetupGeometryDataInterface,
-        TStackIter>;
-
-    using StackWithGeometry = corsika::stack::CombinedStack<
-        typename corsika::stack::nuclear_extension::ParticleDataStack::StackImpl,
-        corsika::stack::node::GeometryData<setup::SetupEnvironment>,
-        StackWithGeometryInterface>;
-
-    // ------------------------------------------
-    // Add [optional] history data to stack, too:
-
-    // combine dummy stack with geometry information for tracking
-    template <typename TStackIter>
-    using StackWithHistoryInterface = corsika::stack::CombinedParticleInterface<
-        StackWithGeometry::MPIType, history::HistoryEventDataInterface, TStackIter>;
-
-    using StackWithHistory =
-        corsika::stack::CombinedStack<typename StackWithGeometry::StackImpl,
-                                      history::HistoryEventData,
-                                      StackWithHistoryInterface>;
-
-  } // namespace detail
-
-  // ---------------------------------------
-  // this is the FINAL stack we use in C8:
-
-#ifdef WITH_HISTORY
-
-  /*
-   * the version with history
-   */
-  using Stack = detail::StackWithHistory;
-  template <typename T1, template <typename> typename M2>
-  using StackViewProducer = corsika::history::HistorySecondaryProducer<T1, M2>;
-
-  namespace detail {
-    /*
-      See Issue 161
-
-      unfortunately clang does not support this in the same way (yet) as
-      gcc, so we have to distinguish here. If clang cataches up, we
-      could remove the clang branch here and also in
-      corsika::Cascade. The gcc code is much more generic and
-      universal. If we could do the gcc version, we won't had to define
-      StackView globally, we could do it with MakeView whereever it is
-      actually needed. Keep an eye on this!
-    */
-#if defined(__clang__)
-    using TheStackView = corsika::stack::SecondaryView<
-        typename corsika::setup::Stack::StackImpl,
-        // CHECK with CLANG: corsika::setup::Stack::MPIType>;
-        corsika::setup::detail::StackWithHistoryInterface, StackViewProducer>;
-#elif defined(__GNUC__) || defined(__GNUG__)
-    using TheStackView =
-        corsika::stack::MakeView<corsika::setup::Stack, StackViewProducer>::type;
-#endif
-  } // namespace detail
-
-#else // WITH_HISTORY
-
-  /*
-   * the version without history
-   */
-  using Stack = detail::StackWithGeometry;
-  template <typename T1, template <typename> typename M2>
-  using StackViewProducer = corsika::stack::DefaultSecondaryProducer<T1, M2>;
-
-  namespace detail {
-    /*
-      See Issue 161
-
-      unfortunately clang does not support this in the same way (yet) as
-      gcc, so we have to distinguish here. If clang cataches up, we
-      could remove the clang branch here and also in
-      corsika::Cascade. The gcc code is much more generic and
-      universal. If we could do the gcc version, we won't had to define
-      StackView globally, we could do it with MakeView whereever it is
-      actually needed. Keep an eye on this!
-    */
-#if defined(__clang__)
-    using TheStackView =
-        corsika::stack::SecondaryView<typename corsika::setup::Stack::StackImpl,
-                                      // CHECK with CLANG:
-                                      // corsika::setup::Stack::MPIType>;
-                                      corsika::setup::detail::StackWithGeometryInterface>;
-#elif defined(__GNUC__) || defined(__GNUG__)
-    using TheStackView = corsika::stack::MakeView<corsika::setup::Stack>::type;
-#endif
-  } // namespace detail
-
-#endif
-
-  // ---------------------------------------
-  // this is the FINAL stackitertor (particle type) we use in C8:
-
-  using StackView = detail::TheStackView;
-
-} // namespace corsika::setup

Setup/SetupTrajectory.h

-/*
- * (c) Copyright 2018 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.
- */
-
-#pragma once
-
-#include <corsika/geometry/Helix.h>
-#include <corsika/geometry/Line.h>
-#include <corsika/geometry/Trajectory.h>
-
-#include <corsika/units/PhysicalUnits.h>
-
-// #include <variant>
-
-namespace corsika::setup {
-
-  /// definition of Trajectory base class, to be used in tracking and cascades
-  typedef corsika::geometry::Trajectory<corsika::geometry::Line> Trajectory;
-
-  /*
-  typedef std::variant<std::monostate, corsika::geometry::Trajectory<Line>,
-                       corsika::geometry::Trajectory<Helix>>
-                       Trajectory;
-
-  /// helper visitor to modify Particle by moving along Trajectory
-  template <typename Particle>
-  class ParticleUpdate {
-
-    Particle& fP;
-
-  public:
-    ParticleUpdate(Particle& p)
-        : fP(p) {}
-    void operator()(std::monostate const&) {}
-
-    template <typename T>
-    void operator()(T const& trajectory) {
-      fP.SetPosition(trajectory.GetPosition(1));
-    }
-  };
-
-  /// helper visitor to modify Particle by moving along Trajectory
-  class GetDuration {
-  public:
-    corsika::units::si::TimeType operator()(std::monostate const&) {
-      return 0 * corsika::units::si::second;
-    }
-    template <typename T>
-    corsika::units::si::TimeType operator()(T const& trajectory) {
-      return trajectory.GetDuration();
-    }
-  };
-  */
-} // namespace corsika::setup

Stack/CMakeLists.txt

-add_subdirectory (DummyStack)
-add_subdirectory (SuperStupidStack)
-add_subdirectory (NuclearStackExtension)
-add_subdirectory (GeometryNodeStackExtension)
-add_subdirectory (History)

Stack/DummyStack/CMakeLists.txt

-set (DummyStack_HEADERS DummyStack.h)
-set (DummyStack_NAMESPACE corsika/stack/dummy)
-
-add_library (DummyStack INTERFACE)
-
-CORSIKA_COPY_HEADERS_TO_NAMESPACE (DummyStack ${DummyStack_NAMESPACE} ${DummyStack_HEADERS})
-
-target_link_libraries (
-  DummyStack
-  INTERFACE
-  CORSIKAstackinterface
-  CORSIKAunits
-  CORSIKAparticles
-  )
-
-target_include_directories (
-  DummyStack
-  INTERFACE
-  $<BUILD_INTERFACE:${PROJECT_BINARY_DIR}/include>
-  $<INSTALL_INTERFACE:include>
-  )
-
-install (
-  FILES
-  ${DummyStack_HEADERS}
-  DESTINATION
-  include/${DummyStack_NAMESPACE}
-  )
-
-# ----------------
-# code unit testing
- CORSIKA_ADD_TEST(testDummyStack)
- target_link_libraries (
-   testDummyStack
-   DummyStack
-   CORSIKAtesting
-   )

Stack/DummyStack/DummyStack.h

-/*
- * (c) Copyright 2018 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.
- */
-
-#pragma once
-
-#include <corsika/logging/Logging.h>
-#include <corsika/particles/ParticleProperties.h>
-#include <corsika/stack/Stack.h>
-#include <corsika/units/PhysicalUnits.h>
-
-#include <string>
-#include <tuple>
-
-namespace corsika::stack {
-
-  namespace dummy {
-
-    /**
-     * Example of a particle object on the stack, with NO DATA.
-     */
-
-    /**
-       however, conceptually we need to provide fake data. A stack without data does not
-       work...
-     */
-
-    struct NoData { /* nothing */
-      int nothing = 0;
-    };
-
-    template <typename StackIteratorInterface>
-    class ParticleInterface
-        : public corsika::stack::ParticleBase<StackIteratorInterface> {
-
-    protected:
-      using corsika::stack::ParticleBase<StackIteratorInterface>::GetStack;
-      using corsika::stack::ParticleBase<StackIteratorInterface>::GetStackData;
-
-    public:
-      using corsika::stack::ParticleBase<StackIteratorInterface>::GetIndex;
-
-    public:
-      void SetParticleData(const std::tuple<NoData>& /*v*/) {}
-      void SetParticleData(ParticleInterface<StackIteratorInterface>& /*parent*/,
-                           const std::tuple<NoData>& /*v*/) {}
-
-      std::string as_string() const { return "dummy-data"; }
-    };
-
-    /**
-     *
-     * Memory implementation of the most simple (no-data) particle stack object.
-     */
-
-    class DummyStackImpl {
-
-    public:
-      void Init() { entries_ = 0; }
-
-      void Clear() { entries_ = 0; }
-
-      int GetSize() const { return entries_; }
-      int GetCapacity() const { return entries_; }
-
-      /**
-       *   Function to copy particle at location i2 in stack to i1
-       */
-      void Copy(const int /*i1*/, const int /*i2*/) {}
-
-      void IncrementSize() { entries_++; }
-      void DecrementSize() { entries_--; }
-
-    private:
-      int entries_ = 0;
-
-    }; // end class DummyStackImpl
-
-    typedef Stack<DummyStackImpl, ParticleInterface> DummyStack;
-
-  } // namespace dummy
-
-} // namespace corsika::stack

Stack/GeometryNodeStackExtension/CMakeLists.txt

-set (GeometryNodeStackExtension_HEADERS GeometryNodeStackExtension.h)
-set (GeometryNodeStackExtension_NAMESPACE corsika/stack/node)
-
-add_library (GeometryNodeStackExtension INTERFACE)
-
-CORSIKA_COPY_HEADERS_TO_NAMESPACE (GeometryNodeStackExtension ${GeometryNodeStackExtension_NAMESPACE} ${GeometryNodeStackExtension_HEADERS})
-
-target_link_libraries (
-  GeometryNodeStackExtension
-  INTERFACE
-  CORSIKAstackinterface
-  CORSIKAunits
-  CORSIKAparticles
-  CORSIKAgeometry
-  )
-
-target_include_directories (
-  GeometryNodeStackExtension
-  INTERFACE
-  $<BUILD_INTERFACE:${PROJECT_BINARY_DIR}/include>
-  $<INSTALL_INTERFACE:include>
-  )
-
-install (
-  FILES
-  ${GeometryNodeStackExtension_HEADERS}
-  DESTINATION
-  include/${GeometryNodeStackExtension_NAMESPACE}
-  )
-
-# ----------------
-# code unit testing
- CORSIKA_ADD_TEST(testGeometryNodeStackExtension)
- target_link_libraries (
-   testGeometryNodeStackExtension
-   GeometryNodeStackExtension
-   CORSIKAtesting
-   )

Stack/GeometryNodeStackExtension/GeometryNodeStackExtension.h

-/*
- * (c) Copyright 2018 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.
- */
-
-#pragma once
-
-#include <corsika/logging/Logging.h>
-#include <corsika/stack/Stack.h>
-
-#include <tuple>
-#include <utility>
-#include <vector>
-
-namespace corsika::stack::node {
-
-  /**
-   * @class GeometryDataInterface
-   *
-   * corresponding defintion of a stack-readout object, the iteractor
-   * dereference operator will deliver access to these function
-  // defintion of a stack-readout object, the iteractor dereference
-  // operator will deliver access to these function
-   */
-  template <typename T, typename TEnvType>
-  class GeometryDataInterface : public T {
-
-  protected:
-    using T::GetStack;
-    using T::GetStackData;
-
-  public:
-    using T::GetIndex;
-    using BaseNodeType = typename TEnvType::BaseNodeType;
-
-  public:
-    // default version for particle-creation from input data
-    void SetParticleData(const std::tuple<BaseNodeType const*> v) {
-      SetNode(std::get<0>(v));
-    }
-    void SetParticleData(GeometryDataInterface& parent,
-                         const std::tuple<BaseNodeType const*>) {
-      SetNode(parent.GetNode()); // copy Node from parent particle!
-    }
-    void SetParticleData() { SetNode(nullptr); }
-    void SetParticleData(GeometryDataInterface& parent) {
-      SetNode(parent.GetNode()); // copy Node from parent particle!
-    }
-
-    std::string as_string() const { return fmt::format("node={}", fmt::ptr(GetNode())); }
-
-    void SetNode(BaseNodeType const* v) { GetStackData().SetNode(GetIndex(), v); }
-    BaseNodeType const* GetNode() const { return GetStackData().GetNode(GetIndex()); }
-  };
-
-  // definition of stack-data object to store geometry information
-  template <typename TEnvType>
-
-  /**
-   * @class GeometryData
-   *
-   * definition of stack-data object to store geometry information
-   */
-  class GeometryData {
-
-  public:
-    using BaseNodeType = typename TEnvType::BaseNodeType;
-
-    // these functions are needed for the Stack interface
-    void Clear() { fNode.clear(); }
-    unsigned int GetSize() const { return fNode.size(); }
-    unsigned int GetCapacity() const { return fNode.size(); }
-    void Copy(const int i1, const int i2) { fNode[i2] = fNode[i1]; }
-    void Swap(const int i1, const int i2) { std::swap(fNode[i1], fNode[i2]); }
-
-    // custom data access function
-    void SetNode(const int i, BaseNodeType const* v) { fNode[i] = v; }
-    BaseNodeType const* GetNode(const int i) const { return fNode[i]; }
-
-    // these functions are also needed by the Stack interface
-    void IncrementSize() { fNode.push_back(nullptr); }
-    void DecrementSize() {
-      if (fNode.size() > 0) { fNode.pop_back(); }
-    }
-
-    // custom private data section
-  private:
-    std::vector<const BaseNodeType*> fNode;
-  };
-
-  template <typename T, typename TEnv>
-  struct MakeGeometryDataInterface {
-    typedef GeometryDataInterface<T, TEnv> type;
-  };
-
-} // namespace corsika::stack::node

Stack/History/CMakeLists.txt

-set (
-  HISTORY_HEADERS
-  EventType.hpp
-  Event.hpp
-  HistorySecondaryProducer.hpp
-  HistoryObservationPlane.hpp
-  HistoryStackExtension.hpp
-  SecondaryParticle.hpp
-  )
-
-set (
-  HISTORY_NAMESPACE
-  corsika/stack/history
-  )
-
-if (WITH_HISTORY)
-  set (
-    HISTORY_SOURCES
-    HistoryObservationPlane.cpp
-    )
-  add_library (CORSIKAhistory STATIC ${HISTORY_SOURCES})
-  set (IS_INTERFACE "")
-else (WITH_HISTORY)
-  add_library (CORSIKAhistory INTERFACE)
-  set (IS_INTERFACE "INTERFACE")
-endif (WITH_HISTORY)
-
-CORSIKA_COPY_HEADERS_TO_NAMESPACE (CORSIKAhistory ${HISTORY_NAMESPACE} ${HISTORY_HEADERS})
-
-target_link_libraries (
-  CORSIKAhistory
-  ${IS_INTERFACE}
-  CORSIKAparticles
-  CORSIKAgeometry
-  CORSIKAprocesssequence # for HistoryObservationPlane
-  CORSIKAsetup # for HistoryObservationPlane
-  CORSIKAlogging
-  SuperStupidStack
-  NuclearStackExtension # for testHistoryView.cc
-  C8::ext::boost # for HistoryObservationPlane
-  )
-
-target_include_directories (
-  CORSIKAhistory
-  INTERFACE 
-  $<BUILD_INTERFACE:${PROJECT_BINARY_DIR}/include>
-  $<INSTALL_INTERFACE:include/include>
-  )
-
-install (
-  FILES ${HISTORY_HEADERS} 
-  DESTINATION include/${HISTORY_NAMESPACE}
-  )
-
-# ----------------
-# code unit testing
-CORSIKA_ADD_TEST(testHistoryStack)
-target_link_libraries (
-  testHistoryStack
-  CORSIKAhistory
-  CORSIKAtesting
-  DummyStack
-  )
-CORSIKA_ADD_TEST(testHistoryView)
-target_link_libraries (
-  testHistoryView
-  CORSIKAhistory
-  CORSIKAtesting
-  )

Stack/History/HistoryObservationPlane.cpp

-/*
- * (c) Copyright 2020 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/logging/Logging.h>
-#include <corsika/stack/history/HistoryObservationPlane.hpp>
-
-#include <boost/histogram/ostream.hpp>
-
-#include <iomanip>
-#include <iostream>
-
-using namespace corsika::units::si;
-using namespace corsika::history;
-using namespace corsika;
-
-HistoryObservationPlane::HistoryObservationPlane(setup::Stack const& stack,
-                                                 geometry::Plane const& obsPlane,
-                                                 bool deleteOnHit)
-    : stack_{stack}
-    , plane_{obsPlane}
-    , deleteOnHit_{deleteOnHit} {}
-
-corsika::process::EProcessReturn HistoryObservationPlane::DoContinuous(
-    setup::Stack::ParticleType const& particle, setup::Trajectory const& trajectory) {
-  TimeType const timeOfIntersection =
-      (plane_.GetCenter() - trajectory.GetR0()).dot(plane_.GetNormal()) /
-      trajectory.GetV0().dot(plane_.GetNormal());
-
-  if (timeOfIntersection < TimeType::zero()) { return process::EProcessReturn::eOk; }
-
-  if (plane_.IsAbove(trajectory.GetR0()) == plane_.IsAbove(trajectory.GetPosition(1))) {
-    return process::EProcessReturn::eOk;
-  }
-
-  C8LOG_DEBUG(fmt::format("HistoryObservationPlane: Particle detected: pid={}",
-                          particle.GetPID()));
-
-  auto const pid = particle.GetPID();
-  if (particles::IsMuon(pid)) { fillHistoryHistogram(particle); }
-
-  if (deleteOnHit_) {
-    return process::EProcessReturn::eParticleAbsorbed;
-  } else {
-    return process::EProcessReturn::eOk;
-  }
-}
-
-LengthType HistoryObservationPlane::MaxStepLength(setup::Stack::ParticleType const&,
-                                                  setup::Trajectory const& trajectory) {
-  TimeType const timeOfIntersection =
-      (plane_.GetCenter() - trajectory.GetR0()).dot(plane_.GetNormal()) /
-      trajectory.GetV0().dot(plane_.GetNormal());
-
-  if (timeOfIntersection < TimeType::zero()) {
-    return std::numeric_limits<double>::infinity() * 1_m;
-  }
-
-  auto const pointOfIntersection = trajectory.GetPosition(timeOfIntersection);
-  return (trajectory.GetR0() - pointOfIntersection).norm() * 1.0001;
-}
-
-void HistoryObservationPlane::fillHistoryHistogram(
-    setup::Stack::ParticleType const& muon) {
-  double const muon_energy = muon.GetEnergy() / 1_GeV;
-
-  int genctr{0};
-  Event const* event = muon.GetEvent().get();
-  while (event) {
-    auto const projectile = stack_.cfirst() + event->projectileIndex();
-    if (event->eventType() == EventType::Interaction) {
-      genctr++;
-      double const projEnergy = projectile.GetEnergy() / 1_GeV;
-      int const pdg = static_cast<int>(particles::GetPDG(projectile.GetPID()));
-
-      histogram_(muon_energy, projEnergy, pdg);
-    }
-    event = event->parentEvent().get(); // projectile.GetEvent().get();
-  }
-}

Stack/History/testHistoryView.cc

-/*
- * (c) Copyright 2018 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/stack/history/Event.hpp>
-#include <corsika/stack/history/HistorySecondaryProducer.hpp>
-#include <corsika/stack/history/HistoryStackExtension.hpp>
-
-#include <corsika/stack/CombinedStack.h>
-#include <corsika/stack/dummy/DummyStack.h>
-#include <corsika/stack/nuclear_extension/NuclearStackExtension.h>
-
-#include <corsika/logging/Logging.h>
-
-#include <catch2/catch.hpp>
-
-using namespace corsika;
-using namespace corsika::geometry;
-using namespace corsika::units::si;
-
-/**
-   Need to replicate setup::SetupStack in a maximally simplified
-   way, but with real particle data
- */
-
-// combine dummy stack with geometry information for tracking
-template <typename TStackIter>
-using StackWithHistoryInterface = corsika::stack::CombinedParticleInterface<
-    stack::nuclear_extension::ParticleDataStack::MPIType,
-    history::HistoryEventDataInterface, TStackIter>;
-
-using TestStack = corsika::stack::CombinedStack<
-    typename stack::nuclear_extension::ParticleDataStack::StackImpl,
-    history::HistoryEventData, StackWithHistoryInterface>;
-
-/*
-    See Issue 161
-
-    unfortunately clang does not support this in the same way (yet) as
-    gcc, so we have to distinguish here. If clang cataches up, we
-    could remove the clang branch here and also in
-    corsika::Cascade. The gcc code is much more generic and
-    universal. If we could do the gcc version, we won't had to define
-    StackView globally, we could do it with MakeView whereever it is
-    actually needed. Keep an eye on this!
-  */
-#if defined(__clang__)
-using TheTestStackView = corsika::stack::SecondaryView<typename TestStack::StackImpl,
-                                                       StackWithHistoryInterface,
-                                                       history::HistorySecondaryProducer>;
-#elif defined(__GNUC__) || defined(__GNUG__)
-using TheTestStackView =
-    corsika::stack::MakeView<TestStack, history::HistorySecondaryProducer>::type;
-#endif
-
-using TestStackView = TheTestStackView;
-
-template <typename Event>
-int count_generations(Event const* event) {
-  int genCounter = 0;
-  while (event) {
-    event = event->parentEvent().get();
-    genCounter++;
-  }
-
-  return genCounter;
-}
-
-TEST_CASE("HistoryStackExtension", "[stack]") {
-
-  logging::SetLevel(logging::level::debug);
-
-  geometry::CoordinateSystem& dummyCS =
-      geometry::RootCoordinateSystem::GetInstance().GetRootCoordinateSystem();
-
-  // in this test we only use one singel stack !
-  TestStack stack;
-
-  // add primary particle
-  auto p0 = stack.AddParticle(
-      std::make_tuple(particles::Code::Electron, 1.5_GeV,
-                      corsika::stack::MomentumVector(dummyCS, {1_GeV, 1_GeV, 1_GeV}),
-                      Point(dummyCS, {1 * meter, 1 * meter, 1 * meter}), 100_s));
-
-  CHECK(stack.getEntries() == 1);
-  corsika::history::EventPtr evt = p0.GetEvent();
-  CHECK(evt == nullptr);
-  CHECK(count_generations(evt.get()) == 0);
-
-  SECTION("interface test, view") {
-
-    // add secondaries, 1st generation
-    TestStackView hview0(p0);
-
-    auto const ev0 = p0.GetEvent();
-    CHECK(ev0 == nullptr);
-
-    C8LOG_DEBUG("loop VIEW");
-
-    // add 5 secondaries
-    for (int i = 0; i < 5; ++i) {
-      auto sec = hview0.AddSecondary(
-          std::make_tuple(particles::Code::Electron, 1.5_GeV,
-                          corsika::stack::MomentumVector(dummyCS, {1_GeV, 1_GeV, 1_GeV}),
-                          Point(dummyCS, {1 * meter, 1 * meter, 1 * meter}), 100_s));
-
-      CHECK(sec.GetParentEventIndex() == i);
-      CHECK(sec.GetEvent() != nullptr);
-      CHECK(sec.GetEvent()->parentEvent() == nullptr);
-      CHECK(count_generations(sec.GetEvent().get()) == 1);
-    }
-
-    // read 1st genertion particle particle
-    auto p1 = stack.GetNextParticle();
-    CHECK(count_generations(p1.GetEvent().get()) == 1);
-
-    TestStackView hview1(p1);
-
-    auto const ev1 = p1.GetEvent();
-
-    // add second generation of secondaries
-    // add 10 secondaries
-    for (int i = 0; i < 10; ++i) {
-      auto sec = hview1.AddSecondary(
-          std::make_tuple(particles::Code::Electron, 1.5_GeV,
-                          corsika::stack::MomentumVector(dummyCS, {1_GeV, 1_GeV, 1_GeV}),
-                          Point(dummyCS, {1 * meter, 1 * meter, 1 * meter}), 100_s));
-
-      CHECK(sec.GetParentEventIndex() == i);
-      CHECK(sec.GetEvent()->parentEvent() == ev1);
-      CHECK(sec.GetEvent()->parentEvent()->parentEvent() == ev0);
-
-      CHECK(count_generations(sec.GetEvent().get()) == 2);
-
-      const auto org_projectile = stack.at(sec.GetEvent()->projectileIndex());
-      CHECK(org_projectile.GetEvent() == sec.GetEvent()->parentEvent());
-    }
-
-    // read 2nd genertion particle particle
-    auto p2 = stack.GetNextParticle();
-
-    TestStackView hview2(p2);
-
-    auto const ev2 = p2.GetEvent();
-
-    // add third generation of secondaries
-    // add 15 secondaries
-    for (int i = 0; i < 15; ++i) {
-      C8LOG_TRACE("loop, view: " + std::to_string(i));
-
-      auto sec = hview2.AddSecondary(
-          std::make_tuple(particles::Code::Electron, 1.5_GeV,
-                          corsika::stack::MomentumVector(dummyCS, {1_GeV, 1_GeV, 1_GeV}),
-                          Point(dummyCS, {1 * meter, 1 * meter, 1 * meter}), 100_s));
-      C8LOG_TRACE("loop, ---- ");
-
-      CHECK(sec.GetParentEventIndex() == i);
-      CHECK(sec.GetEvent()->parentEvent() == ev2);
-      CHECK(sec.GetEvent()->parentEvent()->parentEvent() == ev1);
-      CHECK(sec.GetEvent()->parentEvent()->parentEvent()->parentEvent() == ev0);
-
-      CHECK(count_generations(sec.GetEvent().get()) == 3);
-    }
-  }
-
-  SECTION("also test projectile access") {
-
-    C8LOG_TRACE("projectile test");
-
-    // add secondaries, 1st generation
-    TestStackView hview0(p0);
-    auto proj0 = hview0.GetProjectile();
-    auto const ev0 = p0.GetEvent();
-    CHECK(ev0 == nullptr);
-
-    C8LOG_TRACE("loop");
-
-    // add 5 secondaries
-    for (int i = 0; i < 5; ++i) {
-      C8LOG_TRACE("loop " + std::to_string(i));
-      auto sec = proj0.AddSecondary(
-          std::make_tuple(particles::Code::Electron, 1.5_GeV,
-                          corsika::stack::MomentumVector(dummyCS, {1_GeV, 1_GeV, 1_GeV}),
-                          Point(dummyCS, {1 * meter, 1 * meter, 1 * meter}), 100_s));
-
-      CHECK(sec.GetParentEventIndex() == i);
-      CHECK(sec.GetEvent() != nullptr);
-      CHECK(sec.GetEvent()->parentEvent() == nullptr);
-    }
-    CHECK(stack.getEntries() == 6);
-
-    // read 1st genertion particle particle
-    auto p1 = stack.GetNextParticle();
-
-    TestStackView hview1(p1);
-    auto proj1 = hview1.GetProjectile();
-    auto const ev1 = p1.GetEvent();
-
-    // add second generation of secondaries
-    // add 10 secondaries
-    for (int i = 0; i < 10; ++i) {
-      auto sec = proj1.AddSecondary(
-          std::make_tuple(particles::Code::Electron, 1.5_GeV,
-                          corsika::stack::MomentumVector(dummyCS, {1_GeV, 1_GeV, 1_GeV}),
-                          Point(dummyCS, {1 * meter, 1 * meter, 1 * meter}), 100_s));
-
-      CHECK(sec.GetParentEventIndex() == i);
-      CHECK(sec.GetEvent()->parentEvent() == ev1);
-      CHECK(sec.GetEvent()->secondaries().size() == i + 1);
-      CHECK(sec.GetEvent()->parentEvent()->parentEvent() == ev0);
-
-      const auto org_projectile = stack.at(sec.GetEvent()->projectileIndex());
-      CHECK(org_projectile.GetEvent() == sec.GetEvent()->parentEvent());
-    }
-    CHECK(stack.getEntries() == 16);
-  }
-}

Stack/NuclearStackExtension/CMakeLists.txt

-set (NuclearStackExtension_HEADERS NuclearStackExtension.h)
-set (NuclearStackExtension_NAMESPACE corsika/stack/nuclear_extension)
-
-add_library (NuclearStackExtension INTERFACE)
-
-CORSIKA_COPY_HEADERS_TO_NAMESPACE (NuclearStackExtension ${NuclearStackExtension_NAMESPACE} ${NuclearStackExtension_HEADERS})
-
-target_link_libraries (
-  NuclearStackExtension
-  INTERFACE
-  CORSIKAstackinterface
-  CORSIKAunits
-  CORSIKAparticles
-  CORSIKAgeometry
-  SuperStupidStack
-  )
-
-target_include_directories (
-  NuclearStackExtension
-  INTERFACE
-  $<BUILD_INTERFACE:${PROJECT_BINARY_DIR}/include>
-  $<INSTALL_INTERFACE:include>
-  )
-
-install (
-  FILES
-  ${NuclearStackExtension_HEADERS}
-  DESTINATION
-  include/${NuclearStackExtension_NAMESPACE}
-  )
-
-# ----------------
-# code unit testing
-CORSIKA_ADD_TEST(testNuclearStackExtension)
-target_link_libraries (
-  testNuclearStackExtension
-  NuclearStackExtension
-  CORSIKAparticles
-  CORSIKAgeometry
-  CORSIKAunits
-  CORSIKAtesting
-  )

Stack/NuclearStackExtension/NuclearStackExtension.h

-/*
- * (c) Copyright 2018 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.
- */
-
-#pragma once
-
-#include <corsika/particles/ParticleProperties.h>
-
-#include <corsika/stack/Stack.h>
-#include <corsika/stack/super_stupid/SuperStupidStack.h>
-
-#include <corsika/units/PhysicalUnits.h>
-
-#include <corsika/geometry/Point.h>
-#include <corsika/geometry/Vector.h>
-
-#include <corsika/logging/Logging.h>
-
-#include <algorithm>
-#include <tuple>
-#include <vector>
-
-namespace corsika::stack {
-
-  /**
-   * @namespace nuclear_extension
-   *
-   * Add A and Z data to existing stack (currently SuperStupidStack) of particle
-   * properties. This is done via inheritance, not via CombinedStack since the nuclear
-   * data is stored ONLY when needed (for nuclei) and not for all particles. Thus, this is
-   * a new, derived Stack object.
-   *
-   * Only for Code::Nucleus particles A and Z are stored, not for all
-   * normal elementary particles.
-   *
-   * Thus in your code, make sure to always check <code>
-   * particle.GetPID()==Code::Nucleus </code> before attempting to
-   * read any nuclear information.
-   *
-   *
-   */
-
-  typedef corsika::geometry::Vector<corsika::units::si::hepmomentum_d> MomentumVector;
-
-  namespace nuclear_extension {
-
-    /**
-     * @class NuclearParticleInterface
-     *
-     * Define ParticleInterface for NuclearStackExtension Stack derived from
-     * ParticleInterface of Inner stack class
-     */
-    template <template <typename> typename InnerParticleInterface,
-              typename StackIteratorInterface>
-    class NuclearParticleInterface
-        : public InnerParticleInterface<StackIteratorInterface> {
-
-    protected:
-      using InnerParticleInterface<StackIteratorInterface>::GetStackData;
-      using InnerParticleInterface<StackIteratorInterface>::GetIndex;
-      using InnerParticleInterface<StackIteratorInterface>::as_string;
-
-    public:
-      void SetParticleData(
-          const std::tuple<corsika::particles::Code, corsika::units::si::HEPEnergyType,
-                           corsika::stack::MomentumVector, corsika::geometry::Point,
-                           corsika::units::si::TimeType>& v) {
-        if (std::get<0>(v) == corsika::particles::Code::Nucleus) {
-          std::ostringstream err;
-          err << "NuclearStackExtension: no A and Z specified for new Nucleus!";
-          throw std::runtime_error(err.str());
-        }
-        InnerParticleInterface<StackIteratorInterface>::SetParticleData(v);
-        SetNucleusRef(-1); // this is not a nucleus
-      }
-
-      void SetParticleData(
-          const std::tuple<corsika::particles::Code, corsika::units::si::HEPEnergyType,
-                           corsika::stack::MomentumVector, corsika::geometry::Point,
-                           corsika::units::si::TimeType, unsigned short, unsigned short>&
-              v) {
-        const unsigned short A = std::get<5>(v);
-        const unsigned short Z = std::get<6>(v);
-        if (std::get<0>(v) != corsika::particles::Code::Nucleus || A == 0 || Z == 0) {
-          std::ostringstream err;
-          err << "NuclearStackExtension: no A and Z specified for new Nucleus!";
-          throw std::runtime_error(err.str());
-        }
-        SetNucleusRef(GetStackData().GetNucleusNextRef()); // store this nucleus data ref
-        SetNuclearA(A);
-        SetNuclearZ(Z);
-        InnerParticleInterface<StackIteratorInterface>::SetParticleData(
-            std::tuple<corsika::particles::Code, corsika::units::si::HEPEnergyType,
-                       corsika::stack::MomentumVector, corsika::geometry::Point,
-                       corsika::units::si::TimeType>{std::get<0>(v), std::get<1>(v),
-                                                     std::get<2>(v), std::get<3>(v),
-                                                     std::get<4>(v)});
-      }
-
-      void SetParticleData(
-          InnerParticleInterface<StackIteratorInterface>& p,
-          const std::tuple<corsika::particles::Code, corsika::units::si::HEPEnergyType,
-                           corsika::stack::MomentumVector, corsika::geometry::Point,
-                           corsika::units::si::TimeType>& v) {
-        if (std::get<0>(v) == corsika::particles::Code::Nucleus) {
-          std::ostringstream err;
-          err << "NuclearStackExtension: no A and Z specified for new Nucleus!";
-          throw std::runtime_error(err.str());
-        }
-        InnerParticleInterface<StackIteratorInterface>::SetParticleData(
-            p, std::tuple<corsika::particles::Code, corsika::units::si::HEPEnergyType,
-                          corsika::stack::MomentumVector, corsika::geometry::Point,
-                          corsika::units::si::TimeType>{std::get<0>(v), std::get<1>(v),
-                                                        std::get<2>(v), std::get<3>(v),
-                                                        std::get<4>(v)});
-        SetNucleusRef(-1); // this is not a nucleus
-      }
-
-      void SetParticleData(
-          InnerParticleInterface<StackIteratorInterface>& p,
-          const std::tuple<corsika::particles::Code, corsika::units::si::HEPEnergyType,
-                           corsika::stack::MomentumVector, corsika::geometry::Point,
-                           corsika::units::si::TimeType, unsigned short, unsigned short>&
-              v) {
-        const unsigned short A = std::get<5>(v);
-        const unsigned short Z = std::get<6>(v);
-        if (std::get<0>(v) != corsika::particles::Code::Nucleus || A == 0 || Z == 0) {
-          std::ostringstream err;
-          err << "NuclearStackExtension: no A and Z specified for new Nucleus!";
-          throw std::runtime_error(err.str());
-        }
-        SetNucleusRef(GetStackData().GetNucleusNextRef()); // store this nucleus data ref
-        SetNuclearA(A);
-        SetNuclearZ(Z);
-        InnerParticleInterface<StackIteratorInterface>::SetParticleData(
-            p, std::tuple<corsika::particles::Code, corsika::units::si::HEPEnergyType,
-                          corsika::stack::MomentumVector, corsika::geometry::Point,
-                          corsika::units::si::TimeType>{std::get<0>(v), std::get<1>(v),
-                                                        std::get<2>(v), std::get<3>(v),
-                                                        std::get<4>(v)});
-      }
-
-      std::string as_string() const {
-        return fmt::format(
-            "{}, nuc({})", InnerParticleInterface<StackIteratorInterface>::as_string(),
-            (isNucleus() ? fmt::format("A={}, Z={}", GetNuclearA(), GetNuclearZ())
-                         : "n/a"));
-      }
-
-      /**
-       * @name individual setters
-       * @{
-       */
-      void SetNuclearA(const unsigned short vA) {
-        GetStackData().SetNuclearA(GetIndex(), vA);
-      }
-      void SetNuclearZ(const unsigned short vZ) {
-        GetStackData().SetNuclearZ(GetIndex(), vZ);
-      }
-      /// @}
-
-      /**
-       * @name individual getters
-       * @{
-       */
-      int GetNuclearA() const { return GetStackData().GetNuclearA(GetIndex()); }
-      int GetNuclearZ() const { return GetStackData().GetNuclearZ(GetIndex()); }
-      /// @}
-
-      /**
-       * Overwrite normal GetParticleMass function with nuclear version
-       */
-      corsika::units::si::HEPMassType GetMass() const {
-        if (InnerParticleInterface<StackIteratorInterface>::GetPID() ==
-            corsika::particles::Code::Nucleus)
-          return corsika::particles::GetNucleusMass(GetNuclearA(), GetNuclearZ());
-        return InnerParticleInterface<StackIteratorInterface>::GetMass();
-      }
-      /**
-       * Overwirte normal GetChargeNumber function with nuclear version
-       **/
-      int16_t GetChargeNumber() const {
-        if (InnerParticleInterface<StackIteratorInterface>::GetPID() ==
-            corsika::particles::Code::Nucleus)
-          return GetNuclearZ();
-        return InnerParticleInterface<StackIteratorInterface>::GetChargeNumber();
-      }
-
-      int GetNucleusRef() const { return GetStackData().GetNucleusRef(GetIndex()); }
-
-    protected:
-      void SetNucleusRef(const int vR) { GetStackData().SetNucleusRef(GetIndex(), vR); }
-      bool isNucleus() const { return GetStackData().isNucleus(GetIndex()); }
-    };
-
-    /**
-     * @class NuclearStackExtension
-     *
-     * Memory implementation of adding nuclear inforamtion to the
-     * existing particle stack defined in class InnerStackImpl.
-     *
-     * Inside the NuclearStackExtension class there is a dedicated
-     * fNucleusRef index, where fNucleusRef[i] is referring to the
-     * correct A and Z for a specific particle index i. fNucleusRef[i]
-     * == -1 means that this is not a nucleus, and a subsequent call to
-     * GetNucleusA would produce an exception.
-     */
-    template <typename InnerStackImpl>
-    class NuclearStackExtensionImpl : public InnerStackImpl {
-
-    public:
-      void Init() { InnerStackImpl::Init(); }
-      void Dump() { InnerStackImpl::Dump(); }
-
-      void Clear() {
-        InnerStackImpl::Clear();
-        fNucleusRef.clear();
-        fNuclearA.clear();
-        fNuclearZ.clear();
-      }
-
-      unsigned int GetSize() const { return fNucleusRef.size(); }
-      unsigned int GetCapacity() const { return fNucleusRef.capacity(); }
-
-      void SetNuclearA(const unsigned int i, const unsigned short vA) {
-        fNuclearA[GetNucleusRef(i)] = vA;
-      }
-      void SetNuclearZ(const unsigned int i, const unsigned short vZ) {
-        fNuclearZ[GetNucleusRef(i)] = vZ;
-      }
-      void SetNucleusRef(const unsigned int i, const int v) { fNucleusRef[i] = v; }
-
-      int GetNuclearA(const unsigned int i) const { return fNuclearA[GetNucleusRef(i)]; }
-      int GetNuclearZ(const unsigned int i) const { return fNuclearZ[GetNucleusRef(i)]; }
-      // this function will create new storage for Nuclear Properties, and return the
-      // reference to it
-      int GetNucleusNextRef() {
-        fNuclearA.push_back(0);
-        fNuclearZ.push_back(0);
-        return fNuclearA.size() - 1;
-      }
-
-      int GetNucleusRef(const unsigned int i) const {
-        if (fNucleusRef[i] >= 0) return fNucleusRef[i];
-        std::ostringstream err;
-        err << "NuclearStackExtension: no nucleus at ref=" << i;
-        throw std::runtime_error(err.str());
-      }
-
-      bool isNucleus(const unsigned int i) const { return fNucleusRef[i] >= 0; }
-
-      /**
-       *   Function to copy particle at location i1 in stack to i2
-       */
-      void Copy(const unsigned int i1, const unsigned int i2) {
-        // index range check
-        if (i1 >= GetSize() || i2 >= GetSize()) {
-          std::ostringstream err;
-          err << "NuclearStackExtension: trying to access data beyond size of stack!";
-          throw std::runtime_error(err.str());
-        }
-        // copy internal particle data p[i2] = p[i1]
-        InnerStackImpl::Copy(i1, i2);
-        // check if any of p[i1] or p[i2] was a Code::Nucleus
-        const int ref1 = fNucleusRef[i1];
-        const int ref2 = fNucleusRef[i2];
-        if (ref2 < 0) {
-          if (ref1 >= 0) {
-            // i1 is nucleus, i2 is not
-            fNucleusRef[i2] = GetNucleusNextRef();
-            fNuclearA[fNucleusRef[i2]] = fNuclearA[ref1];
-            fNuclearZ[fNucleusRef[i2]] = fNuclearZ[ref1];
-          } else {
-            // neither i1 nor i2 are nuclei
-          }
-        } else {
-          if (ref1 >= 0) {
-            // both are nuclei, i2 is overwritten with nucleus i1
-            // fNucleusRef stays the same, but A and Z data is overwritten
-            fNuclearA[ref2] = fNuclearA[ref1];
-            fNuclearZ[ref2] = fNuclearZ[ref1];
-          } else {
-            // i2 is overwritten with non-nucleus i1
-            fNucleusRef[i2] = -1;                       // flag as non-nucleus
-            fNuclearA.erase(fNuclearA.cbegin() + ref2); // remove data for i2
-            fNuclearZ.erase(fNuclearZ.cbegin() + ref2); // remove data for i2
-            const int n = fNucleusRef.size(); // update fNucleusRef: indices above ref2
-                                              // must be decremented by 1
-            for (int i = 0; i < n; ++i) {
-              if (fNucleusRef[i] > ref2) { fNucleusRef[i] -= 1; }
-            }
-          }
-        }
-      }
-
-      /**
-       *   Function to copy particle at location i2 in stack to i1
-       */
-      void Swap(const unsigned int i1, const unsigned int i2) {
-        // index range check
-        if (i1 >= GetSize() || i2 >= GetSize()) {
-          std::ostringstream err;
-          err << "NuclearStackExtension: trying to access data beyond size of stack!";
-          throw std::runtime_error(err.str());
-        }
-        // swap original particle data
-        InnerStackImpl::Swap(i1, i2);
-        // swap corresponding nuclear reference data
-        std::swap(fNucleusRef[i2], fNucleusRef[i1]);
-      }
-
-      void IncrementSize() {
-        InnerStackImpl::IncrementSize();
-        fNucleusRef.push_back(-1);
-      }
-
-      void DecrementSize() {
-        InnerStackImpl::DecrementSize();
-        if (fNucleusRef.size() > 0) {
-          const int ref = fNucleusRef.back();
-          fNucleusRef.pop_back();
-          if (ref >= 0) {
-            fNuclearA.erase(fNuclearA.begin() + ref);
-            fNuclearZ.erase(fNuclearZ.begin() + ref);
-            const int n = fNucleusRef.size();
-            for (int i = 0; i < n; ++i) {
-              if (fNucleusRef[i] >= ref) { fNucleusRef[i] -= 1; }
-            }
-          }
-        }
-      }
-
-    private:
-      /// the actual memory to store particle data
-
-      std::vector<int> fNucleusRef;
-      std::vector<unsigned short> fNuclearA;
-      std::vector<unsigned short> fNuclearZ;
-
-    }; // end class NuclearStackExtensionImpl
-
-    template <typename InnerStack, template <typename> typename _PI>
-    using NuclearStackExtension =
-        Stack<NuclearStackExtensionImpl<typename InnerStack::StackImpl>, _PI>;
-
-    //
-    template <typename StackIter>
-    using ExtendedParticleInterfaceType =
-        corsika::stack::nuclear_extension::NuclearParticleInterface<
-            corsika::stack::super_stupid::SuperStupidStack::MPIType, StackIter>;
-
-    // the particle data stack with extra nuclear information:
-    using ParticleDataStack = corsika::stack::nuclear_extension::NuclearStackExtension<
-        corsika::stack::super_stupid::SuperStupidStack, ExtendedParticleInterfaceType>;
-
-  } // namespace nuclear_extension
-} // namespace corsika::stack