diff --git a/corsika/detail/framework/geometry/Helix.inl b/corsika/detail/framework/geometry/Helix.inl
index f023dec0b31713d91437aeffbbed0ce2e433e2a2..1f9d3f94f1b815a3f2fe0420951fd4a7160de430 100644
--- a/corsika/detail/framework/geometry/Helix.inl
+++ b/corsika/detail/framework/geometry/Helix.inl
@@ -17,23 +17,24 @@
namespace corsika {
- Point Helix::GetPosition(TimeType t) const {
- return r0 + vPar * t +
- (vPerp * (cos(omegaC * t) - 1) + uPerp * sin(omegaC * t)) / omegaC;
- }
+ LengthType Helix::getRadius() const { return radius_; }
- Point Helix::PositionFromArclength(LengthType l) const {
- return GetPosition(TimeFromArclength(l));
+ Point Helix::getPosition(TimeType t) const {
+ return r0_ + vPar_ * t +
+ (vPerp_ * (std::cos(omegaC_ * t) - 1) + uPerp_ * std::sin(omegaC_ * t)) /
+ omegaC_;
}
- LengthType Helix::GetRadius() const { return radius; }
+ Point Helix::getPositionFromArclength(LengthType l) const {
+ return getPosition(getTimeFromArclength(l));
+ }
- LengthType Helix::ArcLength(TimeType t1, TimeType t2) const {
- return (vPar + vPerp).norm() * (t2 - t1);
+ LengthType Helix::getArcLength(TimeType t1, TimeType t2) const {
+ return (vPar_ + vPerp_).norm() * (t2 - t1);
}
- TimeType Helix::TimeFromArclength(LengthType l) const {
- return l / (vPar + vPerp).norm();
+ TimeType Helix::getTimeFromArclength(LengthType l) const {
+ return l / (vPar_ + vPerp_).norm();
}
} // namespace corsika
diff --git a/corsika/framework/geometry/Helix.hpp b/corsika/framework/geometry/Helix.hpp
index 85e1e92354064b9d5070a0f51fceeaf921d8b967..55b071586bd8c7d62d4bb5819143d519074eacde 100644
--- a/corsika/framework/geometry/Helix.hpp
+++ b/corsika/framework/geometry/Helix.hpp
@@ -14,8 +14,9 @@ n/*
#include <corsika/framework/geometry/Vector.hpp>
namespace corsika {
+
/*!
- * \class Helix
+ * Defines a helical path
*
* A Helix is defined by the cyclotron frequency \f$ \omega_c \f$, the initial
* Point r0 and
@@ -30,34 +31,36 @@ namespace corsika {
class Helix {
- using VelocityVec = Vector<SpeedType::dimension_type> ;
-
- Point const r0;
- FrequencyType const omegaC;
- VelocityVec const vPar;
- VelocityVec const vPerp, uPerp;
-
- LengthType const radius;
+ ///! \todo move VelocityVec into PhysicalUnits
+ using VelocityVec = Vector<SpeedType::dimension_type>;
public:
Helix(Point const& pR0, FrequencyType pOmegaC, VelocityVec const& pvPar,
VelocityVec const& pvPerp)
- : r0(pR0)
- , omegaC(pOmegaC)
- , vPar(pvPar)
- , vPerp(pvPerp)
- , uPerp(vPerp.cross(vPar.normalized()))
- , radius(pvPar.norm() / abs(pOmegaC)) {}
+ : r0_(pR0)
+ , omegaC_(pOmegaC)
+ , vPar_(pvPar)
+ , vPerp_(pvPerp)
+ , uPerp_(vPerp_.cross(vPar_.normalized()))
+ , radius_(pvPar.norm() / abs(pOmegaC)) {}
- inline Point GetPosition(TimeType t) const;
+ inline LengthType getRadius() const;
- inline Point PositionFromArclength(LengthType l) const;
+ inline Point getPosition(TimeType t) const;
+
+ inline Point getPositionFromArclength(LengthType l) const;
- inline LengthType GetRadius() const;
+ inline LengthType getArcLength(TimeType t1, TimeType t2) const;
- inline LengthType ArcLength(TimeType t1, TimeType t2) const;
+ inline TimeType getTimeFromArclength(LengthType l) const;
- inline TimeType TimeFromArclength(LengthType l) const;
+ private:
+ Point r0_; ///! origin of helix, but this is in the center of the
+ ///! "cylinder" on which the helix rotates
+ FrequencyType omegaC_; ///! speed of angular rotation
+ VelocityVec vPar_; ///! speed along direction of "cylinder"
+ VelocityVec vPerp_, uPerp_;
+ LengthType radius_;
};
} // namespace corsika
diff --git a/tests/framework/CMakeLists.txt b/tests/framework/CMakeLists.txt
index 1eb9c3a75d1c3176920b6ff39acb583b93f70795..81bab8f2273995fa69c6586b9d5a2e78b5e90889 100644
--- a/tests/framework/CMakeLists.txt
+++ b/tests/framework/CMakeLists.txt
@@ -6,6 +6,7 @@ set (test_framework_sources
#testCOMBoost.cpp #FIXME: reenable this
#testCorsikaFenv.cpp # does not work because of use of exceptions in catch2
testFourVector.cpp
+ testHelix.cpp
testFunctionTimer.cpp
testGeometry.cpp
testLogging.cpp
diff --git a/tests/framework/testGeometry.cpp b/tests/framework/testGeometry.cpp
index dbffa55cf3ce0caa3beae0ac9b33df3c54f30dfc..b3eb21e070229726e97bd8f2aeb62e76274e681b 100644
--- a/tests/framework/testGeometry.cpp
+++ b/tests/framework/testGeometry.cpp
@@ -11,7 +11,6 @@
#include <cmath>
#include <corsika/framework/core/PhysicalUnits.hpp>
#include <corsika/framework/geometry/CoordinateSystem.hpp>
-#include <corsika/framework/geometry/Helix.hpp>
#include <corsika/framework/geometry/Line.hpp>
#include <corsika/framework/geometry/Point.hpp>
#include <corsika/framework/geometry/RootCoordinateSystem.hpp>
@@ -19,7 +18,6 @@
#include <corsika/framework/geometry/Trajectory.hpp>
using namespace corsika;
-using namespace corsika::units::si;
double constexpr absMargin = 1.0e-8;
@@ -233,31 +231,4 @@ TEST_CASE("Trajectories") {
.eVector.norm() == Approx(0).margin(absMargin));
}
- SECTION("Helix") {
- Vector<SpeedType::dimension_type> const vPar(
- rootCS, {0_m / second, 0_m / second, 4_m / second});
-
- Vector<SpeedType::dimension_type> const vPerp(
- rootCS, {3_m / second, 0_m / second, 0_m / second});
-
- auto const T = 1_s;
- auto const omegaC = 2 * M_PI / T;
-
- Helix const helix(r0, omegaC, vPar, vPerp);
-
- CHECK((helix.GetPosition(1_s).GetCoordinates() -
- QuantityVector<length_d>(0_m, 0_m, 4_m))
- .norm()
- .magnitude() == Approx(0).margin(absMargin));
-
- CHECK((helix.GetPosition(0.25_s).GetCoordinates() -
- QuantityVector<length_d>(-3_m / (2 * M_PI), -3_m / (2 * M_PI), 1_m))
- .norm()
- .magnitude() == Approx(0).margin(absMargin));
-
- CHECK(
- (helix.GetPosition(7_s) - helix.PositionFromArclength(helix.ArcLength(0_s, 7_s)))
- .norm()
- .magnitude() == Approx(0).margin(absMargin));
- }
}
diff --git a/tests/framework/testHelix.cpp b/tests/framework/testHelix.cpp
new file mode 100644
index 0000000000000000000000000000000000000000..74338e8ab27f2a80c44afd868f912a777f904fd1
--- /dev/null
+++ b/tests/framework/testHelix.cpp
@@ -0,0 +1,66 @@
+/*
+ * (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 <catch2/catch.hpp>
+
+#include <corsika/framework/geometry/Helix.hpp>
+#include <corsika/framework/core/PhysicalUnits.hpp>
+#include <corsika/framework/geometry/CoordinateSystem.hpp>
+#include <corsika/framework/geometry/Line.hpp>
+#include <corsika/framework/geometry/Point.hpp>
+#include <corsika/framework/geometry/RootCoordinateSystem.hpp>
+
+using namespace corsika;
+
+double constexpr absMargin = 1.0e-8;
+
+TEST_CASE("Helix class") {
+
+ CoordinateSystem& rootCS =
+ RootCoordinateSystem::getInstance().GetRootCoordinateSystem();
+ Point r0(rootCS, {0_m, 0_m, 0_m});
+
+ SECTION("Helix") {
+ Vector<SpeedType::dimension_type> const vPar(
+ rootCS, {0_m / second, 0_m / second, 4_m / second});
+
+ Vector<SpeedType::dimension_type> const vPerp(
+ rootCS, {3_m / second, 0_m / second, 0_m / second});
+
+ auto const T = 1_s;
+ auto const omegaC = 2 * M_PI / T;
+
+ Helix const helix(r0, omegaC, vPar, vPerp);
+
+ CHECK((helix.getPosition(1_s).GetCoordinates() -
+ QuantityVector<length_d>(0_m, 0_m, 4_m))
+ .norm()
+ .magnitude() == Approx(0).margin(absMargin));
+
+ CHECK((helix.getPosition(0.25_s).GetCoordinates() -
+ QuantityVector<length_d>(-3_m / (2 * M_PI), -3_m / (2 * M_PI), 1_m))
+ .norm()
+ .magnitude() == Approx(0).margin(absMargin));
+
+ CHECK(
+ (helix.getPosition(7_s) - helix.getPositionFromArclength(helix.getArcLength(0_s, 7_s)))
+ .norm()
+ .magnitude() == Approx(0).margin(absMargin));
+
+ /*
+ // we have to consider this, if we need it
+ auto const t = 1234_s;
+ Trajectory<Helix> const base(helix, t);
+ CHECK(helix.getPosition(t).GetCoordinates() == base.GetPosition(1.).GetCoordinates());
+
+ CHECK(base.ArcLength(0_s, 1_s) / 1_m == Approx(5));
+ */
+ }
+
+}
+