From 2443b75763286e7685e10e38abf4acc5942b370b Mon Sep 17 00:00:00 2001
From: ralfulrich <ralf.ulrich@kit.edu>
Date: Thu, 17 Dec 2020 11:12:43 +0100
Subject: [PATCH] coding guidelines
---
.../detail/framework/stack/SecondaryView.inl | 243 ++++++-------
corsika/detail/media/BaseExponential.inl | 8 +-
corsika/detail/media/Environment.inl | 8 +-
corsika/detail/media/FlatExponential.inl | 2 +-
corsika/detail/media/HomogeneousMedium.inl | 10 +-
.../LayeredSphericalAtmosphereBuilder.inl | 6 +-
.../media/LinearApproximationIntegrator.inl | 16 +-
corsika/detail/media/NuclearComposition.inl | 53 +--
corsika/detail/media/ShowerAxis.inl | 10 +-
.../detail/media/SlidingPlanarExponential.inl | 23 +-
corsika/detail/media/Universe.inl | 6 +-
corsika/detail/media/VolumeTreeNode.inl | 46 +--
.../modules/energy_loss/BetheBlochPDG.inl | 15 +-
corsika/media/BaseExponential.hpp | 6 +-
corsika/media/Environment.hpp | 4 +-
corsika/media/FlatExponential.hpp | 17 +-
corsika/media/HomogeneousMedium.hpp | 20 +-
corsika/media/IMagneticFieldModel.hpp | 4 +-
corsika/media/IMediumModel.hpp | 10 +-
corsika/media/InhomogeneousMedium.hpp | 18 +-
.../LayeredSphericalAtmosphereBuilder.hpp | 2 +-
.../media/LinearApproximationIntegrator.hpp | 11 +-
corsika/media/NuclearComposition.hpp | 45 +--
corsika/media/ShowerAxis.hpp | 32 +-
corsika/media/SlidingPlanarExponential.hpp | 23 +-
corsika/media/Universe.hpp | 4 +-
corsika/media/VolumeTreeNode.hpp | 51 ++-
tests/media/testEnvironment.cpp | 340 +++---------------
tests/media/testMedium.cpp | 7 +-
tests/media/testRefractiveIndex.cpp | 7 +-
tests/media/testShowerAxis.cpp | 31 +-
31 files changed, 382 insertions(+), 696 deletions(-)
diff --git a/corsika/detail/framework/stack/SecondaryView.inl b/corsika/detail/framework/stack/SecondaryView.inl
index a21c8b2d0..cfbd1d35a 100644
--- a/corsika/detail/framework/stack/SecondaryView.inl
+++ b/corsika/detail/framework/stack/SecondaryView.inl
@@ -14,152 +14,150 @@
#include <vector>
namespace corsika {
-/*
-template <typename TStackDataType, template <typename> typename TParticleInterface,
- template <typename T1, template <class> class T2> class MSecondaryProducer>
-SecondaryView< TStackDataType,TParticleInterface, MSecondaryProducer>::
-
-*/
-template <typename TStackDataType, template <typename> typename TParticleInterface,
- template <typename T1, template <class> class T2> class MSecondaryProducer>
-template <typename... Args>
-typename SecondaryView< TStackDataType,TParticleInterface, MSecondaryProducer>::stack_view_iterator
-SecondaryView< TStackDataType,TParticleInterface, MSecondaryProducer>::addSecondary(const Args... v) {
- CORSIKA_LOG_TRACE("SecondaryView::addSecondary(Args&&)");
- stack_view_iterator proj = getProjectile(); // make this const
- return addSecondary(proj, v...);
- }
-
+ /*
+ template <typename TStackDataType, template <typename> typename TParticleInterface,
+ template <typename T1, template <class> class T2> class MSecondaryProducer>
+ SecondaryView< TStackDataType,TParticleInterface, MSecondaryProducer>::
+
+ */
+ template <typename TStackDataType, template <typename> typename TParticleInterface,
+ template <typename T1, template <class> class T2> class MSecondaryProducer>
+ template <typename... Args>
+ typename SecondaryView<TStackDataType, TParticleInterface,
+ MSecondaryProducer>::stack_view_iterator
+ SecondaryView<TStackDataType, TParticleInterface, MSecondaryProducer>::addSecondary(
+ const Args... v) {
+ CORSIKA_LOG_TRACE("SecondaryView::addSecondary(Args&&)");
+ stack_view_iterator proj = getProjectile(); // make this const
+ return addSecondary(proj, v...);
+ }
-template <typename TStackDataType, template <typename> typename TParticleInterface,
- template <typename T1, template <class> class T2> class MSecondaryProducer>
-typename SecondaryView< TStackDataType,TParticleInterface, MSecondaryProducer>::stack_view_iterator
-SecondaryView< TStackDataType,TParticleInterface, MSecondaryProducer>::begin() {
- unsigned int i = 0;
- for (; i < getSize(); ++i) {
+ template <typename TStackDataType, template <typename> typename TParticleInterface,
+ template <typename T1, template <class> class T2> class MSecondaryProducer>
+ typename SecondaryView<TStackDataType, TParticleInterface,
+ MSecondaryProducer>::stack_view_iterator
+ SecondaryView<TStackDataType, TParticleInterface, MSecondaryProducer>::begin() {
+ unsigned int i = 0;
+ for (; i < getSize(); ++i) {
- if (!isErased(i)) break;
- }
- return stack_view_iterator(*this, i + 1);
+ if (!isErased(i)) break;
}
+ return stack_view_iterator(*this, i + 1);
+ }
-template <typename TStackDataType, template <typename> typename TParticleInterface,
- template <typename T1, template <class> class T2> class MSecondaryProducer>
-typename SecondaryView< TStackDataType,TParticleInterface, MSecondaryProducer>::const_stack_view_iterator
-SecondaryView< TStackDataType,TParticleInterface, MSecondaryProducer>::begin() const {
- unsigned int i = 0;
- for (; i < getSize(); ++i) {
- if (!isErased(i)) break;
- }
-
- return const_stack_view_iterator(*this, i + 1);
+ template <typename TStackDataType, template <typename> typename TParticleInterface,
+ template <typename T1, template <class> class T2> class MSecondaryProducer>
+ typename SecondaryView<TStackDataType, TParticleInterface,
+ MSecondaryProducer>::const_stack_view_iterator
+ SecondaryView<TStackDataType, TParticleInterface, MSecondaryProducer>::begin() const {
+ unsigned int i = 0;
+ for (; i < getSize(); ++i) {
+ if (!isErased(i)) break;
}
+ return const_stack_view_iterator(*this, i + 1);
+ }
-template <typename TStackDataType, template <typename> typename TParticleInterface,
- template <typename T1, template <class> class T2> class MSecondaryProducer>
-typename SecondaryView< TStackDataType,TParticleInterface, MSecondaryProducer>::const_stack_view_iterator
-SecondaryView< TStackDataType,TParticleInterface, MSecondaryProducer>::cbegin() const {
- unsigned int i = 0;
- for (; i < getSize(); ++i) {
- if (!isErased(i)) break;
- }
-
- return const_stack_view_iterator(*this, i + 1);
+ template <typename TStackDataType, template <typename> typename TParticleInterface,
+ template <typename T1, template <class> class T2> class MSecondaryProducer>
+ typename SecondaryView<TStackDataType, TParticleInterface,
+ MSecondaryProducer>::const_stack_view_iterator
+ SecondaryView<TStackDataType, TParticleInterface, MSecondaryProducer>::cbegin() const {
+ unsigned int i = 0;
+ for (; i < getSize(); ++i) {
+ if (!isErased(i)) break;
}
-template <typename TStackDataType, template <typename> typename TParticleInterface,
- template <typename T1, template <class> class T2> class MSecondaryProducer>
-typename SecondaryView< TStackDataType,TParticleInterface, MSecondaryProducer>::stack_view_iterator
-SecondaryView< TStackDataType,TParticleInterface, MSecondaryProducer>::last() {
- unsigned int i = 0;
- for (; i < getSize(); ++i) {
- if (!isErased(getSize() - 1 - i)) break;
- }
- return stack_view_iterator(*this, getSize() - 1 - i + 1);
- }
+ return const_stack_view_iterator(*this, i + 1);
+ }
-template <typename TStackDataType, template <typename> typename TParticleInterface,
- template <typename T1, template <class> class T2> class MSecondaryProducer>
-typename SecondaryView< TStackDataType,TParticleInterface, MSecondaryProducer>::const_stack_view_iterator
-SecondaryView< TStackDataType,TParticleInterface, MSecondaryProducer>::last() const {
- unsigned int i = 0;
- for (; i < getSize(); ++i) {
- if (!isErased(getSize() - 1 - i)) break;
- }
- return stack_view_iterator(*this, getSize() - 1 - i + 1);
+ template <typename TStackDataType, template <typename> typename TParticleInterface,
+ template <typename T1, template <class> class T2> class MSecondaryProducer>
+ typename SecondaryView<TStackDataType, TParticleInterface,
+ MSecondaryProducer>::stack_view_iterator
+ SecondaryView<TStackDataType, TParticleInterface, MSecondaryProducer>::last() {
+ unsigned int i = 0;
+ for (; i < getSize(); ++i) {
+ if (!isErased(getSize() - 1 - i)) break;
}
+ return stack_view_iterator(*this, getSize() - 1 - i + 1);
+ }
-template <typename TStackDataType, template <typename> typename TParticleInterface,
- template <typename T1, template <class> class T2> class MSecondaryProducer>
-typename SecondaryView< TStackDataType,TParticleInterface, MSecondaryProducer>::const_stack_view_iterator
-SecondaryView< TStackDataType,TParticleInterface, MSecondaryProducer>::last() const {
- unsigned int i = 0;
- for (; i < getSize(); ++i) {
- if (!isErased(getSize() - 1 - i)) break;
- }
- return stack_view_iterator(*this, getSize() - 1 - i + 1);
+ template <typename TStackDataType, template <typename> typename TParticleInterface,
+ template <typename T1, template <class> class T2> class MSecondaryProducer>
+ typename SecondaryView<TStackDataType, TParticleInterface,
+ MSecondaryProducer>::const_stack_view_iterator
+ SecondaryView<TStackDataType, TParticleInterface, MSecondaryProducer>::last() const {
+ unsigned int i = 0;
+ for (; i < getSize(); ++i) {
+ if (!isErased(getSize() - 1 - i)) break;
}
+ return stack_view_iterator(*this, getSize() - 1 - i + 1);
+ }
-template <typename TStackDataType, template <typename> typename TParticleInterface,
- template <typename T1, template <class> class T2> class MSecondaryProducer>
-typename SecondaryView< TStackDataType,TParticleInterface, MSecondaryProducer>::const_stack_view_iterator
-SecondaryView< TStackDataType,TParticleInterface, MSecondaryProducer>::clast() const {
+ template <typename TStackDataType, template <typename> typename TParticleInterface,
+ template <typename T1, template <class> class T2> class MSecondaryProducer>
+ typename SecondaryView<TStackDataType, TParticleInterface,
+ MSecondaryProducer>::const_stack_view_iterator
+ SecondaryView<TStackDataType, TParticleInterface, MSecondaryProducer>::clast() const {
unsigned int i = 0;
for (; i < getSize(); ++i) {
if (!isErased(getSize() - 1 - i)) break;
}
return const_stack_view_iterator(*this, getSize() - 1 - i + 1);
- }
-
+ }
-template <typename TStackDataType, template <typename> typename TParticleInterface,
- template <typename T1, template <class> class T2> class MSecondaryProducer>
-void SecondaryView< TStackDataType,TParticleInterface, MSecondaryProducer>::swap(stack_view_iterator a, stack_view_iterator b) {
+ template <typename TStackDataType, template <typename> typename TParticleInterface,
+ template <typename T1, template <class> class T2> class MSecondaryProducer>
+ void SecondaryView<TStackDataType, TParticleInterface, MSecondaryProducer>::swap(
+ stack_view_iterator a, stack_view_iterator b) {
CORSIKA_LOG_TRACE("View::swap");
inner_stack_.swap(getIndexFromIterator(a.getIndex()),
getIndexFromIterator(b.getIndex()));
- }
+ }
-template <typename TStackDataType, template <typename> typename TParticleInterface,
- template <typename T1, template <class> class T2> class MSecondaryProducer>
-void SecondaryView< TStackDataType,TParticleInterface, MSecondaryProducer>::copy(stack_view_iterator a, stack_view_iterator b) {
+ template <typename TStackDataType, template <typename> typename TParticleInterface,
+ template <typename T1, template <class> class T2> class MSecondaryProducer>
+ void SecondaryView<TStackDataType, TParticleInterface, MSecondaryProducer>::copy(
+ stack_view_iterator a, stack_view_iterator b) {
+ CORSIKA_LOG_TRACE("View::copy");
+ inner_stack_.copy(getIndexFromIterator(a.getIndex()),
+ getIndexFromIterator(b.getIndex()));
+ }
+
+ template <typename TStackDataType, template <typename> typename TParticleInterface,
+ template <typename T1, template <class> class T2> class MSecondaryProducer>
+ void SecondaryView<TStackDataType, TParticleInterface, MSecondaryProducer>::copy(
+ const_stack_view_iterator a, stack_view_iterator b) {
CORSIKA_LOG_TRACE("View::copy");
inner_stack_.copy(getIndexFromIterator(a.getIndex()),
getIndexFromIterator(b.getIndex()));
- }
+ }
-template <typename TStackDataType, template <typename> typename TParticleInterface,
- template <typename T1, template <class> class T2> class MSecondaryProducer>
-void SecondaryView< TStackDataType,TParticleInterface, MSecondaryProducer>::copy(const_stack_view_iterator a, stack_view_iterator b) {
+ template <typename TStackDataType, template <typename> typename TParticleInterface,
+ template <typename T1, template <class> class T2> class MSecondaryProducer>
+ void SecondaryView<TStackDataType, TParticleInterface, MSecondaryProducer>::erase(
+ stack_view_iterator p) {
- CORSIKA_LOG_TRACE("View::copy");
- inner_stack_.copy(getIndexFromIterator(a.getIndex()),
- getIndexFromIterator(b.getIndex()));
+ CORSIKA_LOG_TRACE("SecondaryView::Delete");
+ if (isEmpty()) { /*error*/
+ throw std::runtime_error("Stack, cannot delete entry since size is zero");
}
-
-template <typename TStackDataType, template <typename> typename TParticleInterface,
- template <typename T1, template <class> class T2> class MSecondaryProducer>
-void SecondaryView< TStackDataType,TParticleInterface, MSecondaryProducer>::erase(stack_view_iterator p) {
-
- CORSIKA_LOG_TRACE("SecondaryView::Delete");
- if (isEmpty()) { /*error*/
- throw std::runtime_error("Stack, cannot delete entry since size is zero");
- }
- if (isErased(p.getIndex() - 1)) { /*error*/
- throw std::runtime_error("Stack, cannot delete entry since already deleted");
- }
- inner_stack_.erase(getIndexFromIterator(p.getIndex()));
- inner_stack_reference_type::nDeleted_++; // also count in SecondaryView
+ if (isErased(p.getIndex() - 1)) { /*error*/
+ throw std::runtime_error("Stack, cannot delete entry since already deleted");
}
+ inner_stack_.erase(getIndexFromIterator(p.getIndex()));
+ inner_stack_reference_type::nDeleted_++; // also count in SecondaryView
+ }
-template <typename TStackDataType, template <typename> typename TParticleInterface,
- template <typename T1, template <class> class T2> class MSecondaryProducer>
-bool SecondaryView< TStackDataType,TParticleInterface, MSecondaryProducer>::purgeLastIfDeleted() {
+ template <typename TStackDataType, template <typename> typename TParticleInterface,
+ template <typename T1, template <class> class T2> class MSecondaryProducer>
+ bool SecondaryView<TStackDataType, TParticleInterface,
+ MSecondaryProducer>::purgeLastIfDeleted() {
CORSIKA_LOG_TRACE("SecondaryView::purgeLastIfDeleted");
if (!isErased(getSize() - 1))
return false; // the last particle is not marked for deletion. Do nothing.
@@ -169,9 +167,9 @@ bool SecondaryView< TStackDataType,TParticleInterface, MSecondaryProducer>::purg
return true;
}
-template <typename TStackDataType, template <typename> typename TParticleInterface,
- template <typename T1, template <class> class T2> class MSecondaryProducer>
-void SecondaryView< TStackDataType,TParticleInterface, MSecondaryProducer>::purge() {
+ template <typename TStackDataType, template <typename> typename TParticleInterface,
+ template <typename T1, template <class> class T2> class MSecondaryProducer>
+ void SecondaryView<TStackDataType, TParticleInterface, MSecondaryProducer>::purge() {
unsigned int iStack = 0;
unsigned int size = getSize();
while (iStack < size) {
@@ -186,9 +184,10 @@ void SecondaryView< TStackDataType,TParticleInterface, MSecondaryProducer>::purg
inner_stack_reference_type::nDeleted_ = 0;
}
-template <typename TStackDataType, template <typename> typename TParticleInterface,
- template <typename T1, template <class> class T2> class MSecondaryProducer>
-std::string SecondaryView< TStackDataType,TParticleInterface, MSecondaryProducer>::as_string() const {
+ template <typename TStackDataType, template <typename> typename TParticleInterface,
+ template <typename T1, template <class> class T2> class MSecondaryProducer>
+ std::string SecondaryView<TStackDataType, TParticleInterface,
+ MSecondaryProducer>::as_string() const {
std::string str(fmt::format("size {}\n", getSize()));
// we make our own begin/end since we want ALL entries
std::string new_line = " ";
@@ -203,11 +202,13 @@ std::string SecondaryView< TStackDataType,TParticleInterface, MSecondaryProducer
return str;
}
-template <typename TStackDataType, template <typename> typename TParticleInterface,
- template <typename T1, template <class> class T2> class MSecondaryProducer>
-template <typename... Args>
-typename SecondaryView< TStackDataType,TParticleInterface, MSecondaryProducer>::stack_view_iterator
-SecondaryView< TStackDataType,TParticleInterface, MSecondaryProducer>::addSecondary(stack_view_iterator& proj, const Args... v) {
+ template <typename TStackDataType, template <typename> typename TParticleInterface,
+ template <typename T1, template <class> class T2> class MSecondaryProducer>
+ template <typename... Args>
+ typename SecondaryView<TStackDataType, TParticleInterface,
+ MSecondaryProducer>::stack_view_iterator
+ SecondaryView<TStackDataType, TParticleInterface, MSecondaryProducer>::addSecondary(
+ stack_view_iterator& proj, const Args... v) {
CORSIKA_LOG_TRACE("SecondaryView::addSecondary(stack_view_iterator&, Args&&)");
// make space on stack
inner_stack_reference_type::getStackData().incrementSize();
diff --git a/corsika/detail/media/BaseExponential.inl b/corsika/detail/media/BaseExponential.inl
index dc9719c3b..c55847cc5 100644
--- a/corsika/detail/media/BaseExponential.inl
+++ b/corsika/detail/media/BaseExponential.inl
@@ -29,8 +29,8 @@ namespace corsika {
Vector<dimensionless_d> const& axis) const {
if (vL == LengthType::zero()) { return GrammageType::zero(); }
- auto const uDotA = line.NormalizedDirection().dot(axis).magnitude();
- auto const rhoStart = getImplementation().getMassDensity(line.getR0());
+ auto const uDotA = line.getNormalizedDirection().dot(axis).magnitude();
+ auto const rhoStart = getImplementation().getMassDensity(line.getStartPoint());
if (uDotA == 0) {
return vL * rhoStart;
@@ -43,8 +43,8 @@ namespace corsika {
LengthType BaseExponential<TDerived>::getArclengthFromGrammage(
Trajectory<Line> const& line, GrammageType grammage,
Vector<dimensionless_d> const& axis) const {
- auto const uDotA = line.NormalizedDirection().dot(axis).magnitude();
- auto const rhoStart = getImplementation().getMassDensity(line.getR0());
+ auto const uDotA = line.getNormalizedDirection().dot(axis).magnitude();
+ auto const rhoStart = getImplementation().getMassDensity(line.getStartPoint());
if (uDotA == 0) {
return grammage / rhoStart;
diff --git a/corsika/detail/media/Environment.inl b/corsika/detail/media/Environment.inl
index b1fc4bc62..855cc0b3b 100644
--- a/corsika/detail/media/Environment.inl
+++ b/corsika/detail/media/Environment.inl
@@ -16,7 +16,7 @@ namespace corsika {
template <typename IEnvironmentModel>
Environment<IEnvironmentModel>::Environment()
- : coordinateSystem_{RootCoordinateSystem::getInstance().GetRootCoordinateSystem()}
+ : coordinateSystem_{get_root_CoordinateSystem()}
, universe_(std::make_unique<BaseNodeType>(
std::make_unique<Universe>(coordinateSystem_))) {}
@@ -31,7 +31,7 @@ namespace corsika {
}
template <typename IEnvironmentModel>
- CoordinateSystem const& Environment<IEnvironmentModel>::getCoordinateSystem() const {
+ CoordinateSystemPtr const& Environment<IEnvironmentModel>::getCoordinateSystem() const {
return coordinateSystem_;
}
@@ -39,9 +39,9 @@ namespace corsika {
template <typename IEnvironmentModel>
template <class TVolumeType, typename... TVolumeArgs>
std::unique_ptr< VolumeTreeNode<IEnvironmentModel> > Environment<IEnvironmentModel>::createNode(TVolumeArgs&&... args) {
- static_assert(std::is_base_of_v<corsika::Volume, TVolumeType>,
+ static_assert(std::is_base_of_v<IVolume, TVolumeType>,
"unusable type provided, needs to be derived from "
- "\"corsika::Volume\"");
+ "\"Volume\"");
return std::make_unique<BaseNodeType>(
std::make_unique<TVolumeType>(std::forward<TVolumeArgs>(args)...));
diff --git a/corsika/detail/media/FlatExponential.inl b/corsika/detail/media/FlatExponential.inl
index e45e63a2d..0835551c9 100644
--- a/corsika/detail/media/FlatExponential.inl
+++ b/corsika/detail/media/FlatExponential.inl
@@ -23,7 +23,7 @@ namespace corsika {
FlatExponential<T>::FlatExponential(Point const& point,
Vector<dimensionless_d> const& axis,
MassDensityType rho, LengthType lambda,
- NuclearComposition nuclComp)
+ NuclearComposition const& nuclComp)
: BaseExponential<FlatExponential<T>>(point, rho, lambda)
, axis_(axis)
, nuclComp_(nuclComp) {}
diff --git a/corsika/detail/media/HomogeneousMedium.inl b/corsika/detail/media/HomogeneousMedium.inl
index 9ceefa380..1952d5f68 100644
--- a/corsika/detail/media/HomogeneousMedium.inl
+++ b/corsika/detail/media/HomogeneousMedium.inl
@@ -18,7 +18,7 @@ namespace corsika {
template <typename T>
HomogeneousMedium<T>::HomogeneousMedium(MassDensityType density,
- NuclearComposition nuclComp)
+ NuclearComposition const& nuclComp)
: density_(density)
, nuclComp_(nuclComp) {}
@@ -32,14 +32,14 @@ namespace corsika {
}
template <typename T>
- GrammageType HomogeneousMedium<T>::integratedGrammage(
- Trajectory<Line> const&, LengthType to) const {
+ GrammageType HomogeneousMedium<T>::getIntegratedGrammage(Trajectory<Line> const&,
+ LengthType to) const {
return to * density_;
}
template <typename T>
- LengthType HomogeneousMedium<T>::getArclengthFromGrammage(
- Trajectory<Line> const&, GrammageType grammage) const {
+ LengthType HomogeneousMedium<T>::getArclengthFromGrammage(Trajectory<Line> const&,
+ GrammageType grammage) const {
return grammage / density_;
}
} // namespace corsika
diff --git a/corsika/detail/media/LayeredSphericalAtmosphereBuilder.inl b/corsika/detail/media/LayeredSphericalAtmosphereBuilder.inl
index fa6cd992d..76ab3bd70 100644
--- a/corsika/detail/media/LayeredSphericalAtmosphereBuilder.inl
+++ b/corsika/detail/media/LayeredSphericalAtmosphereBuilder.inl
@@ -43,7 +43,7 @@ namespace corsika {
auto const rho0 = b / c;
CORSIKA_LOG_INFO("rho0 = {}, c = {}", rho0, c);
- node->SetModelProperties<SlidingPlanarExponential<IMediumModel>>(
+ node->setModelProperties<SlidingPlanarExponential<IMediumModel>>(
center_, rho0, -c, *composition_, seaLevel_);
layers_.push(std::move(node));
@@ -62,7 +62,7 @@ namespace corsika {
auto node = std::make_unique<VolumeTreeNode<IMediumModel>>(
std::make_unique<Sphere>(center_, radius));
- node->SetModelProperties<HomogeneousMedium<IMediumModel>>(rho0, *composition_);
+ node->setModelProperties<HomogeneousMedium<IMediumModel>>(rho0, *composition_);
layers_.push(std::move(node));
}
@@ -80,7 +80,7 @@ namespace corsika {
while (!layers_.empty()) {
auto l = std::move(layers_.top());
auto* tmp = l.get();
- outmost->AddChild(std::move(l));
+ outmost->addChild(std::move(l));
layers_.pop();
outmost = tmp;
}
diff --git a/corsika/detail/media/LinearApproximationIntegrator.inl b/corsika/detail/media/LinearApproximationIntegrator.inl
index d4e90a922..a3420dcc4 100644
--- a/corsika/detail/media/LinearApproximationIntegrator.inl
+++ b/corsika/detail/media/LinearApproximationIntegrator.inl
@@ -20,20 +20,20 @@ namespace corsika {
}
template <typename TDerived>
- auto LinearApproximationIntegrator<TDerived>::integrateGrammage(
+ auto LinearApproximationIntegrator<TDerived>::getIntegrateGrammage(
Trajectory<Line> const& line, LengthType length) const {
auto const c0 = getImplementation().evaluateAt(line.getPosition(0));
- auto const c1 = getImplementation().rho_.FirstDerivative(line.getPosition(0),
- line.NormalizedDirection());
+ auto const c1 = getImplementation().rho_.getFirstDerivative(
+ line.getPosition(0), line.getNormalizedDirection());
return (c0 + 0.5 * c1 * length) * length;
}
template <typename TDerived>
auto LinearApproximationIntegrator<TDerived>::getArclengthFromGrammage(
Trajectory<Line> const& line, GrammageType grammage) const {
- auto const c0 = getImplementation().rho_(line.GetPosition(0));
- auto const c1 = getImplementation().rho_.FirstDerivative(line.GetPosition(0),
- line.NormalizedDirection());
+ auto const c0 = getImplementation().rho_(line.getPosition(0));
+ auto const c1 = getImplementation().rho_.getFirstDerivative(
+ line.getPosition(0), line.getNormalizedDirection());
return (1 - 0.5 * grammage * c1 / (c0 * c0)) * grammage / c0;
}
@@ -41,8 +41,8 @@ namespace corsika {
template <typename TDerived>
auto LinearApproximationIntegrator<TDerived>::getMaximumLength(
Trajectory<Line> const& line, [[maybe_unused]] double relError) const {
- [[maybe_unused]] auto const c1 = getImplementation().rho_.SecondDerivative(
- line.getPosition(0), line.NormalizedDirection());
+ [[maybe_unused]] auto const c1 = getImplementation().rho_.getSecondDerivative(
+ line.getPosition(0), line.getNormalizedDirection());
// todo: provide a real, working implementation
return 1_m * std::numeric_limits<double>::infinity();
diff --git a/corsika/detail/media/NuclearComposition.inl b/corsika/detail/media/NuclearComposition.inl
index feb4607b6..4e4985bfa 100644
--- a/corsika/detail/media/NuclearComposition.inl
+++ b/corsika/detail/media/NuclearComposition.inl
@@ -11,6 +11,8 @@
#include <corsika/framework/core/ParticleProperties.hpp>
#include <corsika/framework/core/PhysicalUnits.hpp>
+#include <corsika/media/WeightProvider.hpp>
+
#include <cassert>
#include <functional>
#include <numeric>
@@ -20,44 +22,6 @@
namespace corsika {
- template <class AConstIterator, class BConstIterator>
- NuclearComposition::WeightProviderIterator<
- AConstIterator, BConstIterator>::WeightProviderIterator(AConstIterator a,
- BConstIterator b)
- : aIter_(a)
- , bIter_(b) {}
-
- template <class AConstIterator, class BConstIterator>
- typename NuclearComposition::WeightProviderIterator<AConstIterator,
- BConstIterator>::value_type
- NuclearComposition::WeightProviderIterator<AConstIterator, BConstIterator>::operator*()
- const {
- return ((*aIter_) * (*bIter_)).magnitude();
- }
-
- template <class AConstIterator, class BConstIterator>
- NuclearComposition::WeightProviderIterator<AConstIterator, BConstIterator>&
- NuclearComposition::WeightProviderIterator<AConstIterator,
- BConstIterator>::operator++() { // prefix ++
- ++aIter_;
- ++bIter_;
- return *this;
- }
-
- template <class AConstIterator, class BConstIterator>
- auto
- NuclearComposition::WeightProviderIterator<AConstIterator, BConstIterator>::operator==(
- WeightProviderIterator other) {
- return aIter_ == other.aIter_;
- }
-
- template <class AConstIterator, class BConstIterator>
- auto
- NuclearComposition::WeightProviderIterator<AConstIterator, BConstIterator>::operator!=(
- WeightProviderIterator other) {
- return !(*this == other);
- }
-
NuclearComposition::NuclearComposition(std::vector<corsika::Code> const& pComponents,
std::vector<float> const& pFractions)
: numberFractions_(pFractions)
@@ -68,8 +32,7 @@ namespace corsika {
if (is_nucleus(compID)) {
return get_nucleus_A(compID) * fraction;
} else {
- return get_mass(compID) / ConvertSIToHEP(constants::u) *
- fraction;
+ return get_mass(compID) / convert_SI_to_HEP(constants::u) * fraction;
}
})) {
assert(pComponents.size() == pFractions.size());
@@ -83,7 +46,7 @@ namespace corsika {
}
template <typename TFunction>
- auto NuclearComposition::weightedSum(TFunction const& func) const {
+ double NuclearComposition::getWeightedSum(TFunction const& func) const {
using ResultQuantity = decltype(func(*components_.cbegin()));
auto const prod = [&](auto const compID, auto const fraction) {
@@ -103,17 +66,17 @@ namespace corsika {
}
}
- auto NuclearComposition::size() const { return numberFractions_.size(); }
+ size_t NuclearComposition::getSize() const { return numberFractions_.size(); }
std::vector<float> const& NuclearComposition::getFractions() const {
return numberFractions_;
}
-
+
std::vector<corsika::Code> const& NuclearComposition::getComponents() const {
return components_;
}
- auto const NuclearComposition::getAverageMassNumber() const { return avgMassNumber_; }
+ double const NuclearComposition::getAverageMassNumber() const { return avgMassNumber_; }
template <class TRNG>
corsika::Code NuclearComposition::sampleTarget(
@@ -135,7 +98,7 @@ namespace corsika {
// Note: when this class ever modifies its internal data, the hash
// must be updated, too!
- size_t NuclearComposition::hash() const { return hash_; }
+ size_t NuclearComposition::getHash() const { return hash_; }
void NuclearComposition::updateHash() {
std::vector<std::size_t> hashes;
diff --git a/corsika/detail/media/ShowerAxis.inl b/corsika/detail/media/ShowerAxis.inl
index aa1a28906..195d05f95 100644
--- a/corsika/detail/media/ShowerAxis.inl
+++ b/corsika/detail/media/ShowerAxis.inl
@@ -78,10 +78,10 @@ namespace corsika {
assert(0 <= lambda && lambda <= 1.);
- CORSIKA_LOG_TRACE("ShowerAxis::X l={} m, lower={}, lambda={}, upper={}", l / 1_m,
+ CORSIKA_LOG_TRACE("ShowerAxis::getX l={} m, lower={}, lambda={}, upper={}", l / 1_m,
lower, lambda, upper);
- // linear interpolation between X[lower] and X[upper]
+ // linear interpolation between getX[lower] and X[upper]
return X_[upper] * lambda + X_[lower] * (1 - lambda);
}
@@ -91,12 +91,12 @@ namespace corsika {
GrammageType ShowerAxis::getMinimumX() const { return GrammageType::zero(); }
- GrammageType ShowerAxis::projectedX(Point const& p) const {
+ GrammageType ShowerAxis::getProjectedX(Point const& p) const {
auto const projectedLength = (p - pointStart_).dot(axis_normalized_);
- return X(projectedLength);
+ return getX(projectedLength);
}
- corsika::Vector<dimensionless_d> const& ShowerAxis::getDirection() const {
+ Vector<dimensionless_d> const& ShowerAxis::getDirection() const {
return axis_normalized_;
}
diff --git a/corsika/detail/media/SlidingPlanarExponential.inl b/corsika/detail/media/SlidingPlanarExponential.inl
index 9aaa3e7ef..05e8749c8 100644
--- a/corsika/detail/media/SlidingPlanarExponential.inl
+++ b/corsika/detail/media/SlidingPlanarExponential.inl
@@ -17,16 +17,15 @@ namespace corsika {
template <typename T>
SlidingPlanarExponential<T>::SlidingPlanarExponential(
Point const& p0, MassDensityType rho0, LengthType lambda,
- NuclearComposition nuclComp, LengthType referenceHeight)
+ NuclearComposition const& nuclComp, LengthType referenceHeight)
: BaseExponential<SlidingPlanarExponential<T>>(p0, rho0, lambda)
, nuclComp_(nuclComp)
, referenceHeight_(referenceHeight) {}
template <typename T>
- MassDensityType SlidingPlanarExponential<T>::getMassDensity(
- Point const& point) const {
+ MassDensityType SlidingPlanarExponential<T>::getMassDensity(Point const& point) const {
auto const height =
- (point - BaseExponential<SlidingPlanarExponential<T>>::point_).norm() -
+ (point - BaseExponential<SlidingPlanarExponential<T>>::point_).getNorm() -
referenceHeight_;
return BaseExponential<SlidingPlanarExponential<T>>::rho0_ *
exp(BaseExponential<SlidingPlanarExponential<T>>::invLambda_ * height);
@@ -38,20 +37,22 @@ namespace corsika {
}
template <typename T>
- GrammageType SlidingPlanarExponential<T>::integratedGrammage(
+ GrammageType SlidingPlanarExponential<T>::getIntegratedGrammage(
Trajectory<Line> const& line, LengthType l) const {
auto const axis =
- (line.getR0() - BaseExponential<SlidingPlanarExponential<T>>::point_).normalized();
- return BaseExponential<SlidingPlanarExponential<T>>::integratedGrammage(line, l,
- axis);
+ (line.getStartPoint() - BaseExponential<SlidingPlanarExponential<T>>::point_)
+ .normalized();
+ return BaseExponential<SlidingPlanarExponential<T>>::getIntegratedGrammage(line, l,
+ axis);
}
template <typename T>
- LengthType SlidingPlanarExponential<T>::arclengthFromGrammage(
+ LengthType SlidingPlanarExponential<T>::getArclengthFromGrammage(
Trajectory<Line> const& line, GrammageType const grammage) const {
auto const axis =
- (line.getR0() - BaseExponential<SlidingPlanarExponential<T>>::point_).normalized();
- return BaseExponential<SlidingPlanarExponential<T>>::arclengthFromGrammage(
+ (line.getStartPoint() - BaseExponential<SlidingPlanarExponential<T>>::point_)
+ .normalized();
+ return BaseExponential<SlidingPlanarExponential<T>>::getArclengthFromGrammage(
line, grammage, axis);
}
diff --git a/corsika/detail/media/Universe.inl b/corsika/detail/media/Universe.inl
index e54d4715a..09c4c101d 100644
--- a/corsika/detail/media/Universe.inl
+++ b/corsika/detail/media/Universe.inl
@@ -14,10 +14,10 @@
namespace corsika {
- Universe::Universe(corsika::CoordinateSystem const& pCS)
- : corsika::Sphere(corsika::Point{pCS, 0 * meter, 0 * meter, 0 * meter},
+ Universe::Universe(CoordinateSystemPtr const& pCS)
+ : corsika::Sphere(Point{pCS, 0 * meter, 0 * meter, 0 * meter},
meter * std::numeric_limits<double>::infinity()) {}
- bool Universe::Contains(corsika::Point const&) const { return true; }
+ bool Universe::isInside(corsika::Point const&) const { return true; }
} // namespace corsika
\ No newline at end of file
diff --git a/corsika/detail/media/VolumeTreeNode.inl b/corsika/detail/media/VolumeTreeNode.inl
index 79e8b4100..409d393f2 100644
--- a/corsika/detail/media/VolumeTreeNode.inl
+++ b/corsika/detail/media/VolumeTreeNode.inl
@@ -15,20 +15,20 @@
namespace corsika {
- //! convenience function equivalent to Volume::Contains
+ //! convenience function equivalent to Volume::isInside
template <typename IModelProperties>
- bool VolumeTreeNode<IModelProperties>::Contains(corsika::Point const& p) const {
- return fGeoVolume->Contains(p);
+ bool VolumeTreeNode<IModelProperties>::isInside(Point const& p) const {
+ return geoVolume_->isInside(p);
}
template <typename IModelProperties>
inline VolumeTreeNode<IModelProperties> const*
- VolumeTreeNode<IModelProperties>::Excludes(corsika::Point const& p) const {
+ VolumeTreeNode<IModelProperties>::isExcluded(Point const& p) const {
auto exclContainsIter =
- std::find_if(fExcludedNodes.cbegin(), fExcludedNodes.cend(),
- [&](auto const& s) { return bool(s->Contains(p)); });
+ std::find_if(excludedNodes_.cbegin(), excludedNodes_.cend(),
+ [&](auto const& s) { return bool(s->isInside(p)); });
- return exclContainsIter != fExcludedNodes.cend() ? *exclContainsIter : nullptr;
+ return exclContainsIter != excludedNodes_.cend() ? *exclContainsIter : nullptr;
}
/** returns a pointer to the sub-VolumeTreeNode which is "responsible" for the given
@@ -36,22 +36,22 @@ namespace corsika {
*/
template <typename IModelProperties>
VolumeTreeNode<IModelProperties> const*
- VolumeTreeNode<IModelProperties>::GetContainingNode(corsika::Point const& p) const {
- if (!Contains(p)) { return nullptr; }
+ VolumeTreeNode<IModelProperties>::getContainingNode(Point const& p) const {
+ if (!isInside(p)) { return nullptr; }
if (auto const childContainsIter =
- std::find_if(fChildNodes.cbegin(), fChildNodes.cend(),
- [&](auto const& s) { return bool(s->Contains(p)); });
- childContainsIter == fChildNodes.cend()) // not contained in any of the children
+ std::find_if(childNodes_.cbegin(), childNodes_.cend(),
+ [&](auto const& s) { return bool(s->isInside(p)); });
+ childContainsIter == childNodes_.cend()) // not contained in any of the children
{
- if (auto const exclContainsIter = Excludes(p)) // contained in any excluded nodes
+ if (auto const exclContainsIter = isExcluded(p)) // contained in any excluded nodes
{
- return exclContainsIter->GetContainingNode(p);
+ return exclContainsIter->getContainingNode(p);
} else {
return this;
}
} else {
- return (*childContainsIter)->GetContainingNode(p);
+ return (*childContainsIter)->getContainingNode(p);
}
}
@@ -60,29 +60,31 @@ namespace corsika {
void VolumeTreeNode<IModelProperties>::walk(TCallable func) {
if constexpr (preorder) { func(*this); }
- std::for_each(fChildNodes.begin(), fChildNodes.end(),
+ std::for_each(childNodes_.begin(), childNodes_.end(),
[&](auto& v) { v->walk(func); });
if constexpr (!preorder) { func(*this); };
}
template <typename IModelProperties>
- void VolumeTreeNode<IModelProperties>::AddChild(typename VolumeTreeNode<IModelProperties>::VTNUPtr pChild) {
- pChild->fParentNode = this;
- fChildNodes.push_back(std::move(pChild));
+ void VolumeTreeNode<IModelProperties>::addChild(
+ typename VolumeTreeNode<IModelProperties>::VTNUPtr pChild) {
+ pChild->parentNode_ = this;
+ childNodes_.push_back(std::move(pChild));
// It is a bad idea to return an iterator to the inserted element
// because it might get invalidated when the vector needs to grow
// later and the caller won't notice.
}
template <typename IModelProperties>
- void VolumeTreeNode<IModelProperties>::ExcludeOverlapWith(typename VolumeTreeNode<IModelProperties>::VTNUPtr const& pNode) {
- fExcludedNodes.push_back(pNode.get());
+ void VolumeTreeNode<IModelProperties>::excludeOverlapWith(
+ typename VolumeTreeNode<IModelProperties>::VTNUPtr const& pNode) {
+ excludedNodes_.push_back(pNode.get());
}
template <typename IModelProperties>
template <class MediumType, typename... Args>
- auto VolumeTreeNode<IModelProperties>::CreateMedium(Args&&... args) {
+ auto VolumeTreeNode<IModelProperties>::createMedium(Args&&... args) {
static_assert(std::is_base_of_v<IMediumModel, MediumType>,
"unusable type provided, needs to be derived from \"IMediumModel\"");
diff --git a/corsika/detail/modules/energy_loss/BetheBlochPDG.inl b/corsika/detail/modules/energy_loss/BetheBlochPDG.inl
index 28ed8be58..da0827cc9 100644
--- a/corsika/detail/modules/energy_loss/BetheBlochPDG.inl
+++ b/corsika/detail/modules/energy_loss/BetheBlochPDG.inl
@@ -161,7 +161,7 @@ namespace corsika::energy_loss {
if (p.GetChargeNumber() == 0) return corsika::EProcessReturn::eOk;
GrammageType const dX =
- p.GetNode()->GetModelProperties().integratedGrammage(t, t.GetLength());
+ p.GetNode()->GetModelProperties().getIntegratedGrammage(t, t.GetLength());
std::cout << "BetheBlochPDG " << p.GetPID() << ", z=" << p.GetChargeNumber()
<< ", dX=" << dX / 1_g * square(1_cm) << "g/cm2" << std::endl;
HEPEnergyType dE = TotalEnergyLoss(p, dX);
@@ -196,8 +196,9 @@ namespace corsika::energy_loss {
auto const maxLoss = 0.01 * vParticle.GetEnergy();
auto const maxGrammage = maxLoss / dE * dX;
- return vParticle.GetNode()->GetModelProperties().arclengthFromGrammage(vTrack,
- maxGrammage) *
+ return vParticle.GetNode()->GetModelProperties().getArclengthFromGrammage(
+ vTrack,
+ maxGrammage) *
1.0001; // to make sure particle gets absorbed when DoContinuous() is called
}
@@ -223,11 +224,11 @@ namespace corsika::energy_loss {
SetupTrack const trajToEndBin(lineToEndBin, 1_s);
GrammageType const grammageStart =
- vP.GetNode()->GetModelProperties().integratedGrammage(trajToStartBin,
- trajToStartBin.GetLength());
+ vP.GetNode()->GetModelProperties().getIntegratedGrammage(
+ trajToStartBin, trajToStartBin.GetLength());
GrammageType const grammageEnd =
- vP.GetNode()->GetModelProperties().integratedGrammage(trajToEndBin,
- trajToEndBin.GetLength());
+ vP.GetNode()->GetModelProperties().getIntegratedGrammage(
+ trajToEndBin, trajToEndBin.GetLength());
const int binStart = grammageStart / dX_;
const int binEnd = grammageEnd / dX_;
diff --git a/corsika/media/BaseExponential.hpp b/corsika/media/BaseExponential.hpp
index 2a85a0b4f..b0afa123a 100644
--- a/corsika/media/BaseExponential.hpp
+++ b/corsika/media/BaseExponential.hpp
@@ -46,7 +46,7 @@ namespace corsika {
* \f]
*/
// clang-format on
- units::si::GrammageType integratedGrammage(
+ units::si::GrammageType getIntegratedGrammage(
Trajectory<Line> const& line, units::si::LengthType vL,
Vector<units::si::dimensionless_d> const& axis) const;
@@ -68,7 +68,7 @@ namespace corsika {
* \f]
*/
// clang-format on
- units::si::LengthType arclengthFromGrammage(
+ units::si::LengthType getArclengthFromGrammage(
Trajectory<Line> const& line, units::si::GrammageType grammage,
Vector<units::si::dimensionless_d> const& axis) const;
@@ -78,6 +78,6 @@ namespace corsika {
}; // class BaseExponential
- } // namespace corsika
+} // namespace corsika
#include <corsika/detail/media/BaseExponential.inl>
diff --git a/corsika/media/Environment.hpp b/corsika/media/Environment.hpp
index f9cf69d82..4a8ddb791 100644
--- a/corsika/media/Environment.hpp
+++ b/corsika/media/Environment.hpp
@@ -45,7 +45,7 @@ namespace corsika {
*
* @retval Retuns a const reference to the CoordinateSystem used
**/
- CoordinateSystem const& getCoordinateSystem() const;
+ CoordinateSystemPtr const& getCoordinateSystem() const;
/** Factory method for creation of VolumeTreeNodes
* @tparam TVolumeType Type of volume to be created
@@ -57,7 +57,7 @@ namespace corsika {
static std::unique_ptr<BaseNodeType> createNode(TVolumeArgs&&... args);
private:
- CoordinateSystem const& coordinateSystem_;
+ CoordinateSystemPtr const coordinateSystem_;
typename BaseNodeType::VTNUPtr universe_;
};
diff --git a/corsika/media/FlatExponential.hpp b/corsika/media/FlatExponential.hpp
index acb456d12..01f16105c 100644
--- a/corsika/media/FlatExponential.hpp
+++ b/corsika/media/FlatExponential.hpp
@@ -30,25 +30,24 @@ namespace corsika {
// clang-format on
template <typename T>
class FlatExponential : public BaseExponential<FlatExponential<T>>, public T {
- Vector<units::si::dimensionless_d> const axis_;
+ Vector<dimensionless_d> const axis_;
NuclearComposition const nuclComp_;
using Base = BaseExponential<FlatExponential<T>>;
public:
- FlatExponential(Point const& point, Vector<units::si::dimensionless_d> const& axis,
- units::si::MassDensityType rho, units::si::LengthType lambda,
- NuclearComposition nuclComp);
+ FlatExponential(Point const& point, Vector<dimensionless_d> const& axis,
+ MassDensityType rho, LengthType lambda,
+ NuclearComposition const& nuclComp);
- units::si::MassDensityType getMassDensity(Point const& point) const override;
+ MassDensityType getMassDensity(Point const& point) const override;
NuclearComposition const& getNuclearComposition() const override;
- units::si::GrammageType integratedGrammage(Trajectory<Line> const& line,
- units::si::LengthType to) const;
+ GrammageType getIntegratedGrammage(Trajectory<Line> const& line, LengthType to) const;
- units::si::LengthType arclengthFromGrammage(Trajectory<Line> const& line,
- units::si::GrammageType grammage) const;
+ LengthType getArclengthFromGrammage(Trajectory<Line> const& line,
+ GrammageType grammage) const;
};
} // namespace corsika
diff --git a/corsika/media/HomogeneousMedium.hpp b/corsika/media/HomogeneousMedium.hpp
index 56ea8a25d..ee79dcb6f 100644
--- a/corsika/media/HomogeneousMedium.hpp
+++ b/corsika/media/HomogeneousMedium.hpp
@@ -22,23 +22,23 @@ namespace corsika {
template <typename T>
class HomogeneousMedium : public T {
- units::si::MassDensityType const density_;
- NuclearComposition const nuclComp_;
public:
- HomogeneousMedium(units::si::MassDensityType density, NuclearComposition nuclComp);
+ HomogeneousMedium(MassDensityType density, NuclearComposition const& nuclComp);
- units::si::MassDensityType getMassDensity(Point const&) const override;
+ MassDensityType getMassDensity(Point const&) const override;
NuclearComposition const& getNuclearComposition() const override;
- units::si::GrammageType integratedGrammage(
- Trajectory<Line> const&,
- units::si::LengthType to) const override;
+ GrammageType getIntegratedGrammage(Trajectory<Line> const&,
+ LengthType to) const override;
- units::si::LengthType arclengthFromGrammage(
- Trajectory<Line> const&,
- units::si::GrammageType grammage) const override;
+ LengthType getArclengthFromGrammage(Trajectory<Line> const&,
+ GrammageType grammage) const override;
+
+ private:
+ MassDensityType const density_;
+ NuclearComposition const nuclComp_;
};
} // namespace corsika
diff --git a/corsika/media/IMagneticFieldModel.hpp b/corsika/media/IMagneticFieldModel.hpp
index f110feb97..91741085f 100644
--- a/corsika/media/IMagneticFieldModel.hpp
+++ b/corsika/media/IMagneticFieldModel.hpp
@@ -25,7 +25,7 @@ namespace corsika {
// a type-alias for a magnetic field vector
using MagneticFieldVector =
- corsika::geometry::Vector<corsika::phys::units::si::magnetic_flux_density_d>;
+ Vector<magnetic_flux_density_d>;
public:
/**
@@ -34,7 +34,7 @@ namespace corsika {
* @param point The location to evaluate the field at.
* @returns The magnetic field vector at that point.
*/
- virtual auto getMagneticField(corsika::geometry::Point const&) const
+ virtual auto getMagneticField(Point const&) const
-> MagneticFieldVector = 0;
/**
diff --git a/corsika/media/IMediumModel.hpp b/corsika/media/IMediumModel.hpp
index c63819251..7fcaafd12 100644
--- a/corsika/media/IMediumModel.hpp
+++ b/corsika/media/IMediumModel.hpp
@@ -20,15 +20,15 @@ namespace corsika {
public:
virtual ~IMediumModel() = default; // LCOV_EXCL_LINE
- virtual units::si::MassDensityType getMassDensity(Point const&) const = 0;
+ virtual 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 units::si::GrammageType integratedGrammage(
- Trajectory<Line> const&, units::si::LengthType) const = 0;
+ virtual GrammageType getIntegratedGrammage(Trajectory<Line> const&,
+ LengthType) const = 0;
- virtual units::si::LengthType arclengthFromGrammage(
- Trajectory<Line> const&, units::si::GrammageType) const = 0;
+ virtual LengthType getArclengthFromGrammage(Trajectory<Line> const&,
+ GrammageType) const = 0;
virtual NuclearComposition const& getNuclearComposition() const = 0;
};
diff --git a/corsika/media/InhomogeneousMedium.hpp b/corsika/media/InhomogeneousMedium.hpp
index 5467e2db3..615fe8910 100644
--- a/corsika/media/InhomogeneousMedium.hpp
+++ b/corsika/media/InhomogeneousMedium.hpp
@@ -23,22 +23,24 @@ namespace corsika {
template <typename T, typename TDensityFunction>
class InhomogeneousMedium : public T {
- NuclearComposition const nuclComp_;
- TDensityFunction const densityFunction_;
public:
template <typename... TArgs>
- InhomogeneousMedium(NuclearComposition nuclComp, TArgs&&... rhoTArgs);
+ InhomogeneousMedium(NuclearComposition const& nuclComp, TArgs&&... rhoTArgs);
- units::si::MassDensityType getMassDensity(Point const& point) const override;
+ MassDensityType getMassDensity(Point const& point) const override;
NuclearComposition const& getNuclearComposition() const override;
- units::si::GrammageType integratedGrammage(Trajectory<Line> const& line,
- units::si::LengthType to) const override;
+ GrammageType getIntegratedGrammage(Trajectory<Line> const& line,
+ LengthType to) const override;
+
+ LengthType getArclengthFromGrammage(Trajectory<Line> const& pLine,
+ GrammageType grammage) const override;
- units::si::LengthType arclengthFromGrammage(
- Trajectory<Line> const& pLine, units::si::GrammageType grammage) const override;
+ private:
+ NuclearComposition const nuclComp_;
+ TDensityFunction const densityFunction_;
};
} // namespace corsika
diff --git a/corsika/media/LayeredSphericalAtmosphereBuilder.hpp b/corsika/media/LayeredSphericalAtmosphereBuilder.hpp
index 5984913ec..60da93956 100644
--- a/corsika/media/LayeredSphericalAtmosphereBuilder.hpp
+++ b/corsika/media/LayeredSphericalAtmosphereBuilder.hpp
@@ -46,7 +46,7 @@ namespace corsika {
void addExponentialLayer(GrammageType, LengthType, LengthType);
- auto size() const { return layers_.size(); }
+ size_t getSize() const { return layers_.size(); }
void addLinearLayer(LengthType, LengthType);
diff --git a/corsika/media/LinearApproximationIntegrator.hpp b/corsika/media/LinearApproximationIntegrator.hpp
index ac49879dd..e14e0e9b3 100644
--- a/corsika/media/LinearApproximationIntegrator.hpp
+++ b/corsika/media/LinearApproximationIntegrator.hpp
@@ -20,14 +20,13 @@ namespace corsika {
auto const& getImplementation() const;
public:
- auto integrateGrammage(Trajectory<Line> const& line,
- units::si::LengthType length) const;
+ auto getIntegrateGrammage(Trajectory<Line> const& line, LengthType length) const;
- auto arclengthFromGrammage(Trajectory<Line> const& line,
- units::si::GrammageType grammage) const;
+ auto getArclengthFromGrammage(Trajectory<Line> const& line,
+ GrammageType grammage) const;
- auto maximumLength(Trajectory<Line> const& line,
- [[maybe_unused]] double relError) const;
+ auto getMaximumLength(Trajectory<Line> const& line,
+ [[maybe_unused]] double relError) const;
};
} // namespace corsika
diff --git a/corsika/media/NuclearComposition.hpp b/corsika/media/NuclearComposition.hpp
index 513ebef7e..c3d2585f6 100644
--- a/corsika/media/NuclearComposition.hpp
+++ b/corsika/media/NuclearComposition.hpp
@@ -24,39 +24,6 @@ namespace corsika {
* Allowes and handles the creation of custom matter compositions
**/
class NuclearComposition {
- private:
- /// TODO: replace with
- /// https://www.boost.org/doc/libs/1_74_0/libs/iterator/doc/zip_iterator.html or
- /// ranges zip
- /** Double Iterator
- * Iterator that allowes the iteration of two individual lists at the same time. The
- *user needs to take care that booth lists have the same length.
- * @tparam AConstIterator Iterator Type of the first list
- * @tparam BConstIterator Iterator Type of the second list
- **/
- template <class AConstIterator, class BConstIterator>
- class WeightProviderIterator {
- AConstIterator aIter_;
- BConstIterator bIter_;
-
- public:
- using value_type = double;
- using iterator_category = std::input_iterator_tag;
- using pointer = value_type*;
- using reference = value_type&;
- using difference_type = ptrdiff_t;
-
- WeightProviderIterator(AConstIterator a, BConstIterator b);
-
- value_type operator*() const;
-
- WeightProviderIterator& operator++();
-
- auto operator==(WeightProviderIterator other);
-
- auto operator!=(WeightProviderIterator other);
- };
-
public:
/** Constructor
* The constructore takes a list of elements and a list which describe the relative
@@ -66,8 +33,8 @@ namespace corsika {
* @param pFractions List of fractions how much each particle contributes. The sum
*needs to add up to 1
**/
- NuclearComposition(std::vector<corsika::Code> pComponents,
- std::vector<float> pFractions);
+ NuclearComposition(std::vector<corsika::Code> const& pComponents,
+ std::vector<float> const& pFractions);
/** Sum all all relative composition weighted by func(element)
* This function sums all relative compositions given during this classes
@@ -79,18 +46,18 @@ namespace corsika {
* @retval returns the weighted sum with the type defined by the return type of func
**/
template <typename TFunction>
- auto weightedSum(TFunction func) const;
+ double getWeightedSum(TFunction const& func) const;
/** Number of elements in the composition array
* @retval returns the number of elements in the composition array
**/
- auto size() const;
+ size_t getSize() const;
/// Returns a const reference to the fraction
std::vector<float> const& getFractions() const;
/// Returns a const reference to the fraction
std::vector<corsika::Code> const& getComponents() const;
- auto const getAverageMassNumber() const;
+ double const getAverageMassNumber() const;
template <class TRNG>
Code sampleTarget(std::vector<units::si::CrossSectionType> const& sigma,
@@ -98,7 +65,7 @@ namespace corsika {
// Note: when this class ever modifies its internal data, the hash
// must be updated, too!
- size_t hash() const;
+ size_t getHash() const;
private:
void updateHash();
diff --git a/corsika/media/ShowerAxis.hpp b/corsika/media/ShowerAxis.hpp
index b58bc9128..a2fea4cd8 100644
--- a/corsika/media/ShowerAxis.hpp
+++ b/corsika/media/ShowerAxis.hpp
@@ -40,42 +40,38 @@ namespace corsika {
*
**/
- ///\todo documentation needs update ...
+ ///\todo documentation needs update ...
class ShowerAxis {
public:
template <typename TEnvModel>
- ShowerAxis(corsika::Point const& pStart,
- corsika::Vector<units::si::length_d> length,
+ ShowerAxis(Point const& pStart, Vector<length_d> length,
Environment<TEnvModel> const& env, int steps = 10'000);
- units::si::LengthType steplength() const;
+ LengthType getSteplength() const;
- units::si::GrammageType maximumX() const;
+ GrammageType getMaximumX() const;
- units::si::GrammageType minimumX() const;
+ GrammageType getMinimumX() const;
- units::si::GrammageType projectedX(Point const& p) const;
+ GrammageType getProjectedX(Point const& p) const;
GrammageType getX(LengthType) const;
- Vector<units::si::dimensionless_d> const& getDirection() const;
+ Vector<dimensionless_d> const& getDirection() const;
Point const& getStart() const;
private:
Point const pointStart_;
- Vector<units::si::length_d> const length_;
- units::si::LengthType const max_length_, steplength_;
- Vector<units::si::dimensionless_d> const axis_normalized_;
- std::vector<units::si::GrammageType> X_;
+ Vector<length_d> const length_;
+ LengthType const max_length_, steplength_;
+ Vector<dimensionless_d> const axis_normalized_;
+ std::vector<GrammageType> X_;
// for storing the lengths corresponding to equidistant X values
- units::si::GrammageType const X_binning_ = std::invoke([]() {
- using namespace units::si;
- return 1_g / 1_cm / 1_cm;
- });
- std::vector<units::si::LengthType> d_;
+ GrammageType const X_binning_ = 1_g / 1_cm / 1_cm;
+ std::vector<LengthType> d_;
};
} // namespace corsika
-#include <corsika/detail/media/ShowerAxis.inl>
\ No newline at end of file
+#include <corsika/detail/media/ShowerAxis.inl>
diff --git a/corsika/media/SlidingPlanarExponential.hpp b/corsika/media/SlidingPlanarExponential.hpp
index 013135491..66e51eea7 100644
--- a/corsika/media/SlidingPlanarExponential.hpp
+++ b/corsika/media/SlidingPlanarExponential.hpp
@@ -34,26 +34,27 @@ namespace corsika {
template <typename T>
class SlidingPlanarExponential : public BaseExponential<SlidingPlanarExponential<T>>,
public T {
- NuclearComposition const nuclComp_;
- units::si::LengthType const referenceHeight_;
using Base = BaseExponential<SlidingPlanarExponential<T>>;
public:
- SlidingPlanarExponential(
- Point const& p0, units::si::MassDensityType rho0, units::si::LengthType lambda,
- NuclearComposition nuclComp,
- units::si::LengthType referenceHeight = units::si::LengthType::zero());
+ SlidingPlanarExponential(Point const& p0, MassDensityType rho0, LengthType lambda,
+ NuclearComposition const& nuclComp,
+ LengthType referenceHeight = LengthType::zero());
- units::si::MassDensityType getMassDensity(Point const& point) const override;
+ MassDensityType getMassDensity(Point const& point) const override;
NuclearComposition const& getNuclearComposition() const override;
- units::si::GrammageType integratedGrammage(Trajectory<Line> const& line,
- units::si::LengthType l) const override;
+ GrammageType getIntegratedGrammage(Trajectory<Line> const& line,
+ LengthType l) const override;
+
+ LengthType getArclengthFromGrammage(Trajectory<Line> const& line,
+ GrammageType grammage) const override;
- units::si::LengthType arclengthFromGrammage(
- Trajectory<Line> const& line, units::si::GrammageType grammage) const override;
+ private:
+ NuclearComposition const nuclComp_;
+ LengthType const referenceHeight_;
};
} // namespace corsika
diff --git a/corsika/media/Universe.hpp b/corsika/media/Universe.hpp
index 753526edd..af1ead8cb 100644
--- a/corsika/media/Universe.hpp
+++ b/corsika/media/Universe.hpp
@@ -14,8 +14,8 @@
namespace corsika {
struct Universe : public corsika::Sphere {
- inline Universe(corsika::CoordinateSystem const& pCS);
- inline bool Contains(corsika::Point const&) const override;
+ inline Universe(corsika::CoordinateSystemPtr const& pCS);
+ inline bool isInside(corsika::Point const&) const override;
};
} // namespace corsika
diff --git a/corsika/media/VolumeTreeNode.hpp b/corsika/media/VolumeTreeNode.hpp
index 051f8c96c..98ba1d442 100644
--- a/corsika/media/VolumeTreeNode.hpp
+++ b/corsika/media/VolumeTreeNode.hpp
@@ -25,22 +25,21 @@ namespace corsika {
using VTN_type = VolumeTreeNode<IModelProperties>;
using VTNUPtr = std::unique_ptr<VolumeTreeNode<IModelProperties>>;
using IMPSharedPtr = std::shared_ptr<IModelProperties>;
- using VolUPtr = std::unique_ptr<corsika::Volume>;
+ using VolUPtr = std::unique_ptr<IVolume>;
VolumeTreeNode(VolUPtr pVolume = nullptr)
- : fGeoVolume(std::move(pVolume)) {}
+ : geoVolume_(std::move(pVolume)) {}
- //! convenience function equivalent to Volume::Contains
- inline bool Contains(corsika::Point const& p) const;
+ //! convenience function equivalent to Volume::isInside
+ inline bool isInside(Point const& p) const;
- inline VolumeTreeNode<IModelProperties> const* Excludes(
- corsika::Point const& p) const;
+ inline VolumeTreeNode<IModelProperties> const* isExcluded(Point const& p) const;
/** returns a pointer to the sub-VolumeTreeNode which is "responsible" for the given
* \class Point \p p, or nullptr iff \p p is not contained in this volume.
*/
- inline VolumeTreeNode<IModelProperties> const* GetContainingNode(
- corsika::Point const& p) const;
+ inline VolumeTreeNode<IModelProperties> const* getContainingNode(
+ Point const& p) const;
/**
* Traverses the VolumeTree pre- or post-order and calls the functor \p func for each
@@ -50,43 +49,43 @@ namespace corsika {
template <typename TCallable, bool preorder = true>
inline void walk(TCallable func);
- inline void AddChild(VTNUPtr pChild);
+ inline void addChild(VTNUPtr pChild);
- inline void ExcludeOverlapWith(VTNUPtr const& pNode);
+ inline void excludeOverlapWith(VTNUPtr const& pNode);
- inline auto* GetParent() const { return fParentNode; };
+ inline auto* getParent() const { return parentNode_; };
- inline auto const& GetChildNodes() const { return fChildNodes; }
+ inline auto const& getChildNodes() const { return childNodes_; }
- inline auto const& GetExcludedNodes() const { return fExcludedNodes; }
+ inline auto const& getExcludedNodes() const { return excludedNodes_; }
- inline auto const& GetVolume() const { return *fGeoVolume; }
+ inline auto const& getVolume() const { return *geoVolume_; }
- inline auto const& GetModelProperties() const { return *fModelProperties; }
+ inline auto const& getModelProperties() const { return *modelProperties_; }
- inline bool HasModelProperties() const { return fModelProperties.get() != nullptr; }
+ inline bool hasModelProperties() const { return modelProperties_.get() != nullptr; }
template <typename ModelProperties, typename... Args>
- inline auto SetModelProperties(Args&&... args) {
+ inline auto setModelProperties(Args&&... args) {
static_assert(std::is_base_of_v<IModelProperties, ModelProperties>,
"unusable model properties type provided");
- fModelProperties = std::make_shared<ModelProperties>(std::forward<Args>(args)...);
- return fModelProperties;
+ modelProperties_ = std::make_shared<ModelProperties>(std::forward<Args>(args)...);
+ return modelProperties_;
}
- inline void SetModelProperties(IMPSharedPtr ptr) { fModelProperties = ptr; }
+ inline void setModelProperties(IMPSharedPtr ptr) { modelProperties_ = ptr; }
/*
template <class MediumType, typename... Args>
- static auto CreateMedium(Args&&... args);
+ static auto createMedium(Args&&... args);
private:
- std::vector<VTNUPtr> fChildNodes;
- std::vector<VolumeTreeNode<IModelProperties> const*> fExcludedNodes;
- VolumeTreeNode<IModelProperties> const* fParentNode = nullptr;
- VolUPtr fGeoVolume;
- IMPSharedPtr fModelProperties;
+ std::vector<VTNUPtr> childNodes_;
+ std::vector<VolumeTreeNode<IModelProperties> const*> excludedNodes_;
+ VolumeTreeNode<IModelProperties> const* parentNode_ = nullptr;
+ VolUPtr geoVolume_;
+ IMPSharedPtr modelProperties_;
};
} // namespace corsika
diff --git a/tests/media/testEnvironment.cpp b/tests/media/testEnvironment.cpp
index 989958a45..43ed6524f 100644
--- a/tests/media/testEnvironment.cpp
+++ b/tests/media/testEnvironment.cpp
@@ -27,10 +27,8 @@
#include <catch2/catch.hpp>
using namespace corsika;
-using namespace corsika::units::si;
-CoordinateSystem const& gCS =
- RootCoordinateSystem::getInstance().GetRootCoordinateSystem();
+CoordinateSystemPtr const& gCS = get_root_CoordinateSystem();
Point const gOrigin(gCS, {0_m, 0_m, 0_m});
@@ -56,8 +54,8 @@ TEST_CASE("FlatExponential") {
gCS, {20_cm / second, 0_m / second, 0_m / second}));
Trajectory<Line> const trajectory(line, tEnd);
- REQUIRE((medium.integratedGrammage(trajectory, 2_m) / (rho0 * 2_m)) == Approx(1));
- REQUIRE((medium.arclengthFromGrammage(trajectory, rho0 * 5_m) / 5_m) == Approx(1));
+ CHECK((medium.getIntegratedGrammage(trajectory, 2_m) / (rho0 * 2_m)) == Approx(1));
+ CHECK((medium.getArclengthFromGrammage(trajectory, rho0 * 5_m) / 5_m) == Approx(1));
}
SECTION("vertical") {
@@ -67,8 +65,8 @@ TEST_CASE("FlatExponential") {
LengthType const length = 2 * lambda;
GrammageType const exact = rho0 * lambda * (exp(length / lambda) - 1);
- REQUIRE((medium.integratedGrammage(trajectory, length) / exact) == Approx(1));
- REQUIRE((medium.arclengthFromGrammage(trajectory, exact) / length) == Approx(1));
+ CHECK((medium.getIntegratedGrammage(trajectory, length) / exact) == Approx(1));
+ CHECK((medium.getArclengthFromGrammage(trajectory, exact) / length) == Approx(1));
}
SECTION("escape grammage") {
@@ -78,9 +76,9 @@ TEST_CASE("FlatExponential") {
GrammageType const escapeGrammage = rho0 * lambda;
- REQUIRE(trajectory.NormalizedDirection().dot(axis).magnitude() < 0);
- REQUIRE(medium.arclengthFromGrammage(trajectory, 1.2 * escapeGrammage) ==
- std::numeric_limits<typename GrammageType::value_type>::infinity() * 1_m);
+ CHECK(trajectory.getNormalizedDirection().dot(axis).magnitude() < 0);
+ CHECK(medium.getArclengthFromGrammage(trajectory, 1.2 * escapeGrammage) ==
+ std::numeric_limits<typename GrammageType::value_type>::infinity() * 1_m);
}
SECTION("inclined") {
@@ -91,8 +89,8 @@ TEST_CASE("FlatExponential") {
LengthType const length = 2 * lambda;
GrammageType const exact =
rho0 * lambda * (exp(cosTheta * length / lambda) - 1) / cosTheta;
- REQUIRE((medium.integratedGrammage(trajectory, length) / exact) == Approx(1));
- REQUIRE((medium.arclengthFromGrammage(trajectory, exact) / length) == Approx(1));
+ CHECK((medium.getIntegratedGrammage(trajectory, length) / exact) == Approx(1));
+ CHECK((medium.getArclengthFromGrammage(trajectory, exact) / length) == Approx(1));
}
}
@@ -124,31 +122,31 @@ TEST_CASE("SlidingPlanarExponential") {
CHECK(medium.getMassDensity({gCS, {0_mm, 0_m, 3_m}}).magnitude() ==
flat.getMassDensity({gCS, {0_mm, 0_m, 3_m}}).magnitude());
- CHECK(medium.integratedGrammage(trajectory, 2_m).magnitude() ==
- flat.integratedGrammage(trajectory, 2_m).magnitude());
- CHECK(medium.arclengthFromGrammage(trajectory, rho0 * 5_m).magnitude() ==
- flat.arclengthFromGrammage(trajectory, rho0 * 5_m).magnitude());
+ CHECK(medium.getIntegratedGrammage(trajectory, 2_m).magnitude() ==
+ flat.getIntegratedGrammage(trajectory, 2_m).magnitude());
+ CHECK(medium.getArclengthFromGrammage(trajectory, rho0 * 5_m).magnitude() ==
+ flat.getArclengthFromGrammage(trajectory, rho0 * 5_m).magnitude());
}
}
auto constexpr rho0 = 1_kg / 1_m / 1_m / 1_m;
struct Exponential {
- auto operator()(corsika::Point const& p) const {
- return exp(p.GetCoordinates()[0] / 1_m) * rho0;
+ auto operator()(Point const& p) const {
+ return exp(p.getCoordinates()[0] / 1_m) * rho0;
}
template <int N>
- auto Derivative(Point const& p, Vector<dimensionless_d> const& v) const {
- return v.GetComponents()[0] * (*this)(p) / static_pow<N>(1_m);
+ auto getDerivative(Point const& p, Vector<dimensionless_d> const& v) const {
+ return v.getComponents()[0] * (*this)(p) / static_pow<N>(1_m);
}
- auto FirstDerivative(Point const& p, Vector<dimensionless_d> const& v) const {
- return Derivative<1>(p, v);
+ auto getFirstDerivative(Point const& p, Vector<dimensionless_d> const& v) const {
+ return getDerivative<1>(p, v);
}
- auto SecondDerivative(Point const& p, Vector<dimensionless_d> const& v) const {
- return Derivative<2>(p, v);
+ auto getSecondDerivative(Point const& p, Vector<dimensionless_d> const& v) const {
+ return getDerivative<2>(p, v);
}
};
@@ -165,9 +163,9 @@ TEST_CASE("InhomogeneousMedium") {
DensityFunction<decltype(e), LinearApproximationIntegrator> const rho(e);
SECTION("DensityFunction") {
- REQUIRE(e.Derivative<1>(gOrigin, direction) / (1_kg / 1_m / 1_m / 1_m / 1_m) ==
- Approx(1));
- REQUIRE(rho.evaluateAt(gOrigin) == e(gOrigin));
+ CHECK(e.getDerivative<1>(gOrigin, direction) / (1_kg / 1_m / 1_m / 1_m / 1_m) ==
+ Approx(1));
+ CHECK(rho.evaluateAt(gOrigin) == e(gOrigin));
}
auto const exactGrammage = [](auto l) { return 1_m * rho0 * (exp(l / 1_m) - 1); };
@@ -179,18 +177,18 @@ TEST_CASE("InhomogeneousMedium") {
InhomogeneousMedium<IMediumModel, decltype(rho)> const inhMedium(composition, rho);
SECTION("Integration") {
- REQUIRE(rho.integrateGrammage(trajectory, l) / exactGrammage(l) ==
- Approx(1).epsilon(1e-2));
- REQUIRE(rho.arclengthFromGrammage(trajectory, exactGrammage(l)) /
- exactLength(exactGrammage(l)) ==
- Approx(1).epsilon(1e-2));
- REQUIRE(rho.maximumLength(trajectory, 1e-2) >
- l); // todo: write reasonable test when implementation is working
-
- REQUIRE(rho.integrateGrammage(trajectory, l) ==
- inhMedium.integratedGrammage(trajectory, l));
- REQUIRE(rho.arclengthFromGrammage(trajectory, 20_g / (1_cm * 1_cm)) ==
- inhMedium.arclengthFromGrammage(trajectory, 20_g / (1_cm * 1_cm)));
+ CHECK(rho.getIntegrateGrammage(trajectory, l) / exactGrammage(l) ==
+ Approx(1).epsilon(1e-2));
+ CHECK(rho.getArclengthFromGrammage(trajectory, exactGrammage(l)) /
+ exactLength(exactGrammage(l)) ==
+ Approx(1).epsilon(1e-2));
+ CHECK(rho.getMaximumLength(trajectory, 1e-2) >
+ l); // todo: write reasonable test when implementation is working
+
+ CHECK(rho.getIntegrateGrammage(trajectory, l) ==
+ inhMedium.getIntegratedGrammage(trajectory, l));
+ CHECK(rho.getArclengthFromGrammage(trajectory, 20_g / (1_cm * 1_cm)) ==
+ inhMedium.getArclengthFromGrammage(trajectory, 20_g / (1_cm * 1_cm)));
}
}
@@ -202,264 +200,26 @@ TEST_CASE("LayeredSphericalAtmosphereBuilder") {
builder.addLinearLayer(2_km, 20_km);
builder.addExponentialLayer(540.1778_g / (1_cm * 1_cm), 772170.16_cm, 30_km);
- CHECK(builder.size() == 3);
-
- auto const builtEnv = builder.assemble();
- auto const& univ = builtEnv.GetUniverse();
-
- CHECK(builder.size() == 0);
-
- auto const R = builder.getEarthRadius();
-
- CHECK(univ->GetChildNodes().size() == 1);
-
- CHECK(univ->GetContainingNode(Point(gCS, 0_m, 0_m, R + 35_km)) == univ.get());
- CHECK(dynamic_cast<Sphere const&>(
- univ->GetContainingNode(Point(gCS, 0_m, 0_m, R + 8_km))->GetVolume())
- .GetRadius() == R + 10_km);
- CHECK(dynamic_cast<Sphere const&>(
- univ->GetContainingNode(Point(gCS, 0_m, 0_m, R + 12_km))->GetVolume())
- .GetRadius() == R + 20_km);
- CHECK(dynamic_cast<Sphere const&>(
- univ->GetContainingNode(Point(gCS, 0_m, 0_m, R + 24_km))->GetVolume())
- .GetRadius() == R + 30_km);
-}
-
-TEST_CASE("UniformMagneticField w/ Homogeneous Medium") {
-
- // setup our interface types
- using IModelInterface = IMagneticFieldModel<IMediumModel>;
- using AtmModel = UniformMagneticField<HomogeneousMedium<IModelInterface>>;
-
- // the composition we use for the homogenous medium
- NuclearComposition const protonComposition(std::vector<Code>{Code::Proton},
- std::vector<float>{1.f});
-
- // create a magnetic field vector
- Vector B0(gCS, 0_T, 0_T, 0_T);
-
- // the constant density
- const auto density{19.2_g / cube(1_cm)};
-
- // create our atmospheric model
- AtmModel medium(B0, density, protonComposition);
-
- // and test at several locations
- CHECK(B0.GetComponents(gCS) ==
- medium.GetMagneticField(Point(gCS, -10_m, 4_m, 35_km)).GetComponents(gCS));
- CHECK(
- B0.GetComponents(gCS) ==
- medium.GetMagneticField(Point(gCS, 1000_km, -1000_km, 1000_km)).GetComponents(gCS));
- CHECK(B0.GetComponents(gCS) ==
- medium.GetMagneticField(Point(gCS, 0_m, 0_m, 0_m)).GetComponents(gCS));
-
- // create a new magnetic field vector
- Vector B1(gCS, 23_T, 57_T, -4_T);
-
- // and update this atmospheric model
- medium.SetMagneticField(B1);
-
- // and test at several locations
- CHECK(B1.GetComponents(gCS) ==
- medium.GetMagneticField(Point(gCS, -10_m, 4_m, 35_km)).GetComponents(gCS));
- CHECK(
- B1.GetComponents(gCS) ==
- medium.GetMagneticField(Point(gCS, 1000_km, -1000_km, 1000_km)).GetComponents(gCS));
- CHECK(B1.GetComponents(gCS) ==
- medium.GetMagneticField(Point(gCS, 0_m, 0_m, 0_m)).GetComponents(gCS));
-
- // 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
- setup::Trajectory const trajectory = setup::testing::make_track<setup::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));
-}
-
-TEST_CASE("LayeredSphericalAtmosphereBuilder w/ magnetic field") {
-
- // setup our interface types
- using ModelInterface = IMagneticFieldModel<IMediumModel>;
-
- // the composition we use for the homogenous medium
- NuclearComposition const protonComposition(std::vector<Code>{Code::Proton},
- std::vector<float>{1.f});
-
- // create magnetic field vectors
- Vector B0(gCS, 0_T, 0_T, 1_T);
-
- LayeredSphericalAtmosphereBuilder builder =
- environment::make_layered_spherical_atmosphere_builder<
- ModelInterface,
- UniformMagneticField>::create(gOrigin, units::constants::EarthRadius::Mean, B0);
-
- builder.setNuclearComposition(
- {{{particles::Code::Nitrogen, particles::Code::Oxygen}}, {{.6, .4}}});
- builder.addLinearLayer(1_km, 10_km);
- builder.addExponentialLayer(1222.6562_g / (1_cm * 1_cm), 994186.38_cm, 20_km);
-
- CHECK(builder.size() == 2);
+ CHECK(builder.getSize() == 3);
auto const builtEnv = builder.assemble();
auto const& univ = builtEnv.getUniverse();
- CHECK(builder.size() == 0);
- CHECK(univ->GetChildNodes().size() == 1);
- auto const R = builder.getEarthRadius();
-
- // check magnetic field at several locations
- const Point pTest(gCS, -10_m, 4_m, R + 35_m);
- CHECK(B0.GetComponents(gCS) == univ->GetContainingNode(pTest)
- ->GetModelProperties()
- .GetMagneticField(pTest)
- .GetComponents(gCS));
- const Point pTest2(gCS, 10_m, -4_m, R + 15_km);
- CHECK(B0.GetComponents(gCS) == univ->GetContainingNode(pTest2)
- ->GetModelProperties()
- .GetMagneticField(pTest2)
- .GetComponents(gCS));
-}
-
-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);
-
- // and require that it is constant
- CHECK(n == medium.GetRefractiveIndex(Point(gCS, -10_m, 4_m, 35_km)));
- CHECK(n == medium.GetRefractiveIndex(Point(gCS, +210_m, 0_m, 7_km)));
- CHECK(n == medium.GetRefractiveIndex(Point(gCS, 0_m, 0_m, 0_km)));
- CHECK(n == medium.GetRefractiveIndex(Point(gCS, 100_km, 400_km, 350_km)));
-
- // a new refractive index
- const double n2{2.3472123};
-
- // update the refractive index of this atmospheric model
- medium.SetRefractiveIndex(n2);
-
- // check that the returned refractive index is correct
- CHECK(n2 == medium.GetRefractiveIndex(Point(gCS, -10_m, 4_m, 35_km)));
- CHECK(n2 == medium.GetRefractiveIndex(Point(gCS, +210_m, 0_m, 7_km)));
- CHECK(n2 == medium.GetRefractiveIndex(Point(gCS, 0_m, 0_m, 0_km)));
- CHECK(n2 == medium.GetRefractiveIndex(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}));
+ CHECK(builder.getSize() == 0);
// the end time of our line
auto const tEnd = 1_s;
- // and the associated trajectory
- setup::Trajectory const trajectory = setup::testing::make_track<setup::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));
-}
-
-TEST_CASE("MediumProperties") {
-
- // test access of medium properties via enum and class types
-
- const Medium type = Medium::AirDry1Atm;
- const MediumData& air = mediumData(type);
- CHECK(air.Ieff() == 85.7);
- CHECK(air.Cbar() == 10.5961);
- CHECK(air.x0() == 1.7418);
- CHECK(air.x1() == 4.2759);
- CHECK(air.aa() == 0.10914);
- CHECK(air.sk() == 3.3994);
- CHECK(air.dlt0() == 0.0);
-}
-
-TEST_CASE("MediumPropertyModel w/ Homogeneous") {
-
- // setup our interface types
- using IModelInterface = IMediumPropertyModel<IMediumModel>;
- using AtmModel = MediumPropertyModel<HomogeneousMedium<IModelInterface>>;
+ CHECK(univ->getChildNodes().size() == 1);
- // 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
- setup::Trajectory const trajectory = setup::testing::make_track<setup::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));
+ CHECK(univ->getContainingNode(Point(gCS, 0_m, 0_m, R + 35_km)) == univ.get());
+ CHECK(dynamic_cast<Sphere const&>(
+ univ->getContainingNode(Point(gCS, 0_m, 0_m, R + 8_km))->getVolume())
+ .getRadius() == R + 10_km);
+ CHECK(dynamic_cast<Sphere const&>(
+ univ->getContainingNode(Point(gCS, 0_m, 0_m, R + 12_km))->getVolume())
+ .getRadius() == R + 20_km);
+ CHECK(dynamic_cast<Sphere const&>(
+ univ->getContainingNode(Point(gCS, 0_m, 0_m, R + 24_km))->getVolume())
+ .getRadius() == R + 30_km);
}
diff --git a/tests/media/testMedium.cpp b/tests/media/testMedium.cpp
index d9dcf0c86..d658d7c31 100644
--- a/tests/media/testMedium.cpp
+++ b/tests/media/testMedium.cpp
@@ -25,8 +25,7 @@ using namespace corsika::units::si;
TEST_CASE("MediumPropertyModel w/ Homogeneous") {
- CoordinateSystem const& gCS =
- RootCoordinateSystem::getInstance().GetRootCoordinateSystem();
+ CoordinateSystemPtr gCS = get_root_CoordinateSystem();
Point const gOrigin(gCS, {0_m, 0_m, 0_m});
@@ -83,6 +82,6 @@ TEST_CASE("MediumPropertyModel w/ Homogeneous") {
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));
+ CHECK((medium.getIntegratedGrammage(trajectory, 3_m) / (density * 3_m)) == Approx(1));
+ CHECK((medium.getArclengthFromGrammage(trajectory, density * 5_m) / 5_m) == Approx(1));
}
diff --git a/tests/media/testRefractiveIndex.cpp b/tests/media/testRefractiveIndex.cpp
index 79215e7ac..65f086c88 100644
--- a/tests/media/testRefractiveIndex.cpp
+++ b/tests/media/testRefractiveIndex.cpp
@@ -24,8 +24,7 @@ using namespace corsika::units::si;
TEST_CASE("UniformRefractiveIndex w/ Homogeneous") {
- CoordinateSystem const& gCS =
- RootCoordinateSystem::getInstance().GetRootCoordinateSystem();
+ CoordinateSystemPtr const& gCS = get_root_CoordinateSystem();
Point const gOrigin(gCS, {0_m, 0_m, 0_m});
@@ -82,6 +81,6 @@ TEST_CASE("UniformRefractiveIndex w/ Homogeneous") {
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));
+ CHECK((medium.getIntegratedGrammage(trajectory, 3_m) / (density * 3_m)) == Approx(1));
+ CHECK((medium.getArclengthFromGrammage(trajectory, density * 5_m) / 5_m) == Approx(1));
}
diff --git a/tests/media/testShowerAxis.cpp b/tests/media/testShowerAxis.cpp
index fcf2b920d..cba73d6cc 100644
--- a/tests/media/testShowerAxis.cpp
+++ b/tests/media/testShowerAxis.cpp
@@ -20,7 +20,7 @@
#include <catch2/catch.hpp>
-//using namespace
+// using namespace
using namespace corsika;
const auto density = 1_kg / (1_m * 1_m * 1_m);
@@ -31,15 +31,13 @@ auto setupEnvironment(Code vTargetCode) {
auto& universe = *(env->getUniverse());
const CoordinateSystem& cs = env->getCoordinateSystem();
- auto theMedium =
- Environment<IMediumModel>::createNode<Sphere>(
- Point{cs, 0_m, 0_m, 0_m},
- 1_km * std::numeric_limits<double>::infinity());
+ auto theMedium = Environment<IMediumModel>::createNode<Sphere>(
+ Point{cs, 0_m, 0_m, 0_m}, 1_km * std::numeric_limits<double>::infinity());
using MyHomogeneousModel = HomogeneousMedium<IMediumModel>;
theMedium->SetModelProperties<MyHomogeneousModel>(
- density, NuclearComposition(std::vector<Code>{vTargetCode},
- std::vector<float>{1.}));
+ density,
+ NuclearComposition(std::vector<Code>{vTargetCode}, std::vector<float>{1.}));
auto const* nodePtr = theMedium.get();
universe.AddChild(std::move(theMedium));
@@ -59,22 +57,21 @@ TEST_CASE("Homogeneous Density") {
Point const showerCore{cs, 0_m, 0_m, observationHeight};
Point const injectionPos = showerCore + Vector<dimensionless_d>{cs, {0, 0, 1}} * t;
- ShowerAxis const showerAxis{injectionPos, (showerCore - injectionPos),
- *env, 20};
+ ShowerAxis const showerAxis{injectionPos, (showerCore - injectionPos), *env, 20};
- CHECK(showerAxis.steplength() == 500_m);
+ CHECK(showerAxis.getSteplength() == 500_m);
- CHECK(showerAxis.maximumX() / (10_km * density) == Approx(1).epsilon(1e-8));
+ CHECK(showerAxis.getMaximumX() / (10_km * density) == Approx(1).epsilon(1e-8));
- CHECK(showerAxis.minimumX() == 0_g / square(1_cm));
+ CHECK(showerAxis.getMinimumX() == 0_g / square(1_cm));
const Point p{cs, 10_km, 20_km, 8.3_km};
- CHECK(showerAxis.projectedX(p) / (1.7_km * density) == Approx(1).epsilon(1e-8));
+ CHECK(showerAxis.getProjectedX(p) / (1.7_km * density) == Approx(1).epsilon(1e-8));
- const units::si::LengthType d = 6.789_km;
- CHECK(showerAxis.X(d) / (d * density) == Approx(1).epsilon(1e-8));
+ const LengthType d = 6.789_km;
+ CHECK(showerAxis.getX(d) / (d * density) == Approx(1).epsilon(1e-8));
const Vector<dimensionless_d> dir{cs, {0, 0, -1}};
- CHECK(showerAxis.getDirection().GetComponents(cs) == dir.GetComponents(cs));
- CHECK(showerAxis.getStart().GetCoordinates() == injectionPos.GetCoordinates());
+ CHECK(showerAxis.getDirection().getComponents(cs) == dir.getComponents(cs));
+ CHECK(showerAxis.getStart().getCoordinates() == injectionPos.getCoordinates());
}
--
GitLab