corsika/media/FlatExponential.hpp

+/*
+ * (c) Copyright 2020 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/geometry/Line.hpp>
+#include <corsika/framework/geometry/Point.hpp>
+#include <corsika/media/BaseExponential.hpp>
+#include <corsika/media/NuclearComposition.hpp>
+#include <corsika/framework/geometry/BaseTrajectory.hpp>
+
+namespace corsika {
+
+  // clang-format off
+  /**
+   * flat exponential density distribution with
+   * \f[
+   *  \varrho(\vec{r}) = \varrho_0 \exp\left( \frac{1}{\lambda} (\vec{r} - \vec{p}) \cdot
+   *    \vec{a} \right).
+   * \f]
+   * \f$ \vec{a} \f$ denotes the axis and should be normalized to avoid degeneracy
+   * with the scale parameter \f$ \lambda \f$, \f$ \vec{r} \f$ is the location of 
+   * the evaluation, \f$ \vec{p} \f$ is the anchor point at which \f$ \varrho_0 \f$ 
+   * is given. Thus, the unit vector \f$ \vec{a} \f$ specifies the direction of
+   * <em>decreasing</em> <b>height/altitude</b>.
+   */
+  // clang-format on
+  template <typename T>
+  class FlatExponential : public BaseExponential<FlatExponential<T>>, public T {
+    using base_type = BaseExponential<FlatExponential<T>>;
+
+  public:
+    FlatExponential(Point const& point, Vector<dimensionless_d> const& axis,
+                    MassDensityType const rho, LengthType const lambda,
+                    NuclearComposition const& nuclComp);
+
+    MassDensityType getMassDensity(Point const& point) const override;
+
+    NuclearComposition const& getNuclearComposition() const override;
+
+    GrammageType getIntegratedGrammage(BaseTrajectory const& line) const override;
+
+    LengthType getArclengthFromGrammage(BaseTrajectory const& line,
+                                        GrammageType const grammage) const override;
+
+  private:
+    DirectionVector const axis_;
+    NuclearComposition const nuclComp_;
+  };
+
+} // namespace corsika
+
+#include <corsika/detail/media/FlatExponential.inl>

corsika/media/GeomagneticModel.hpp

+/*
+ * (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 <corsika/framework/core/PhysicalUnits.hpp>
+#include <corsika/framework/geometry/PhysicalGeometry.hpp>
+#include <corsika/framework/utility/CorsikaData.hpp>
+
+#include <boost/filesystem.hpp>
+#include <fstream>
+#include <string>
+
+namespace corsika {
+
+  /**
+   * A magnetic field calculated with the WMM or IGRF model.
+   * WMM: https://geomag.bgs.ac.uk/documents/WMM2020_Report.pdf
+   * IGRF: https://www.ngdc.noaa.gov/IAGA/vmod/igrf.html
+   */
+
+  class GeomagneticModel {
+
+    /**
+     * Internal data structure for a single shell of the spherical harmonic
+     * parameterization.
+     */
+    struct ParameterLine {
+      int n;
+      int m;
+      double g;
+      double h;
+      double dg; //< time dependence of g in "per year"
+      double dh; //< time dependence of h in "per year"
+    };
+
+  public:
+    /**
+     * Construct a new World Magnetic Model object.
+     *
+     * @param center Center of Earth.
+     * @param data Data table to read.
+     */
+    GeomagneticModel(Point const& center, boost::filesystem::path const path =
+                                              corsika::corsika_data("GeoMag/WMM.COF"));
+
+    /**
+     * Calculates the value of the magnetic field.
+     *
+     * @param  year        Year of the evaluation, between 2020 and 2025.
+     * @param  altitude    Height of the location to evaluate the field at,
+     *                    in km between -1 and 850.
+     * @param  latitude    Latitude of the location to evaluate the field at,
+     *                   in degrees between -90 and 90 (negative for southern hemisphere).
+     * @param  longitute   Longitude of the location to evaluate the field at,
+     *                      in degrees between -180 and 180 (negative for western
+     *   hemisphere).
+     *
+     * @returns    The magnetic field vector in nT.
+     */
+    MagneticFieldVector getField(double const year, LengthType const altitude,
+                                 double const latitude, double const longitude);
+
+  private:
+    Point center_;                                         //< Center of Earth.
+    std::map<int, std::vector<ParameterLine>> parameters_; //< epoch to Parameters map
+  };
+} // namespace corsika
+
+#include <corsika/detail/media/GeomagneticModel.inl>
\ No newline at end of file

corsika/media/GladstoneDaleRefractiveIndex.hpp

