IAP GITLAB

Skip to content
Snippets Groups Projects
Commit 32d03c35 authored by Remy Prechelt's avatar Remy Prechelt
Browse files

Change interface to use raw double for rindex.

parent 2dc06e42
No related branches found
No related tags found
No related merge requests found
......@@ -23,9 +23,6 @@ namespace corsika::environment {
template <typename Model>
class IRefractiveIndexModel : public Model {
// a type-alias for a dimensionless refractive index
using RefractiveIndex = corsika::units::si::RefractiveIndexType;
public:
/**
* Evaluate the refractive index at a given location.
......@@ -33,7 +30,7 @@ namespace corsika::environment {
* @param point The location to evaluate at.
* @returns The refractive index at this point.
*/
virtual RefractiveIndex GetRefractiveIndex(corsika::geometry::Point const&) const = 0;
virtual double GetRefractiveIndex(corsika::geometry::Point const&) const = 0;
/**
* A virtual default destructor.
......
......@@ -23,10 +23,7 @@ namespace corsika::environment {
template <typename T>
class UniformRefractiveIndex : public T {
// a type-alias for a dimensionless refractive index
using RefractiveIndex = corsika::units::si::RefractiveIndexType;
RefractiveIndex n_; ///< The constant refractive index that we use.
double n_; ///< The constant refractive index that we use.
public:
/**
......@@ -38,7 +35,7 @@ namespace corsika::environment {
* @param field The refractive index to return everywhere.
*/
template <typename... Args>
UniformRefractiveIndex(RefractiveIndex const n, Args&&... args)
UniformRefractiveIndex(double const n, Args&&... args)
: T(std::forward<Args>(args)...)
, n_(n) {}
......@@ -48,8 +45,7 @@ namespace corsika::environment {
* @param point The location to evaluate at.
* @returns The refractive index at this point.
*/
RefractiveIndex GetRefractiveIndex(
corsika::geometry::Point const&) const final override {
double GetRefractiveIndex(corsika::geometry::Point const&) const final override {
return n_;
}
......@@ -59,7 +55,7 @@ namespace corsika::environment {
* @param point The location to evaluate at.
* @returns The refractive index at this location.
*/
void SetRefractiveIndex(RefractiveIndex const& n) { n_ = n; }
void SetRefractiveIndex(double const& n) { n_ = n; }
}; // END: class RefractiveIndex
......
......@@ -19,11 +19,8 @@
#include <corsika/environment/LinearApproximationIntegrator.h>
#include <corsika/environment/NuclearComposition.h>
#include <corsika/environment/SlidingPlanarExponential.h>
<<<<<<< HEAD
#include <corsika/environment/UniformMagneticField.h>
=======
#include <corsika/environment/UniformRefractiveIndex.h>
>>>>>>> a9fca8f... Initial refractive index interface and models with unittests.
#include <corsika/environment/VolumeTreeNode.h>
#include <corsika/geometry/Line.h>
#include <corsika/geometry/RootCoordinateSystem.h>
......@@ -253,12 +250,7 @@ TEST_CASE("UniformMagneticField w/ Homogeneous Medium") {
// the constant density
const auto density{19.2_g / cube(1_cm)};
// the constant density
const auto density{19.2_g / cube(1_cm)};
// create our atmospheric model
<<<<<<< HEAD
<<<<<<< HEAD
AtmModel medium(B0, density, protonComposition);
// and test at several locations
......@@ -294,11 +286,26 @@ TEST_CASE("UniformMagneticField w/ Homogeneous Medium") {
// and check the integrated grammage
REQUIRE((medium.IntegratedGrammage(trajectory, 3_m) / (density * 3_m)) == Approx(1));
REQUIRE((medium.ArclengthFromGrammage(trajectory, density * 5_m) / 5_m) == Approx(1));
=======
AtmModel medium(n, 19.2_g / cube(1_cm), protonComposition);
=======
}
TEST_CASE("UniformRefractiveIndex w/ Homogeneous") {
// setup our interface types
using IModelInterface = IRefractiveIndexModel<IMediumModel>;
using AtmModel = UniformRefractiveIndex<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 double n{1.000327};
// create the atmospheric model
AtmModel medium(n, density, protonComposition);
>>>>>>> a61a703... Add additional density checks for homogeneous medium.
// and require that it is constant
REQUIRE(n == medium.GetRefractiveIndex(Point(gCS, -10_m, 4_m, 35_km)));
......@@ -307,7 +314,7 @@ TEST_CASE("UniformMagneticField w/ Homogeneous Medium") {
REQUIRE(n == medium.GetRefractiveIndex(Point(gCS, 100_km, 400_km, 350_km)));
// a new refractive index
RefractiveIndexType n2 = 2.3472123__;
const double n2{2.3472123};
// update the refractive index of this atmospheric model
medium.SetRefractiveIndex(n2);
......@@ -317,54 +324,10 @@ TEST_CASE("UniformMagneticField w/ Homogeneous Medium") {
REQUIRE(n2 == medium.GetRefractiveIndex(Point(gCS, +210_m, 0_m, 7_km)));
REQUIRE(n2 == medium.GetRefractiveIndex(Point(gCS, 0_m, 0_m, 0_km)));
REQUIRE(n2 == medium.GetRefractiveIndex(Point(gCS, 100_km, 400_km, 350_km)));
<<<<<<< HEAD
>>>>>>> 12a0be6... Add test for SetRefractiveIndex
}
TEST_CASE("UniformRefractiveIndex w/ FlatExponential") {
// setup our interface types
using IModelInterface = IRefractiveIndexModel<IMediumModel>;
using AtmModel = UniformRefractiveIndex<FlatExponential<IModelInterface>>;
// the composition we use for the homogenous medium
NuclearComposition const protonComposition(std::vector<Code>{Code::Proton},
std::vector<float>{1.f});
// define our axis vector
Vector const axis(gCS, QuantityVector<dimensionless_d>(0, 0, 1));
// the parameters of our exponential model
LengthType const lambda = 3_m;
auto const rho0 = 1_g / units::si::detail::static_pow<3>(1_cm);
// the refractive index we want returned everywhere
RefractiveIndexType n = 1.000327__;
=======
>>>>>>> a61a703... Add additional density checks for homogeneous medium.
// check the density and nuclear composition
REQUIRE(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}));
<<<<<<< HEAD
<<<<<<< HEAD
// the refractive index we want returned everywhere
RefractiveIndexType n = 1.1234__;
// create our atmospheric model
AtmModel const medium(n, 19.2_g / cube(1_cm), protonComposition);
// and require that it is constant
REQUIRE(n == medium.GetRefractiveIndex(Point(gCS, -10_m, 4_m, 35_km)));
REQUIRE(n == medium.GetRefractiveIndex(Point(gCS, +210_m, 0_m, 7_km)));
REQUIRE(n == medium.GetRefractiveIndex(Point(gCS, 0_m, 0_m, 0_km)));
REQUIRE(n == medium.GetRefractiveIndex(Point(gCS, 100_km, 400_km, 350_km)));
// check the density and nuclear composition
REQUIRE(density == medium.GetMassDensity(Point(gCS, 0_m, 0_m, 0_m)));
medium.GetNuclearComposition();
......@@ -382,26 +345,4 @@ TEST_CASE("UniformRefractiveIndex w/ FlatExponential") {
// and check the integrated grammage
REQUIRE((medium.IntegratedGrammage(trajectory, 3_m) / (density * 3_m)) == Approx(1));
REQUIRE((medium.ArclengthFromGrammage(trajectory, density * 5_m) / 5_m) == Approx(1));
=======
// a new refractive index
RefractiveIndexType n2 = 1.123456__;
=======
// the end time of our line
auto const tEnd = 1_s;
>>>>>>> a61a703... Add additional density checks for homogeneous medium.
// and the associated trajectory
Trajectory<Line> const trajectory(line, tEnd);
<<<<<<< HEAD
// check that the returned refractive index is correct
REQUIRE(n2 == medium.GetRefractiveIndex(Point(gCS, -10_m, 4_m, 35_km)));
REQUIRE(n2 == medium.GetRefractiveIndex(Point(gCS, +210_m, 0_m, 7_km)));
REQUIRE(n2 == medium.GetRefractiveIndex(Point(gCS, 0_m, 0_m, 0_km)));
REQUIRE(n2 == medium.GetRefractiveIndex(Point(gCS, 100_km, 400_km, 350_km)));
>>>>>>> 12a0be6... Add test for SetRefractiveIndex
=======
// and check the integrated grammage
REQUIRE((medium.IntegratedGrammage(trajectory, 3_m) / (density * 3_m)) == Approx(1));
>>>>>>> a61a703... Add additional density checks for homogeneous medium.
}
......@@ -69,12 +69,6 @@ namespace corsika::units::si {
phys::units::quantity<phys::units::dimensions<0, 0, -1>, double>;
using InverseGrammageType =
phys::units::quantity<phys::units::dimensions<2, -1, 0>, double>;
using RefractiveIndexType = phys::units::quantity<dimensionless_d, double>;
// define a new literal to construct dimensionless types
// This is a 'double underscore' literal, So 3__ is a
// dimensionless quantity = 3
QUANTITY_DEFINE_LITERALS(_, dimensionless_d)
namespace detail {
template <int N, typename T>
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment