corsika/detail/output/ParquetStreamer.inl

 /*
  * (c) Copyright 2021 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.
+ * This software is distributed under the terms of the 3-clause BSD license.
+ * See file LICENSE for a full version of the license.
  */
 
 #pragma once
@@ -23,7 +22,7 @@ namespace corsika {
 
     // add run and event tags to the file
     addField("shower", parquet::Repetition::REQUIRED, parquet::Type::INT32,
-             parquet::ConvertedType::INT_32);
+             parquet::ConvertedType::UINT_32);
   }
 
   template <typename... TArgs>
@@ -31,9 +30,9 @@ namespace corsika {
     fields_.push_back(parquet::schema::PrimitiveNode::Make(args...));
   }
 
-  inline void ParquetStreamer::enableCompression(int const /*level*/) {
-    // builder_.compression(parquet::Compression::ZSTD);
-    // builder_.compression_level(level);
+  inline void ParquetStreamer::enableCompression(int const level) {
+    builder_.compression(parquet::Compression::LZ4);
+    builder_.compression_level(level);
   }
 
   inline void ParquetStreamer::buildStreamer() {

corsika/detail/output/YAMLStreamer.inl

+/*
+ * (c) Copyright 2021 CORSIKA Project, corsika-project@lists.kit.edu
+ *
+ * This software is distributed under the terms of the 3-clause BSD license.
+ * See file LICENSE for a full version of the license.
+ */
+
+#pragma once
+
+#include <boost/filesystem/fstream.hpp>
+
+namespace corsika {
+
+  inline void YAMLStreamer::writeYAML(YAML::Node const& node,
+                                      boost::filesystem::path const& path) const {
+
+    // construct a YAML emitter for this config file
+    YAML::Emitter out;
+
+    // and write the node to the output
+    out << node;
+
+    // open the output file - this is <output name>.yaml
+    boost::filesystem::ofstream file(path);
+
+    // dump the YAML to the file
+    file << out.c_str() << std::endl;
+  }
+
+} // namespace corsika

corsika/detail/setup/SetupEnvironment.inl

-/*
- * (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.
- */
-
-#pragma once
-
-#include <corsika/framework/geometry/Point.hpp>
-#include <corsika/framework/geometry/CoordinateSystem.hpp>
-
-#include <limits>

corsika/detail/setup/SetupStack.hpp

 /*
  * (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.
+ * This software is distributed under the terms of the 3-clause BSD license.
+ * See file LICENSE for a full version of the license.
  */
 
 #pragma once
@@ -14,46 +13,84 @@
 #include <corsika/stack/WeightStackExtension.hpp>
 #include <corsika/stack/history/HistorySecondaryProducer.hpp>
 #include <corsika/stack/history/HistoryStackExtension.hpp>
-
-#include <corsika/setup/SetupEnvironment.hpp>
+#include <corsika/media/Environment.hpp>
+#include <corsika/media/IMagneticFieldModel.hpp>
+#include <corsika/media/IMediumModel.hpp>
+#include <corsika/media/IMediumPropertyModel.hpp>
 
 namespace corsika {
 
   namespace setup::detail {
+    template <typename TEnvironment>
+    class StackGenerator {
+    private:
+      using env_type = TEnvironment;
+
+      // ------------------------------------------
+      // add geometry node data to stack. This is fundamentally needed
+      // for robust tracking through multiple volumes.
+
+      // the GeometryNode stack needs to know the type of geometry-nodes from the
+      // environment:
+
+      template <typename TStackIter>
+      using SetupGeometryDataInterface =
+          typename node::MakeGeometryDataInterface<TStackIter, env_type>::type;
+
+      // combine particle data stack with geometry information for tracking
+      template <typename TStackIter>
+      using StackWithGeometryInterface =
+          CombinedParticleInterface<VectorStack::pi_type, SetupGeometryDataInterface,
+                                    TStackIter>;
+
+      using StackWithGeometry =
+          CombinedStack<typename VectorStack::stack_data_type,
+                        node::GeometryData<env_type>, StackWithGeometryInterface,
+                        DefaultSecondaryProducer>;
+
+      template <class T>
+      using StackWithGeometry_PI_type = typename StackWithGeometry::template pi_type<T>;
+
+      // ------------------------------------------
+      // add weight data to stack. This is fundamentally needed
+      // for thinning.
+
+      // the "pure" weight stack (interface)
+      template <typename TStackIter>
+      using SetupWeightDataInterface =
+          typename weights::MakeWeightDataInterface<TStackIter>::type;
+
+      // combine geometry-node-vector data stack with weight information for tracking
+      template <typename TStackIter>
+      using StackWithWeightInterface =
+          CombinedParticleInterface<StackWithGeometry_PI_type, SetupWeightDataInterface,
+                                    TStackIter>;
+
+    public:
+      // the combined stack data: particle + geometry + weight
+      using StackWithWeight =
+          CombinedStack<typename StackWithGeometry::stack_data_type, weights::WeightData,
+                        StackWithWeightInterface, DefaultSecondaryProducer>;
+
+    private:
+      template <typename T>
+      using StackWithWeight_PI_type = typename StackWithWeight::template pi_type<T>;
+
+      // ------------------------------------------
+      // Add [OPTIONAL] history data to stack, too.
+      // This keeps the entire lineage of particles in memory.
+
+      template <typename TStackIter>
+      using StackWithHistoryInterface =
+          CombinedParticleInterface<StackWithWeight_PI_type,
+                                    history::HistoryEventDataInterface, TStackIter>;
 
-    // ------------------------------------------
-    // 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 node::MakeGeometryDataInterface<TStackIter, setup::Environment>::type;
-
-    // combine particle data stack with geometry information for tracking
-    template <typename TStackIter>
-    using StackWithGeometryInterface =
-        CombinedParticleInterface<VectorStack::pi_type, SetupGeometryDataInterface,
-                                  TStackIter>;
-
-    using StackWithGeometry =
-        CombinedStack<typename VectorStack::stack_data_type,
-                      node::GeometryData<setup::Environment>, StackWithGeometryInterface,
-                      DefaultSecondaryProducer>;
-
-    // ------------------------------------------
-    // Add [optional] history data to stack, too:
-
-    // combine dummy stack with geometry information for tracking
-    template <typename TStackIter>
-    using StackWithHistoryInterface =
-        CombinedParticleInterface<StackWithGeometry::pi_type,
-                                  history::HistoryEventDataInterface, TStackIter>;
-
-    using StackWithHistory =
-        CombinedStack<typename StackWithGeometry::stack_data_type,
-                      history::HistoryEventData, StackWithHistoryInterface,
-                      history::HistorySecondaryProducer>;
+    public:
+      using StackWithHistory =
+          CombinedStack<typename StackWithWeight::stack_data_type,
+                        history::HistoryEventData, StackWithHistoryInterface,
+                        history::HistorySecondaryProducer>;
+    };
 
   } // namespace setup::detail
 

corsika/detail/setup/SetupStack.inl

-/*
- * (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.
- */

corsika/detail/stack/VectorStack.inl

 /*
  * (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.
+ * This software is distributed under the terms of the 3-clause BSD license.
+ * See file LICENSE for a full version of the license.
  */
 
 #pragma once

corsika/detail/stack/WeightStackExtension.inl

 /*
  * (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.
+ * This software is distributed under the terms of the 3-clause BSD license.
+ * See file LICENSE for a full version of the license.
  */
 
 #pragma once

corsika/detail/stack/history/HistoryObservationPlane.hpp

 /*
  * (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.
+ * This software is distributed under the terms of the 3-clause BSD license.
+ * See file LICENSE for a full version of the license.
  */
 
 #pragma once

corsika/detail/stack/history/HistoryObservationPlane.inl

 /*
  * (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.
+ * This software is distributed under the terms of the 3-clause BSD license.
+ * See file LICENSE for a full version of the license.
  */
 
 #include <corsika/logging/Logging.h>

corsika/framework/analytics/ClassTimer.hpp

 /*
  * (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.
+ * This software is distributed under the terms of the 3-clause BSD license.
+ * See file LICENSE for a full version of the license.
  */
 
 // Another possibility:

corsika/framework/analytics/FunctionTimer.hpp

 /*
  * (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.
+ * This software is distributed under the terms of the 3-clause BSD license.
+ * See file LICENSE for a full version of the license.
  */
 #pragma once
 

corsika/framework/analytics/Timer.hpp

 /*
  * (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.
+ * This software is distributed under the terms of the 3-clause BSD license.
+ * See file LICENSE for a full version of the license.
  */
 
 #pragma once

corsika/framework/core/Cascade.hpp

 /*
  * (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.
+ * This software is distributed under the terms of the 3-clause BSD license.
+ * See file LICENSE for a full version of the license.
  */
 
 #pragma once
@@ -55,7 +54,7 @@ namespace corsika {
 
     typedef typename TStack::particle_type particle_type;
 
-    typedef std::remove_pointer_t<decltype(((particle_type*)nullptr)->getNode())>
+    typedef std::remove_pointer_t<decltype(std::declval<particle_type>().getNode())>
         volume_tree_node_type;
 
     typedef typename volume_tree_node_type::IModelProperties medium_interface_type;
@@ -73,21 +72,7 @@ namespace corsika {
     ~Cascade() = default;
     Cascade& operator=(Cascade const&) = default;
     Cascade(Environment<medium_interface_type> const& env, TTracking& tr,
-            TProcessList& pl, TOutput& out, TStack& stack)
-        : environment_(env)
-        , tracking_(tr)
-        , sequence_(pl)
-        , output_(out)
-        , stack_(stack) {
-      CORSIKA_LOG_INFO(c8_ascii_);
-      CORSIKA_LOG_INFO("This is CORSIKA {}.{}.{}.{}", CORSIKA_RELEASE_NUMBER,
-                       CORSIKA_MAJOR_NUMBER, CORSIKA_MINOR_NUMBER, CORSIKA_PATCH_NUMBER);
-      CORSIKA_LOG_INFO("Tracking algorithm: {} (version {})", TTracking::getName(),
-                       TTracking::getVersion());
-      if constexpr (stack_view_type::has_event) {
-        CORSIKA_LOG_INFO("Stack - with full cascade HISTORY.");
-      }
-    }
+            TProcessList& pl, TOutput& out, TStack& stack);
     //! @}
 
     /**
@@ -106,9 +91,18 @@ namespace corsika {
      * Force an interaction of the top particle of the stack at its current position.
      * Note that setNodes() or an equivalent procedure needs to be called first if you
      * want to call forceInteraction() for the primary interaction.
+     * Incompatible with forceDecay()
      */
     void forceInteraction();
 
+    /**
+     * Force an decay of the top particle of the stack at its current position.
+     * Note that setNodes() or an equivalent procedure needs to be called first if you
+     * want to call forceDecay() for the primary interaction.
+     * Incompatible with forceInteraction()
+     */
+    void forceDecay();
+
   private:
     /**
      * The Step function is executed for each particle from the
@@ -135,6 +129,8 @@ namespace corsika {
     TOutput& output_;
     TStack& stack_;
     default_prng_type& rng_ = RNGManager<>::getInstance().getRandomStream("cascade");
+    bool forceInteraction_;
+    bool forceDecay_;
     unsigned int count_ = 0;
 
     // but this here temporarily. Should go into dedicated file later:

corsika/framework/core/EnergyMomentumOperations.hpp

 /*
  * (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.
+ * This software is distributed under the terms of the 3-clause BSD license.
+ * See file LICENSE for a full version of the license.
  */
 
 #pragma once
@@ -117,4 +116,17 @@ namespace corsika {
     return (sqrtS_sqr - static_pow<2>(m_proj) - static_pow<2>(m_targ)) / (2 * m_targ);
   }
 
+  /**
+   * \f[E_{com}=sqrt{2 * m_{proj} * m_{targ} * E_{lab} + m_{proj}^2 + m_{targ}^2} \f]
+   *
+   * @param E lab. energy.
+   * @param m particle mass.
+   * @return HEPEnergyType
+   */
+  HEPEnergyType constexpr calculate_com_energy(HEPEnergyType Elab,
+                                               HEPMassType const m_proj,
+                                               HEPMassType const m_targ) {
+    return sqrt(2 * Elab * m_targ + static_pow<2>(m_proj) + static_pow<2>(m_targ));
+  }
+
 } // namespace corsika

corsika/framework/core/Logging.hpp

 /*
  * (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.
+ * This software is distributed under the terms of the 3-clause BSD license.
+ * See file LICENSE for a full version of the license.
  */
 
 /**
@@ -22,23 +21,31 @@
 
 // use the coarse system clock. This is *much* faster
 // but introduces a timestamp error of O(10 ms) which is fine for us.
+#ifndef SPDLOG_CLOCK_COARSE
 #define SPDLOG_CLOCK_COARSE
-
+#endif
 // do not create a default logger (we provide our own "corsika" logger)
+#ifndef SPDLOG_DISABLE_DEFAULT_LOGGER
 #define SPDLOG_DISABLE_DEFAULT_LOGGER
-
+#endif
 // use __PRETTY_FUNCTION__ instead of __FUNCTION__ where
 // printing function names in trace statements. This is much
 // nicer than __FUNCTION__ under GCC/clang.
+#ifndef SPDLOG_FUNCTION
 #define SPDLOG_FUNCTION __PRETTY_FUNCTION__
+#endif
 
 // if this is a Debug build, include debug messages in objects
-#ifdef DEBUG
+#ifdef _C8_DEBUG_
 // trace is the highest level of logging (ALL messages will be printed)
+#ifndef SPDLOG_ACTIVE_LEVEL
 #define SPDLOG_ACTIVE_LEVEL SPDLOG_LEVEL_TRACE
+#endif
 #else // otherwise, remove everything but "error" and worse messages
+#ifndef SPDLOG_ACTIVE_LEVEL
 #define SPDLOG_ACTIVE_LEVEL SPDLOG_LEVEL_DEBUG
 #endif
+#endif
 
 #include <spdlog/fmt/ostr.h> // will output whenerver a streaming operator is found
 #include <spdlog/sinks/stdout_color_sinks.h>
@@ -149,3 +156,4 @@ namespace corsika {
 } // namespace corsika
 
 #include <corsika/detail/framework/core/Logging.inl>
+#include <corsika/detail/framework/core/SpdlogSpecializations.inl>

corsika/framework/core/ParticleProperties.hpp

 /*
  * (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.
+ * This software is distributed under the terms of the 3-clause BSD license.
+ * See file LICENSE for a full version of the license.
  */
 
 /**
@@ -31,7 +30,7 @@ namespace corsika {
    *
    * The properties of all particles are saved in static and flat
    * arrays. There is a enum corsika::Code to identify each
-   * particles, and each individual particles has its own static class,
+   * particle, and each individual particle has its own static class,
    * which can be used to retrieve its physical properties.
    *
    * The properties of all elementary particles are accessible here. The data
@@ -54,6 +53,25 @@ namespace corsika {
    * The names, relations and properties of all particles known to CORSIKA 8 are listed
    * below.
    *
+   * **Note** on energy threshold on particle production as well as particle propagation.
+   * The functions:
+   * @code {.cpp}
+   * HEPEnergyType constexpr get_energy_production_threshold(Code const);
+   * void constexpr set_energy_production_threshold(Code const, HEPEnergyType const);
+   * @endcode
+   * can be used to tune the transition where explicit production of new particles, e.g.
+   * in Bremsstrahlung, is simulated versus a continuous handling of low-energy particles
+   * as generic energy losses. The default value for all particle types is 1 MeV.
+   *
+   * Furthermore, the functions:
+   * @code {.cpp}
+   * HEPEnergyType constexpr get_kinetic_energy_propagation_threshold(Code const);
+   * void constexpr set_kinetic_energy_propagation_threshold(Code const, HEPEnergyType
+   *                                                         const);
+   * @endcode
+   * are used to discard low energy particle during tracking. The default value for all
+   * particle types is 1 GeV.
+   *
    * @addtogroup Particles
    * @{
    */
@@ -63,15 +81,14 @@ namespace corsika {
    *
    * The Code enum is the actual place to define CORSIKA 8 particle codes.
    */
-  enum class Code : int32_t;
+  enum class Code : std::int32_t;
 
   /**
    * @enum PDGCode
    *
    * Specifically for PDG ids.
    */
-  enum class PDGCode : int32_t;
-
+  enum class PDGCode : std::int32_t;
   /**
    * Internal integer type for enum Code.
    */
@@ -89,34 +106,53 @@ namespace corsika {
 namespace corsika {
 
   // forward declarations to be used in GeneratedParticleProperties
+  struct full_name {}; //!< tag class for get_name()
+
   int16_t constexpr get_charge_number(Code const);     //!< electric charge in units of e
   ElectricChargeType constexpr get_charge(Code const); //!< electric charge
   HEPMassType constexpr get_mass(Code const);          //!< mass
-  HEPEnergyType constexpr calculate_kinetic_energy_threshold(
-      Code const); //!< get kinetic energy threshold below which the particle is
-                   //!< discarded, by default set to zero
-  void constexpr set_kinetic_energy_threshold(
-      Code const, HEPEnergyType const); //!< set kinetic energy threshold below which the
-                                        //!< particle is discarded
-
-  inline void set_kinetic_energy_threshold(std::pair<Code const, HEPEnergyType const> p) {
-    set_kinetic_energy_threshold(p.first, p.second);
-  }
-  inline void set_kinetic_energy_thresholds(
-      std::unordered_map<Code const, HEPEnergyType const> const& eCuts) {
-    for (auto v : eCuts) set_kinetic_energy_threshold(v);
-  }
+
+  /**
+   * Get the kinetic energy propagation threshold.
+   *
+   * Particles are tracked only above the kinetic energy propagation threshold. Below
+   * this, they are discarded and removed. Sensible default values must be configured for
+   * a simulation.
+   */
+  HEPEnergyType get_kinetic_energy_propagation_threshold(Code const);
+
+  /**
+   * Set the kinetic energy propagation threshold object.
+   */
+  void set_kinetic_energy_propagation_threshold(Code const, HEPEnergyType const);
+
+  /**
+   * Get the particle production energy threshold.
+   *
+   * The (total) energy below which a particle is only handled stoachastically (no
+   * production below this energy). This is for example important for stochastic discrete
+   * Bremsstrahlung versus low-energy Bremsstrahlung as part of continuous energy losses.
+   */
+  HEPEnergyType get_energy_production_threshold(Code const); //!<
+
+  /**
+   * Set the particle production energy threshold in total energies.
+   */
+  void set_energy_production_threshold(Code const, HEPEnergyType const);
 
   //! Particle code according to PDG, "Monte Carlo Particle Numbering Scheme"
   PDGCode constexpr get_PDG(Code const);
-
+  PDGCode constexpr get_PDG(unsigned int const A, unsigned int const Z);
   std::string_view constexpr get_name(Code const); //!< name of the particle as string
-  TimeType constexpr get_lifetime(Code const);     //!< lifetime
+  std::string get_name(Code,
+                       full_name); //!< get name of particle, including (A,Z) for nuclei
+  TimeType constexpr get_lifetime(Code const); //!< lifetime
 
   bool constexpr is_hadron(Code const); //!< true if particle is hadron
   bool constexpr is_em(Code const); //!< true if particle is electron, positron or photon
   bool constexpr is_muon(Code const);     //!< true if particle is mu+ or mu-
   bool constexpr is_neutrino(Code const); //!< true if particle is (anti-) neutrino
+  bool constexpr is_charged(Code const);  //!< true if particle is charged
 
   /**
    * @brief Creates the Code for a nucleus of type 10LZZZAAAI.
@@ -169,7 +205,7 @@ namespace corsika {
    * @param code
    * @return std::string_view
    */
-  inline std::string_view get_nucleus_name(Code const code);
+  inline std::string get_nucleus_name(Code const code);
 
   /**
    * @brief convert PDG code to CORSIKA 8 internal code.

corsika/framework/core/PhysicalConstants.hpp

@@ -2,9 +2,8 @@
  * (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.
+ * This software is distributed under the terms of the 3-clause BSD license.
+ * See file LICENSE for a full version of the license.
  */
 
 #pragma once
@@ -35,6 +34,10 @@ namespace corsika::constants {
   // elementary charge
   constexpr quantity<electric_charge_d> e{Rep(1.6021766208e-19L) * coulomb};
 
+  // vacuum permittivity
+  constexpr quantity<dimensions<-3, -1, 4, 2>> epsilonZero{Rep(8.8541878128e-12L) *
+                                                           farad / meter};
+
   // electronvolt
   // constexpr quantity<hepenergy_d> eV{e / coulomb * joule};
 
@@ -65,7 +68,8 @@ namespace corsika::constants {
    */
   namespace EarthRadius {
     static constexpr auto Mean{6'371'000 * meter};
-    static constexpr auto Eqautorial{6'378'137 * meter};
+    static constexpr auto Geomagnetic_reference{6'371'200 * meter};
+    static constexpr auto Equatorial{6'378'137 * meter};
     static constexpr auto Polar{6'356'752 * meter};
     static constexpr auto PolarCurvature{6'399'593 * meter};
   } // namespace EarthRadius

corsika/framework/core/PhysicalGeometry.hpp

 /*
  * (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.
+ * This software is distributed under the terms of the 3-clause BSD license.
+ * See file LICENSE for a full version of the license.
  */
 
 #pragma once
@@ -23,18 +22,34 @@ namespace corsika {
   /**
    * A 3D vector defined in a specific coordinate system with units HEPMomentumType
    **/
-  typedef Vector<hepmomentum_d> MomentumVector;
+  using MomentumVector = Vector<hepmomentum_d>;
 
   /**
    * A 3D vector defined in a specific coordinate system with no units. But, note, this is
    * not automatically normaliyed! It is not a "NormalVector".
    **/
-  typedef Vector<dimensionless_d> DirectionVector;
+  using DirectionVector = Vector<dimensionless_d>;
 
   /**
    * A 3D vector defined in a specific coordinate system with units "velocity_t".
    *
    **/
-  typedef Vector<SpeedType::dimension_type> VelocityVector;
+  using VelocityVector = Vector<SpeedType::dimension_type>;
+
+  /**
+   * A 3D vector defined in a specific coordinate system with units "length_t".
+   *
+   **/
+  using LengthVector = Vector<length_d>;
+
+  /**
+   * A 3D vector defined in a specific coordinate system with units ElectricFieldType
+   **/
+  typedef Vector<ElectricFieldType::dimension_type> ElectricFieldVector;
+
+  /**
+   * A 3D vector defined in a specific coordinate system with units VectorPotentialType
+   **/
+  typedef Vector<VectorPotentialType::dimension_type> VectorPotential;
 
 } // namespace corsika

corsika/framework/core/PhysicalUnits.hpp

 /*
  * (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.
+ * This software is distributed under the terms of the 3-clause BSD license.
+ * See file LICENSE for a full version of the license.
  */
 
 #pragma once
@@ -92,10 +91,16 @@ namespace corsika::units::si {
       phys::units::quantity<phys::units::dimensions<-1, 0, 0>, double>;
   using InverseTimeType =
       phys::units::quantity<phys::units::dimensions<0, 0, -1>, double>;
+  using InverseMassDensityType =
+      phys::units::quantity<phys::units::dimensions<3, -1, 0>, double>;
   using InverseGrammageType =
       phys::units::quantity<phys::units::dimensions<2, -1, 0>, double>;
   using MagneticFluxType =
       phys::units::quantity<phys::units::magnetic_flux_density_d, double>;
+  using ElectricFieldType =
+      phys::units::quantity<phys::units::dimensions<1, 1, -3, -1>, double>;
+  using VectorPotentialType =
+      phys::units::quantity<phys::units::dimensions<1, 1, -2, -1>, double>;
 
   template <typename DimFrom, typename DimTo>
   auto constexpr conversion_factor_HEP_to_SI() {

corsika/framework/core/Step.hpp

+/*
+ * (c) Copyright 2022 CORSIKA Project, corsika-project@lists.kit.edu
+ *
+ * This software is distributed under the terms of the 3-clause BSD license.
+ * See file LICENSE for a full version of the license.
+ */
+
+#pragma once
+
+#include <corsika/framework/core/PhysicalUnits.hpp>
+#include <corsika/framework/core/PhysicalGeometry.hpp>
+#include <corsika/framework/geometry/Vector.hpp>
+#include <corsika/framework/geometry/Point.hpp>
+#include <corsika/framework/geometry/StraightTrajectory.hpp>
+
+namespace corsika {
+
+  struct DeltaParticleState {
+  public:
+    DeltaParticleState(HEPEnergyType dEkin, TimeType dT, LengthVector const& ds,
+                       DirectionVector const& du)
+        : delta_Ekin_{dEkin}
+        , delta_time_{dT}
+        , displacement_{ds}
+        , delta_direction_{du} {}
+
+    HEPEnergyType delta_Ekin_ = HEPEnergyType::zero();
+    TimeType delta_time_ = TimeType::zero();
+    LengthVector displacement_;
+    DirectionVector delta_direction_;
+  };
+
+  template <typename TParticle>
+  class Step {
+  public:
+    template <typename TTrajectory>
+    Step(TParticle const& particle, TTrajectory const& track)
+        : particlePreStep_{particle} //~ , track_{track}
+        , diff_{HEPEnergyType::zero(), track.getDuration(1),
+                track.getPosition(1) - particle.getPosition(),
+                track.getDirection(1) - particle.getDirection()}
+
+    {}
+
+    HEPEnergyType const& add_dEkin(HEPEnergyType dEkin) {
+      diff_.delta_Ekin_ += dEkin;
+      return diff_.delta_Ekin_;
+    }
+
+    TimeType const& add_dt(TimeType dt) {
+      diff_.delta_time_ += dt;
+      return diff_.delta_time_;
+    }
+
+    LengthVector const& add_displacement(LengthVector const& dis) {
+      diff_.displacement_ += dis;
+      return diff_.displacement_;
+    }
+
+    DirectionVector const& add_dU(DirectionVector const& du) {
+      diff_.delta_direction_ += du;
+      return diff_.delta_direction_;
+    }
+
+    // getters for difference
+
+    HEPEnergyType getDiffEkin() const { return diff_.delta_Ekin_; }
+
+    TimeType getDiffT() const { return diff_.delta_time_; };
+
+    DirectionVector const& getDiffDirection() const { return diff_.delta_direction_; }
+
+    LengthVector const& getDisplacement() const { return diff_.displacement_; }
+
+    //! alias for getDisplacement()
+    LengthVector const& getDiffPosition() const { return getDisplacement(); }
+
+    // getters for absolute
+    TParticle const& getParticlePre() const { return particlePreStep_; }
+
+    HEPEnergyType getEkinPre() const { return getParticlePre().getKineticEnergy(); }
+
+    HEPEnergyType getEkinPost() const { return getEkinPre() + getDiffEkin(); }
+
+    TimeType getTimePre() const { return getParticlePre().getTime(); }
+
+    TimeType getTimePost() const { return getTimePre() + getDiffT(); }
+
+    DirectionVector const getDirectionPre() const {
+      return getParticlePre().getDirection();
+    }
+
+    VelocityVector getVelocityVector() const { return getDisplacement() / getDiffT(); }
+
+    StraightTrajectory getStraightTrack() const {
+      Line const line(getPositionPre(), getVelocityVector());
+      StraightTrajectory track(line, getDiffT());
+      return track;
+    }
+
+    DirectionVector getDirectionPost() const {
+      return (getDirectionPre() + getDiffDirection()).normalized();
+    }
+
+    Point const& getPositionPre() const { return getParticlePre().getPosition(); }
+
+    Point getPositionPost() const { return getPositionPre() + getDisplacement(); }
+
+  private:
+    TParticle const& particlePreStep_;
+    //~ TTrajectory const& track_; // TODO: perhaps remove
+    DeltaParticleState diff_;
+  };
+
+} // namespace corsika