+/*
+ * (c) Copyright 2023 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/media/IRefractiveIndexModel.hpp>
+
+namespace corsika {
+
+  /**
+   * A tabulated refractive index.
+   *
+   * This class returns the value of refractive index
+   * for a specific height bin.
+   */
+  template <typename T>
+  class GladstoneDaleRefractiveIndex : public T {
+
+    double const
+        referenceRefractivity_; ///< The constant reference refractive index minus one.
+    InverseMassDensityType const
+        referenceInvDensity_; ///< The constant inverse mass density at reference point.
+
+  public:
+    /**
+     * Construct a GladstoneDaleRefractiveIndex.
+     *
+     * This is initialized with a fixed refractive index
+     * at sea level and uses the Gladstone-Dale law to
+     * calculate the refractive index at a given point.
+     *
+     * @param seaLevelRefractiveIndex  The refractive index at sea level.
+     * @param point AA point at earth's surface.
+     */
+    template <typename... Args>
+    GladstoneDaleRefractiveIndex(double const seaLevelRefractiveIndex, Point const point,
+                                 Args&&... args);
+
+    /**
+     * Evaluate the refractive index at a given location.
+     *
+     * @param  point    The location to evaluate at.
+     * @returns    The refractive index at this point.
+     */
+    double getRefractiveIndex(Point const& point) const override;
+
+  }; // END: class GladstoneDaleRefractiveIndex
+
+} // namespace corsika
+
+#include <corsika/detail/media/GladstoneDaleRefractiveIndex.inl>

corsika/media/HomogeneousMedium.hpp

+/*
+ * (c) Copyright 2020 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/geometry/Line.hpp>
+#include <corsika/framework/geometry/Point.hpp>
+#include <corsika/media/NuclearComposition.hpp>
+#include <corsika/framework/geometry/BaseTrajectory.hpp>
+
+namespace corsika {
+
+  /**
+   * a homogeneous medium
+   */
+
+  template <typename T>
+  class HomogeneousMedium : public T {
+
+  public:
+    HomogeneousMedium(MassDensityType density, NuclearComposition const& nuclComp);
+
+    MassDensityType getMassDensity(Point const&) const override;
+
+    NuclearComposition const& getNuclearComposition() const override;
+
+    GrammageType getIntegratedGrammage(BaseTrajectory const&) const override;
+
+    LengthType getArclengthFromGrammage(BaseTrajectory const&,
+                                        GrammageType grammage) const override;
+
+  private:
+    MassDensityType const density_;
+    NuclearComposition const nuclComp_;
+  };
+
+} // namespace corsika
+
+#include <corsika/detail/media/HomogeneousMedium.inl>

corsika/media/IEmpty.hpp

+/*
+ * (c) Copyright 2020 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/geometry/BaseTrajectory.hpp>
+
+namespace corsika {
+
+  /**
+   * Intended for usage as default template argument for environments
+   * with no properties.  For now, the ArclengthFromGrammage is
+   * mandatory, since it is used even in the most simple Cascade code.
+   *
+   * - IEmpty is the interface definition.
+   * - Empty<IEmpty> is a possible model implementation
+   *
+   **/
+
+  class IEmpty {
+  public:
+    virtual LengthType getArclengthFromGrammage(BaseTrajectory const&,
+                                                GrammageType) const = 0;
+
+    virtual NuclearComposition const& getNuclearComposition() const = 0;
+
+    virtual ~IEmpty() {}
+  };
+
+  template <typename TModel = IEmpty>
+  class Empty : public TModel {
+  public:
+    LengthType getArclengthFromGrammage(BaseTrajectory const&, GrammageType) const {
+      return 0. * meter;
+    }
+    NuclearComposition const& getNuclearComposition() const {
+      return NuclearComposition(std::vector<Code>{}, {});
+    };
+  };
+
+} // namespace corsika

corsika/media/IMagneticFieldModel.hpp

+/*
+ * (c) Copyright 2020 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/geometry/Point.hpp>
+#include <corsika/framework/geometry/Vector.hpp>
+#include <corsika/framework/core/PhysicalUnits.hpp>
+
+namespace corsika {
+
+  /**
+   * An interface for magnetic field models.
+   *
+   * This is the base interface for magnetic field model mixins.
+   *
+   */
+  template <typename Model>
+  class IMagneticFieldModel : public Model {
+
+    // a type-alias for a magnetic field vector
+    using MagneticFieldVector = Vector<magnetic_flux_density_d>;
+
+  public:
+    /**
+     * Evaluate the magnetic field at a given location.
+     *
+     * @param  point    The location to evaluate the field at.
+     * @returns    The magnetic field vector at that point.
+     */
+    virtual auto getMagneticField(Point const&) const -> MagneticFieldVector = 0;
+
+    /**
+     * A virtual default destructor.
+     */
+    virtual ~IMagneticFieldModel() = default; // LCOV_EXCL_LINE
+
+  }; // END: class MagneticField
+
+} // namespace corsika

corsika/media/IMediumModel.hpp

+/*
+ * (c) Copyright 2018 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/media/NuclearComposition.hpp>
+#include <corsika/framework/geometry/Point.hpp>
+#include <corsika/framework/core/PhysicalUnits.hpp>
+#include <corsika/framework/geometry/BaseTrajectory.hpp>
+
+namespace corsika {
+
+  class IMediumModel {
+  public:
+    virtual ~IMediumModel() = default; // LCOV_EXCL_LINE
+
+    virtual MassDensityType getMassDensity(Point const&) const = 0;
+
+    /**
+     * Integrate the matter density along trajectory.
+     *
+     * @return GrammageType as integrated matter density along the BaseTrajectory
+     *
+     * @todo think about the mixin inheritance of the trajectory vs the BaseTrajectory
+     *       approach; for now, only lines are supported (?).
+     */
+    virtual GrammageType getIntegratedGrammage(BaseTrajectory const&) const = 0;
+
+    /**
+     * Calculates the length along the trajectory.
+     *
+     * The length along the trajectory is determined at which the integrated matter
+     * density is reached. If the specified matter density cannot be reached (is too
+     * large) the result becomes meaningless and could be "infinity" (discuss this).
+     *
+     * @return LengthType the length corresponding to grammage.
+     */
+    virtual LengthType getArclengthFromGrammage(BaseTrajectory const&,
+                                                GrammageType) const = 0;
+
+    virtual NuclearComposition const& getNuclearComposition() const = 0;
+  };
+
+} // namespace corsika

corsika/media/IMediumPropertyModel.hpp

+/*
+ * (c) Copyright 2020 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/media/MediumProperties.hpp>
+
+#include <corsika/framework/geometry/Point.hpp>
+#include <corsika/framework/core/PhysicalUnits.hpp>
+
+namespace corsika {
+
+  /**
+   * An interface for type of media, needed e.g. to determine energy losses.
+   *
+   * This is the base interface for media types.
+   */
+  template <typename TModel>
+  class IMediumPropertyModel : public TModel {
+
+  public:
+    /**
+     * Evaluate the medium type at a given location.
+     *
+     * @param  point    The location to evaluate at.
+     * @returns    The media type
+     */
+    virtual Medium getMedium() const = 0;
+
+    /**
+     * A virtual default destructor.
+     */
+    virtual ~IMediumPropertyModel() = default;
+
+  }; // END: class IMediumTypeModel
+
+} // namespace corsika

corsika/media/IRefractiveIndexModel.hpp

+/*
+ * (c) Copyright 2020 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/geometry/Point.hpp>
+
+namespace corsika {
+
+  /**
+   * An interface for refractive index models.
+   *
+   * This is the base interface for refractive index mixins.
+   *
+   */
+  template <typename TModel>
+  class IRefractiveIndexModel : public TModel {
+
+  public:
+    /**
+     * Evaluate the refractive index at a given location.
+     *
+     * @param  point    The location to evaluate at.
+     * @returns    The refractive index at this point.
+     */
+    virtual double getRefractiveIndex(Point const&) const = 0;
+
+    /**
+     * A virtual default destructor.
+     */
+    virtual ~IRefractiveIndexModel() = default;
+
+  }; // END: class IRefractiveIndexModel
+
+} // namespace corsika

corsika/media/InhomogeneousMedium.hpp

+/*
+ * (c) Copyright 2020 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/geometry/Line.hpp>
+#include <corsika/framework/geometry/Point.hpp>
+#include <corsika/media/NuclearComposition.hpp>
+#include <corsika/framework/geometry/BaseTrajectory.hpp>
+
+namespace corsika {
+
+  /**
+   * A general inhomogeneous medium. The mass density distribution TDensityFunction must
+   * be a \f$C^2\f$-function.
+   */
+
+  template <typename T, typename TDensityFunction>
+  class InhomogeneousMedium : public T {
+
+  public:
+    template <typename... TArgs>
+    InhomogeneousMedium(NuclearComposition const& nuclComp, TArgs&&... rhoTArgs);
+
+    MassDensityType getMassDensity(Point const& point) const override;
+
+    NuclearComposition const& getNuclearComposition() const override;
+
+    GrammageType getIntegratedGrammage(BaseTrajectory const& line) const override;
+
+    LengthType getArclengthFromGrammage(BaseTrajectory const& pLine,
+                                        GrammageType grammage) const override;
+
+  private:
+    NuclearComposition const nuclComp_;
+    TDensityFunction const densityFunction_;
+  };
+
+} // namespace corsika
+
+#include <corsika/detail/media/InhomogeneousMedium.inl>

corsika/media/LayeredSphericalAtmosphereBuilder.hpp

