From acc57b53dd6902de67f11f619a0d467c84c7e56e Mon Sep 17 00:00:00 2001
From: Remy Prechelt <prechelt@hawaii.edu>
Date: Mon, 16 Nov 2020 23:18:59 -1000
Subject: [PATCH] Port over medium property models.
---
corsika/detail/media/MediumPropertyModel.inl | 31 +++++++
corsika/media/IMediumModel.hpp | 19 ++---
corsika/media/IMediumPropertyModel.hpp | 43 ++++++++++
corsika/media/MediumPropertyModel.hpp | 53 ++++++++++++
tests/media/testMedium.cpp | 87 ++++++++++++++++++++
5 files changed, 223 insertions(+), 10 deletions(-)
create mode 100644 corsika/detail/media/MediumPropertyModel.inl
create mode 100644 corsika/media/IMediumPropertyModel.hpp
create mode 100644 corsika/media/MediumPropertyModel.hpp
create mode 100644 tests/media/testMedium.cpp
diff --git a/corsika/detail/media/MediumPropertyModel.inl b/corsika/detail/media/MediumPropertyModel.inl
new file mode 100644
index 000000000..3e584252b
--- /dev/null
+++ b/corsika/detail/media/MediumPropertyModel.inl
@@ -0,0 +1,31 @@
+/*
+ * (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/media/IMediumPropertyModel.hpp>
+
+namespace corsika {
+
+ template <typename T>
+ template <typename... Args>
+ MediumPropertyModel<T>::MediumPropertyModel(Medium const medium, Args&&... args)
+ : T(std::forward<Args>(args)...)
+ , medium_(medium) {}
+
+ template <typename T>
+ Medium MediumPropertyModel<T>::getMedium(Point const&) const final override {
+ return medium_;
+ }
+
+ template <typename T>
+ void MediumPropertyModel<T>::setMedium(Medium const medium) {
+ medium_ = medium;
+ }
+
+} // namespace corsika
diff --git a/corsika/media/IMediumModel.hpp b/corsika/media/IMediumModel.hpp
index 538657cb6..c63819251 100644
--- a/corsika/media/IMediumModel.hpp
+++ b/corsika/media/IMediumModel.hpp
@@ -1,5 +1,5 @@
/*
- * (c) Copyright 2020 CORSIKA Project, corsika-project@lists.kit.edu
+ * (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
@@ -8,30 +8,29 @@
#pragma once
-#include <corsika/framework/core/PhysicalUnits.hpp>
+#include <corsika/media/NuclearComposition.hpp>
#include <corsika/framework/geometry/Line.hpp>
#include <corsika/framework/geometry/Point.hpp>
#include <corsika/framework/geometry/Trajectory.hpp>
-#include <corsika/media/NuclearComposition.hpp>
+#include <corsika/framework/core/PhysicalUnits.hpp>
namespace corsika {
class IMediumModel {
-
public:
virtual ~IMediumModel() = default; // LCOV_EXCL_LINE
- virtual MassDensityType GetMassDensity(Point const&) const = 0;
+ virtual units::si::MassDensityType getMassDensity(Point const&) const = 0;
// todo: think about the mixin inheritance of the trajectory vs the BaseTrajectory
// approach; for now, only lines are supported
- virtual GrammageType IntegratedGrammage(Trajectory<Line> const&,
- LengthType) const = 0;
+ virtual units::si::GrammageType integratedGrammage(
+ Trajectory<Line> const&, units::si::LengthType) const = 0;
- virtual LengthType ArclengthFromGrammage(Trajectory<Line> const&,
- GrammageType) const = 0;
+ virtual units::si::LengthType arclengthFromGrammage(
+ Trajectory<Line> const&, units::si::GrammageType) const = 0;
- virtual NuclearComposition const& GetNuclearComposition() const = 0;
+ virtual NuclearComposition const& getNuclearComposition() const = 0;
};
} // namespace corsika
diff --git a/corsika/media/IMediumPropertyModel.hpp b/corsika/media/IMediumPropertyModel.hpp
new file mode 100644
index 000000000..4455a6242
--- /dev/null
+++ b/corsika/media/IMediumPropertyModel.hpp
@@ -0,0 +1,43 @@
+/*
+ * (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/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(Point const&) const = 0;
+
+ /**
+ * A virtual default destructor.
+ */
+ virtual ~IMediumPropertyModel() = default;
+
+ }; // END: class IMediumTypeModel
+
+} // namespace corsika
diff --git a/corsika/media/MediumPropertyModel.hpp b/corsika/media/MediumPropertyModel.hpp
new file mode 100644
index 000000000..bf47b9514
--- /dev/null
+++ b/corsika/media/MediumPropertyModel.hpp
@@ -0,0 +1,53 @@
+/*
+ * (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/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(Point const&) 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) override;
+
+ }; // END: class MediumPropertyModel
+
+} // namespace corsika
+
+#include <corsika/detail/media/MediumPropertyModel.inl>
diff --git a/tests/media/testMedium.cpp b/tests/media/testMedium.cpp
new file mode 100644
index 000000000..d963c5d1e
--- /dev/null
+++ b/tests/media/testMedium.cpp
@@ -0,0 +1,87 @@
+/*
+ * (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/framework/core/PhysicalUnits.hpp>
+#include <corsika/framework/geometry/Line.hpp>
+#include <corsika/framework/geometry/RootCoordinateSystem.hpp>
+#include <corsika/framework/geometry/Vector.hpp>
+#include <corsika/media/FlatExponential.hpp>
+#include <corsika/media/HomogeneousMedium.hpp>
+#include <corsika/media/IMediumModel.hpp>
+#include <corsika/media/InhomogeneousMedium.hpp>
+#include <corsika/media/NuclearComposition.hpp>
+#include <corsika/media/MediumPropertyModel.hpp>
+
+#include <catch2/catch.hpp>
+
+using namespace corsika;
+using namespace corsika::units::si;
+
+TEST_CASE("MediumPropertyModel w/ Homogeneous") {
+
+ CoordinateSystem const& gCS =
+ RootCoordinateSystem::GetInstance().GetRootCoordinateSystem();
+
+ Point const gOrigin(gCS, {0_m, 0_m, 0_m});
+
+ // setup our interface types
+ using IModelInterface = IMediumPropertyModel<IMediumModel>;
+ using AtmModel = MediumPropertyModel<HomogeneousMedium<IModelInterface>>;
+
+ // the constant density
+ const auto density{19.2_g / cube(1_cm)};
+
+ // the composition we use for the homogenous medium
+ NuclearComposition const protonComposition(std::vector<Code>{Code::Proton},
+ std::vector<float>{1.f});
+
+ // the refrative index that we use
+ const Medium type = Medium::AirDry1Atm;
+
+ // create the atmospheric model
+ AtmModel medium(type, density, protonComposition);
+
+ // and require that it is constant
+ CHECK(type == medium.medium(Point(gCS, -10_m, 4_m, 35_km)));
+ CHECK(type == medium.medium(Point(gCS, +210_m, 0_m, 7_km)));
+ CHECK(type == medium.medium(Point(gCS, 0_m, 0_m, 0_km)));
+ CHECK(type == medium.medium(Point(gCS, 100_km, 400_km, 350_km)));
+
+ // a new refractive index
+ const Medium type2 = Medium::StandardRock;
+
+ // update the refractive index of this atmospheric model
+ medium.set_medium(type2);
+
+ // check that the returned refractive index is correct
+ CHECK(type2 == medium.medium(Point(gCS, -10_m, 4_m, 35_km)));
+ CHECK(type2 == medium.medium(Point(gCS, +210_m, 0_m, 7_km)));
+ CHECK(type2 == medium.medium(Point(gCS, 0_m, 0_m, 0_km)));
+ CHECK(type2 == medium.medium(Point(gCS, 100_km, 400_km, 350_km)));
+
+ // define our axis vector
+ Vector const axis(gCS, QuantityVector<dimensionless_d>(0, 0, 1));
+
+ // check the density and nuclear composition
+ CHECK(density == medium.GetMassDensity(Point(gCS, 0_m, 0_m, 0_m)));
+ medium.GetNuclearComposition();
+
+ // create a line of length 1 m
+ Line const line(gOrigin, Vector<SpeedType::dimension_type>(
+ gCS, {1_m / second, 0_m / second, 0_m / second}));
+
+ // the end time of our line
+ auto const tEnd = 1_s;
+
+ // and the associated trajectory
+ Trajectory<Line> const trajectory(line, tEnd);
+
+ // and check the integrated grammage
+ CHECK((medium.IntegratedGrammage(trajectory, 3_m) / (density * 3_m)) == Approx(1));
+ CHECK((medium.ArclengthFromGrammage(trajectory, density * 5_m) / 5_m) == Approx(1));
+}
--
GitLab