+/*
+ * (c) Copyright 2020 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/geometry/Point.hpp>
+#include <corsika/media/Environment.hpp>
+#include <corsika/media/IMediumModel.hpp>
+#include <corsika/media/NuclearComposition.hpp>
+#include <corsika/media/VolumeTreeNode.hpp>
+
+// for detail namespace, NoExtraModelInner, NoExtraModel and traits
+#include <corsika/detail/media/LayeredSphericalAtmosphereBuilder.hpp>
+
+#include <memory>
+#include <stack>
+#include <tuple>
+#include <type_traits>
+#include <functional>
+
+/**
+ * @file LayeredSphericalAtmosphereBuilder.hpp
+ */
+
+namespace corsika {
+
+  /**
+   * make_layered_spherical_atmosphere_builder.
+   *
+   * Helper class to create LayeredSphericalAtmosphereBuilder, the
+   * extra environment models have to be passed as template-template
+   * argument to make_layered_spherical_atmosphere_builder, the member
+   * function `create` does then take an unspecified number of extra
+   * parameters to internalize those models for all layers later
+   * produced.
+   */
+  template <typename TMediumInterface = IMediumModel,
+            template <typename> typename MExtraEnvirnoment = detail::NoExtraModel>
+  struct make_layered_spherical_atmosphere_builder;
+
+  /**
+   * Helper class to setup concentric spheres of layered atmosphere
+   * with spcified density profiles (exponential, linear, ...).
+   *
+   * This can be used most importantly to replicate CORSIKA7
+   * atmospheres.
+   *
+   * Each layer by definition has a density profile and a (constant)
+   * nuclear composition model.
+   */
+
+  template <typename TMediumInterface = IMediumModel,
+            template <typename> typename TMediumModelExtra = detail::NoExtraModel,
+            typename... TModelArgs>
+  class LayeredSphericalAtmosphereBuilder {
+
+    LayeredSphericalAtmosphereBuilder() = delete;
+    LayeredSphericalAtmosphereBuilder(const LayeredSphericalAtmosphereBuilder&) = delete;
+    LayeredSphericalAtmosphereBuilder(const LayeredSphericalAtmosphereBuilder&&) = delete;
+    LayeredSphericalAtmosphereBuilder& operator=(
+        const LayeredSphericalAtmosphereBuilder&) = delete;
+
+    // friend, to allow construction
+    template <typename, template <typename> typename>
+    friend struct make_layered_spherical_atmosphere_builder;
+
+  protected:
+    LayeredSphericalAtmosphereBuilder(TModelArgs... args, Point const& center,
+                                      LengthType const planetRadius)
+        : center_(center)
+        , planetRadius_(planetRadius)
+        , additionalModelArgs_{args...} {}
+
+  public:
+    typedef typename VolumeTreeNode<TMediumInterface>::VTN_type volume_tree_node;
+    typedef typename VolumeTreeNode<TMediumInterface>::VTNUPtr volume_tree_node_uptr;
+
+    void setNuclearComposition(NuclearComposition const& composition);
+    volume_tree_node* addExponentialLayer(GrammageType const b,
+                                          LengthType const scaleHeight,
+                                          LengthType const upperBoundary);
+    void addLinearLayer(GrammageType const b, LengthType const scaleHeight,
+                        LengthType const upperBoundary);
+
+    void addTabularLayer(std::function<MassDensityType(LengthType)> const& funcRho,
+                         unsigned int const nBins, LengthType const deltaHeight,
+                         LengthType const upperBoundary);
+
+    int getSize() const { return layers_.size(); }
+
+    void assemble(Environment<TMediumInterface>& env);
+    Environment<TMediumInterface> assemble();
+
+    /**
+     * Get the current planet radius.
+     */
+    LengthType getPlanetRadius() const { return planetRadius_; }
+
+  private:
+    void checkRadius(LengthType const r) const;
+
+    std::unique_ptr<NuclearComposition> composition_;
+    Point center_;
+    LengthType previousRadius_{LengthType::zero()};
+    LengthType planetRadius_;
+    std::tuple<TModelArgs...> const additionalModelArgs_;
+
+    std::stack<volume_tree_node_uptr> layers_; // innermost layer first
+
+  }; // end class LayeredSphericalAtmosphereBuilder
+
+} // namespace corsika
+
+#include <corsika/detail/media/LayeredSphericalAtmosphereBuilder.inl>

corsika/media/LinearApproximationIntegrator.hpp

+/*
+ * (c) Copyright 2020 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 <limits>
+
+#include <corsika/framework/geometry/Line.hpp>
+#include <corsika/framework/geometry/BaseTrajectory.hpp>
+
+namespace corsika {
+
+  /**
+   * Helper class to integrate 1D density functions.
+   *
+   * Integrated density (grammage) can be calculated along a Trajectory. Also the inverse
+   * (arclength) for a specified grammage can be calculated. In addition, also a maximum
+   * possible length for a maximum allowed uncertainty can be evaluated.
+   *
+   * @tparam TDerived Must provide the **evaluateAt** and **getFirstDerivative**
+   * methods.
+   */
+
+  template <typename TDerived>
+  class LinearApproximationIntegrator {
+
+    /**
+     * Return the type TDerived of this object.
+     *
+     * @return TDerived const&
+     */
+    TDerived const& getImplementation() const;
+
+  public:
+    /**
+     * Get the integrated Grammage along line.
+     *
+     * Calculated as \f[ X=(\varrho_0 + 0.5* d\varrho/dl \cdot l ) l \f].
+     *
+     * @param line
+     * @return GrammageType
+     */
+    GrammageType getIntegrateGrammage(BaseTrajectory const& line) const;
+
+    /**
+     * Get the arclength from Grammage along line.
+     *
+     * Calculate as \f[ (1 - 0.5 * X * d\varrho/dl / (\varrho_0^2)) * X / \varrho_0 \f].
+     *
+     * @param line
+     * @param grammage
+     * @return LengthType
+     */
+    LengthType getArclengthFromGrammage(BaseTrajectory const& line,
+                                        GrammageType grammage) const;
+
+    /**
+     * Get the maximum length l to keep uncertainty below \f[ relError \f].
+     *
+     * @param line
+     * @param relError
+     * @return LengthType
+     */
+    LengthType getMaximumLength(BaseTrajectory const& line, double relError) const;
+  };
+
+} // namespace corsika
+
+#include <corsika/detail/media/LinearApproximationIntegrator.inl>

corsika/media/MediumProperties.hpp

+/*
+ * (c) Copyright 2020 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
+
+namespace corsika {
+
+  /**
+   * @defgroup MediaProperties Media Properties
+   *
+   * Medium types are useful most importantly for effective models
+   * like energy losses. a particular medium (mixture of components)
+   * may have specif properties not reflected by its mixture of
+   * components.
+   *
+   * The data provided here is automatically parsed from the file
+   * properties8.dat from NIST.
+   *
+   * The data of each known medium can be access via the global functions in namespace
+   * corsika, or via a static class object with the following interface (here at the
+   * example of the class HydrogenGas):
+   *
+   * @code
+   * static constexpr Medium medium() { return Medium::HydrogenGas; }
+   *
+   * static std::string const getName() { return data_.getName(); }
+   * static std::string const getPrettyName() { return data_.getPrettyName(); }
+   * static double getWeight() { return data_.getWeight (); }
+   * static int weight_significant_figure() { return data_.weight_significant_figure (); }
+   * static int weight_error_last_digit() { return data_.weight_error_last_digit(); }
+   * static double Z_over_A() { return data_.Z_over_A(); }
+   * static double getSternheimerDensity() { return data_.getSternheimerDensity(); }
+   * static double getCorrectedDensity() { return data_.getCorrectedDensity(); }
+   * static StateOfMatter getStateOfMatter() { return data_.getStateOfMatter(); }
+   * static MediumType getType() { return data_.getType(); }
+   * static std::string const getSymbol() { return data_.getSymbol(); }
+   *
+   * static double getIeff() { return data_.getIeff(); }
+   * static double getCbar() { return data_.getCbar(); }
+   * static double getX0() { return data_.getX0(); }
+   * static double getX1() { return data_.getX1(); }
+   * static double getAA() { return data_.getAA(); }
+   * static double getSK() { return data_.getSK(); }
+   * static double getDlt0() { return data_.getDlt0(); }
+   *
+   * inline static const MediumData data_ { "hydrogen_gas", "hydrogen gas (H%2#)", 1.008,
+   * 3, 7, 0.99212,
+   * 8.3748e-05, 8.3755e-05, StateOfMatter::DiatomicGas,
+   * MediumType::Element, "H", 19.2, 9.5835, 1.8639, 3.2718, 0.14092, 5.7273, 0.0 };
+   * @endcode
+   *
+   * The numeric data known to CORSIKA 8 (and obtained from NIST) can be browsed in
+   * @ref MediaPropertiesClasses.
+   *
+   * @ingroup MediaProperties
+   * @{
+   */
+
+  //! General type of medium
+  enum class MediumType {
+    Unknown,
+    Element,
+    RadioactiveElement,
+    InorganicCompound,
+    OrganicCompound,
+    Polymer,
+    Mixture,
+    BiologicalDosimetry
+  };
+
+  //! Physical state of medium
+  enum class StateOfMatter { Unknown, Solid, Liquid, Gas, DiatomicGas };
+
+  /**
+   * Enum for all known Media types.
+   */
+  enum class Medium : int16_t;
+
+  /**
+   * The internal integer type of a medium enum.
+   */
+  typedef std::underlying_type<Medium>::type MediumIntType;
+
+  /**
+   *
+   * Simple object to group together the properties of a medium.
+   */
+  struct MediumData {
+    std::string name_;
+    std::string pretty_name_;
+    double weight_;
+    int weight_significant_figure_;
+    int weight_error_last_digit_;
+    double Z_over_A_;
+    double sternheimer_density_;
+    double corrected_density_;
+    StateOfMatter state_;
+    MediumType type_;
+    std::string symbol_;
+    double Ieff_;
+    double Cbar_;
+    double x0_;
+    double x1_;
+    double aa_;
+    double sk_;
+    double dlt0_;
+
+    //! @name MediumDataInterface Interface methods
+    //! Interface functions for MediumData
+    //! @{
+    std::string getName() const { return name_; }              /// returns name
+    std::string getPrettyName() const { return pretty_name_; } /// returns pretty name
+    double getWeight() const { return weight_; }               /// return weight
+    const int& weight_significant_figure() const {
+      return weight_significant_figure_;
+    } /// return significnat figures of weight
+    const int& weight_error_last_digit() const {
+      return weight_error_last_digit_;
+    }                                                    /// return error of weight
+    const double& Z_over_A() const { return Z_over_A_; } /// Z_over_A_
+    double getSternheimerDensity() const {
+      return sternheimer_density_;
+    } /// Sternheimer density
+    double getCorrectedDensity() const {
+      return corrected_density_;
+    }                                                         /// corrected density
+    StateOfMatter getStateOfMatter() const { return state_; } /// state
+    MediumType getType() const { return type_; }              /// type
+    std::string getSymbol() const { return symbol_; }         /// symbol
+    double getIeff() const { return Ieff_; }                  /// Ieff
+    double getCbar() const { return Cbar_; }                  /// Cbar
+    double getX0() const { return x0_; }                      /// X0
+    double getX1() const { return x1_; }                      /// X1
+    double getAA() const { return aa_; }                      /// AA
+    double getSK() const { return sk_; }                      /// Sk
+    double getDlt0() const { return dlt0_; }                  /// Delta0
+    //! @}
+  };
+
+  //! @}
+
+} // namespace corsika
+
+#include <corsika/media/GeneratedMediaProperties.inc>
+
+namespace corsika {
+
+  /**
+   * @file MediaProperties.hpp
+   *
+   * @ingroup MediaProperties
+   * @{
+   *
+   * Returns MediumData object for medium identifed by enum Medium.
+   */
+
+  constexpr MediumData const& mediumData(Medium const m) {
+    return corsika::detail::medium_data[static_cast<MediumIntType>(m)];
+  }
+
+  //! @}
+
+} // namespace corsika

corsika/media/MediumPropertyModel.hpp

+/*
+ * (c) Copyright 2020 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/media/IMediumPropertyModel.hpp>
+
+namespace corsika {
+
+  /**
+   * A model for the energy loss property of a medium.
+   */
+  template <typename T>
+  class MediumPropertyModel : public T {
+
+    Medium medium_; ///< The medium that is used for this medium.
+
+  public:
+    /**
+     * Construct a MediumPropertyModel.
+     *
+     * @param field    The refractive index to return everywhere.
+     */
+    template <typename... Args>
+    MediumPropertyModel(Medium const medium, Args&&... args);
+
+    /**
+     * Evaluate the medium type at a given location.
+     *
+     * @param  point    The location to evaluate at.
+     * @returns    The medium type as enum environment::Medium
+     */
+    Medium getMedium() const override;
+
+    /**
+     * Set the medium type.
+     *
+     * @param  medium    The medium to store.
+     * @returns    The medium type as enum environment::Medium
+     */
+    void setMedium(Medium const medium);
+
+  }; // END: class MediumPropertyModel
+
+} // namespace corsika
+
+#include <corsika/detail/media/MediumPropertyModel.inl>

corsika/media/NameModel.hpp

+/*
+ * (c) Copyright 2020 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 <string>
+
+namespace corsika {
+
+  template <typename T>
+  struct NameModel : public T {
+    virtual std::string const& GetName() const = 0;
+    virtual ~NameModel() = default;
+  };
+
+} // namespace corsika

corsika/media/NuclearComposition.hpp

+/*
+ * (c) Copyright 2020 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/ParticleProperties.hpp>
+#include <corsika/framework/core/PhysicalUnits.hpp>
+
+#include <cassert>
+#include <functional>
+#include <numeric>
+#include <random>
+#include <stdexcept>
+#include <vector>
+
+namespace corsika {
+
+  /**
+   * Describes the composition of matter
+   * Allowes and handles the creation of custom matter compositions.
+   */
+
+  class NuclearComposition {
+  public:
+    /**
+     * Constructor
+     *  The constructore takes a list of elements and a list which describe the relative
+     *  amount. Booth lists need to have the same length and the sum all of fractions
+     *  should be 1. Otherwise an exception is thrown.
+     *  Example for air: Composition (N2,O2,Ar) = (78.084, 20.946, 0.934)
+     *  Pfraction(Ar) = Ar/(2*N2 + 2*O2 + Ar) = 0.00469
+     *  {Code::Nitrogen, Code::Oxygen, Code::Argon}, {0.78479, 0.21052, 0.00469}}
+     *
+     *  @param pComponents List of particle types.
+     *  @param pFractions List of fractions how much each particle contributes. The sum
+     *         needs to add up to 1.
+     */
+    NuclearComposition(std::vector<Code> const& pComponents,
+                       std::vector<double> const& pFractions);
+
+    /**
+     * Returns a vector of the same length as elements in the material with the weighted
+     * return of "func". The typical default application is for cross section weighted
+     * with fraction in the material.
+     *
+     *  @tparam TFunction Type of functions for the weights. The type should be
+     *          Code -> CrossSectionType.
+     *  @param func Functions for reweighting specific elements.
+     *  @retval returns the vector with weighted return types of func.
+     */
+    template <typename TFunction>
+    auto getWeighted(TFunction func) const;
+
+    /**
+     * Sum all all relative composition weighted by func(element)
+     *  This function sums all relative compositions given during this classes
+     * construction. Each entry is weighted by the user defined function func given to
+     * this function.
+     *
+     *  @tparam TFunction Type of functions for the weights. The type should be
+     *          Code -> double.
+     *  @param func Functions for reweighting specific elements.
+     *  @retval returns the weighted sum with the type defined by the return type of func.
+     */
+    template <typename TFunction>
+    auto getWeightedSum(TFunction func) const -> decltype(func(std::declval<Code>()));
+
+    /**
+     * Number of elements in the composition array
+     *  @retval returns the number of elements in the composition array.
+     */
+    size_t getSize() const;
+
+    //! Returns a const reference to the fraction
+    std::vector<double> const& getFractions() const;
+    //! Returns a const reference to the fraction
+    std::vector<Code> const& getComponents() const;
+    double const getAverageMassNumber() const;
+
+    template <class TRNG>
+    Code sampleTarget(std::vector<CrossSectionType> const& sigma,
+                      TRNG&& randomStream) const;
+
+    // Note: when this class ever modifies its internal data, the hash
+    // must be updated, too!
+    size_t getHash() const;
+
+    //! based on hash value
+    bool operator==(NuclearComposition const& v) const;
+
+  private:
+    void updateHash();
+
+    std::vector<double> const numberFractions_; //!< relative fractions of number density
+    std::vector<Code> const components_;        //!< particle codes of consitutents
+
+    double const avgMassNumber_;
+
+    std::size_t hash_;
+  };
+
+} // namespace corsika
+
+#include <corsika/detail/media/NuclearComposition.inl>

corsika/media/ShowerAxis.hpp

+/*
+ * (c) Copyright 2020 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/media/Environment.hpp>
+#include <corsika/framework/geometry/Point.hpp>
+#include <corsika/framework/geometry/Vector.hpp>
+#include <corsika/framework/core/PhysicalUnits.hpp>
+
+#include <cassert>
+#include <cstdlib>
+#include <fstream>
+#include <functional>
+#include <iterator>
+#include <memory>
+#include <stdexcept>
+#include <vector>
+
+#include <iostream>
+
+#include <boost/math/quadrature/gauss_kronrod.hpp>
+
+namespace corsika {
+
+  /**
+   * \class ShowerAxis
+   *
+   * The environment::ShowerAxis is created from a Point and
+   * a Vector and inside an Environment. It internally uses
+   * a table with steps=10000 (default) rows for interpolation.
+   *
+   * The shower axis can convert location in the shower into a
+   * projected grammage along the shower axis.
+   *
+   **/
+
+  ///\todo documentation needs update ...
+  class ShowerAxis {
+  public:
+    template <typename TEnvModel>
+    ShowerAxis(Point const& pStart, Vector<length_d> const& length,
+               Environment<TEnvModel> const& env, bool const doThrow = false,
+               int const steps = 10'000);
+
+    LengthType getSteplength() const;
+
+    GrammageType getMaximumX() const;
+
+    GrammageType getMinimumX() const;
+
+    /**
+     * Returns the grammage along the shower axis of the projection of a point p
+     * onto the shower axis.
+     * Will return either getMinimumX() or getMaximumX() in case the projection is outside
+     * the shower axis.
+     *
+     *  @param p Point to project onto the shower axis.
+     *  @retval Grammage along shower axis for projection of point p.
+     */
+    GrammageType getProjectedX(Point const& p) const;
+
+    /**
+     * Returns the grammage along the shower axis for a given length along the shower
+     * axis. Will return either getMinimumX() or getMaximumX() in case the length is
+     * outside the shower axis.
+     *
+     *  @param l Length along shower axis.
+     *  @retval Grammage along shower axis for length l.
+     */
+    GrammageType getX(LengthType l) const;
+
+    DirectionVector const& getDirection() const;
+
+    Point const& getStart() const;
+
+  private:
+    Point const pointStart_;
+    Vector<length_d> const length_;
+    bool throw_ = false;
+    LengthType const max_length_, steplength_;
+    DirectionVector const axis_normalized_;
+    std::vector<GrammageType> X_;
+
+    // for storing the lengths corresponding to equidistant X values
+    GrammageType const X_binning_ = 1_g / 1_cm / 1_cm;
+    std::vector<LengthType> d_;
+  };
+} // namespace corsika
+
+#include <corsika/detail/media/ShowerAxis.inl>

corsika/media/SlidingPlanarExponential.hpp

+/*
+ * (c) Copyright 2020 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/ParticleProperties.hpp>
+#include <corsika/framework/core/PhysicalUnits.hpp>
+#include <corsika/framework/geometry/Line.hpp>
+#include <corsika/framework/geometry/Point.hpp>
+#include <corsika/framework/random/RNGManager.hpp>
+#include <corsika/media/NuclearComposition.hpp>
+#include <corsika/framework/geometry/BaseTrajectory.hpp>
+
+namespace corsika {
+
+  // clang-format off
+  /**
+   * The SlidingPlanarExponential models mass density as:
+   * \f[
+   *   \varrho(r) = \varrho_0 \exp\left( \frac{|p_0 - r|}{\lambda} \right).
+   * \f]
+   * For grammage/length conversion, the density distribution is approximated as
+   * locally flat at the starting point \f$ r_0 \f$ of the trajectory with the axis pointing
+   * from \f$ p_0 \f$ to \f$ r_0 \f$.
+   */
+  // clang-format on
+
+  template <typename T>
+  class SlidingPlanarExponential : public BaseExponential<SlidingPlanarExponential<T>>,
+                                   public T {
+
+    using Base = BaseExponential<SlidingPlanarExponential<T>>;
+
+  public:
+    SlidingPlanarExponential(Point const& p0, MassDensityType const rho0,
+                             LengthType const lambda, NuclearComposition const& nuclComp,
+                             LengthType const referenceHeight = LengthType::zero());
+
+    MassDensityType getMassDensity(Point const& point) const override;
+
+    NuclearComposition const& getNuclearComposition() const override;
+
+    GrammageType getIntegratedGrammage(BaseTrajectory const& line) const override;
+
+    LengthType getArclengthFromGrammage(BaseTrajectory const& line,
+                                        GrammageType const grammage) const override;
+
+  private:
+    NuclearComposition const nuclComp_;
+  };
+
+} // namespace corsika
+
+#include <corsika/detail/media/SlidingPlanarExponential.inl>

corsika/media/SlidingPlanarTabular.hpp

+/*
+ * (c) Copyright 2020 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/ParticleProperties.hpp>
+#include <corsika/framework/core/PhysicalUnits.hpp>
+#include <corsika/framework/geometry/Line.hpp>
+#include <corsika/framework/geometry/Point.hpp>
+#include <corsika/framework/geometry/BaseTrajectory.hpp>
+#include <corsika/media/NuclearComposition.hpp>
+#include <corsika/media/BaseTabular.hpp>
+
+namespace corsika {
+
+  // clang-format off
+  /**
+   * The SlidingPlanarTabular models mass density as
+   * \f[
+   *   \varrho(r) = \varrho_0 \rho\left( |p_0 - r| \right).
+   * \f]
+   * For grammage/length conversion, the density distribution is approximated as
+   * locally flat at the starting point \f$ r_0 \f$ of the trajectory with the 
+   * axis pointing rom \f$ p_0 \f$ to \f$ r_0 \f$ defining the local height.
+   */
+  // clang-format on
+
+  template <typename TDerived>
+  class SlidingPlanarTabular : public BaseTabular<SlidingPlanarTabular<TDerived>>,
+                               public TDerived {
+
+    using Base = BaseTabular<SlidingPlanarTabular<TDerived>>;
+
+  public:
+    SlidingPlanarTabular(Point const& p0,
+                         std::function<MassDensityType(LengthType)> const& rho,
+                         unsigned int const nBins, LengthType const deltaHeight,
+                         NuclearComposition const& nuclComp,
+                         LengthType referenceHeight = LengthType::zero());
+
+    MassDensityType getMassDensity(Point const& point) const override;
+
+    NuclearComposition const& getNuclearComposition() const override;
+
+    GrammageType getIntegratedGrammage(BaseTrajectory const& line) const override;
+
+    LengthType getArclengthFromGrammage(BaseTrajectory const& line,
+                                        GrammageType grammage) const override;
+
+  private:
+    NuclearComposition const nuclComp_;
+  };
+
+} // namespace corsika
+
+#include <corsika/detail/media/SlidingPlanarTabular.inl>

corsika/media/UniformMagneticField.hpp

+/*
+ * (c) Copyright 2020 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/media/IMagneticFieldModel.hpp>
+#include <corsika/framework/geometry/Vector.hpp>
+#include <corsika/framework/geometry/PhysicalGeometry.hpp>
+
+namespace corsika {
+
+  /**
+   * A uniform (constant) magnetic field.
+   *
+   * This class returns the same magnetic field vector
+   * for all evaluated locations.
+   */
+  template <typename T>
+  class UniformMagneticField : public T {
+
+  public:
+    /**
+     * Construct a UniformMagneticField.
+     *
+     * This is initialized with a fixed magnetic field
+     * and returns this magnetic field at all locations.
+     *
+     * @param field   The fixed magnetic field to return.
+     */
+    template <typename... Args>
+    UniformMagneticField(MagneticFieldVector const& field, Args&&... args)
+        : T(std::forward<Args>(args)...)
+        , B_(field) {}
+
+    /**
+     * Evaluate the magnetic field at a given location.
+     *
+     * @param  point    The location to evaluate the field at (not used internally).
+     * @returns    The magnetic field vector.
+     */
+    MagneticFieldVector getMagneticField([
+        [maybe_unused]] Point const& point) const final override {
+      return B_;
+    }
+
+    /**
+     * Set the magnetic field returned by this instance.
+     *
+     * @param  Bfield    The new vaue of the global magnetic field.
+     */
+    void setMagneticField(MagneticFieldVector const& Bfield) { B_ = Bfield; }
+
+  private:
+    MagneticFieldVector B_; ///< The constant magnetic field we use.
+
+  }; // END: class MagneticField
+
+} // namespace corsika