IAP GITLAB

Skip to content
Snippets Groups Projects
Commit 03c803c6 authored by ralfulrich's avatar ralfulrich
Browse files

reworked Trajectory interface

parent e60dab45
No related branches found
No related tags found
No related merge requests found
...@@ -12,11 +12,11 @@ ...@@ -12,11 +12,11 @@
#ifndef _include_Cascade_h_ #ifndef _include_Cascade_h_
#define _include_Cascade_h_ #define _include_Cascade_h_
#include <corsika/geometry/LineTrajectory.h> // to be removed. for dummy trajectory only
#include <corsika/geometry/Point.h> // to be removed. for dummy trajectory only
#include <corsika/process/ProcessReturn.h> #include <corsika/process/ProcessReturn.h>
#include <corsika/units/PhysicalUnits.h> #include <corsika/units/PhysicalUnits.h>
#include <corsika/setup/SetupTrajectory.h>
using namespace corsika::units::si; using namespace corsika::units::si;
namespace corsika::cascade { namespace corsika::cascade {
...@@ -62,16 +62,17 @@ namespace corsika::cascade { ...@@ -62,16 +62,17 @@ namespace corsika::cascade {
// DoCascadeEquations(); // // DoCascadeEquations(); //
} }
} }
void Step(Particle& particle) { void Step(Particle& particle) {
[[maybe_unused]] double nextStep = fProcesseList.MinStepLength(particle); [[maybe_unused]] double nextStep = fProcesseList.MinStepLength(particle);
// corsika::utls::ignore(nextStep); // corsika::utls::ignore(nextStep);
auto const root = corsika::geometry::CoordinateSystem::CreateRootCS(); auto const root = corsika::geometry::CoordinateSystem::CreateRootCS();
corsika::geometry::LineTrajectory corsika::geometry::Trajectory<corsika::geometry::Line>
trajectory( // trajectory is not yet used. this is a dummy. trajectory( // trajectory is not yet used. this is a dummy.
corsika::geometry::Point(root, {0_m, 0_m, 0_m}), corsika::geometry::Line(corsika::geometry::Point(root, {0_m, 0_m, 0_m}),
corsika::geometry::Vector<corsika::units::si::SpeedType::dimension_type>( corsika::geometry::Vector<corsika::units::si::SpeedType::dimension_type>(
root, 0 * 1_m / second, 0 * 1_m / second, 1 * 1_m / second)); root, 0 * 1_m / second, 0 * 1_m / second, 1 * 1_m / second)),
0_s, 1_s);
corsika::process::EProcessReturn status = corsika::process::EProcessReturn status =
fProcesseList.DoContinuous(particle, trajectory, fStack); fProcesseList.DoContinuous(particle, trajectory, fStack);
if (status == corsika::process::EProcessReturn::eParticleAbsorbed) { if (status == corsika::process::EProcessReturn::eParticleAbsorbed) {
......
...@@ -10,7 +10,6 @@ ...@@ -10,7 +10,6 @@
*/ */
#include <corsika/cascade/Cascade.h> #include <corsika/cascade/Cascade.h>
#include <corsika/geometry/LineTrajectory.h>
#include <corsika/process/ProcessSequence.h> #include <corsika/process/ProcessSequence.h>
#include <corsika/process/stack_inspector/StackInspector.h> #include <corsika/process/stack_inspector/StackInspector.h>
......
...@@ -24,17 +24,37 @@ namespace corsika::geometry { ...@@ -24,17 +24,37 @@ namespace corsika::geometry {
*/ */
class BaseTrajectory { class BaseTrajectory {
BaseTrajectory() = delete;
public: public:
BaseTrajectory(corsika::units::si::TimeType start, corsika::units::si::TimeType end)
: fTStart(start)
, fTEnd(end) {}
//!< for \f$ t = 0 \f$, the starting Point shall be returned. //!< for \f$ t = 0 \f$, the starting Point shall be returned.
virtual Point GetPosition(corsika::units::si::TimeType) const = 0; virtual Point GetPosition(corsika::units::si::TimeType) const = 0;
//!< the Point is return from u=0 (start) to u=1 (end)
virtual Point GetPosition(double u) const = 0;
/*! /*!
* returns the arc length between two points of the trajectory * returns the length between two points of the trajectory
* parameterized by \arg t1 and \arg t2. Requires \arg t2 > \arg t1. * parameterized by \arg t1 and \arg t2. Requires \arg t2 > \arg t1.
*/ */
virtual LengthType DistanceBetween(corsika::units::si::TimeType t1, virtual LengthType GetDistance(corsika::units::si::TimeType t1,
corsika::units::si::TimeType t2) const = 0; corsika::units::si::TimeType t2) const = 0;
virtual corsika::units::si::TimeType GetDuration(
corsika::units::si::TimeType t1, corsika::units::si::TimeType t2) const {
return t2 - t1;
}
virtual Point GetEndpoint() const { return GetPosition(fTEnd); }
virtual Point GetStartpoint() const { return GetPosition(fTStart); }
protected:
corsika::units::si::TimeType const fTStart, fTEnd;
}; };
} // namespace corsika::geometry } // namespace corsika::geometry
......
...@@ -8,14 +8,14 @@ set ( ...@@ -8,14 +8,14 @@ set (
GEOMETRY_HEADERS GEOMETRY_HEADERS
Vector.h Vector.h
Point.h Point.h
Line.h
Sphere.h Sphere.h
CoordinateSystem.h CoordinateSystem.h
Helix.h Helix.h
BaseVector.h BaseVector.h
QuantityVector.h QuantityVector.h
BaseTrajectory.h
LineTrajectory.h
Trajectory.h Trajectory.h
BaseTrajectory.h
) )
set ( set (
......
...@@ -12,7 +12,6 @@ ...@@ -12,7 +12,6 @@
#ifndef _include_HELIX_H_ #ifndef _include_HELIX_H_
#define _include_HELIX_H_ #define _include_HELIX_H_
#include <corsika/geometry/BaseTrajectory.h>
#include <corsika/geometry/Point.h> #include <corsika/geometry/Point.h>
#include <corsika/geometry/Vector.h> #include <corsika/geometry/Vector.h>
#include <corsika/units/PhysicalUnits.h> #include <corsika/units/PhysicalUnits.h>
...@@ -31,7 +30,8 @@ namespace corsika::geometry { ...@@ -31,7 +30,8 @@ namespace corsika::geometry {
\f} \f}
*/ */
class Helix : public BaseTrajectory { class Helix {
using VelocityVec = Vector<corsika::units::si::SpeedType::dimension_type>; using VelocityVec = Vector<corsika::units::si::SpeedType::dimension_type>;
Point const r0; Point const r0;
...@@ -51,7 +51,7 @@ namespace corsika::geometry { ...@@ -51,7 +51,7 @@ namespace corsika::geometry {
, uPerp(vPerp.cross(vPar.normalized())) , uPerp(vPerp.cross(vPar.normalized()))
, radius(pvPar.norm() / abs(pOmegaC)) {} , radius(pvPar.norm() / abs(pOmegaC)) {}
Point GetPosition(corsika::units::si::TimeType t) const override { Point GetPosition(corsika::units::si::TimeType t) const {
return r0 + vPar * t + return r0 + vPar * t +
(vPerp * (cos(omegaC * t) - 1) + uPerp * sin(omegaC * t)) / omegaC; (vPerp * (cos(omegaC * t) - 1) + uPerp * sin(omegaC * t)) / omegaC;
} }
...@@ -59,7 +59,7 @@ namespace corsika::geometry { ...@@ -59,7 +59,7 @@ namespace corsika::geometry {
auto GetRadius() const { return radius; } auto GetRadius() const { return radius; }
LengthType DistanceBetween(corsika::units::si::TimeType t1, LengthType DistanceBetween(corsika::units::si::TimeType t1,
corsika::units::si::TimeType t2) const override { corsika::units::si::TimeType t2) const {
return (vPar + vPerp).norm() * (t2 - t1); return (vPar + vPerp).norm() * (t2 - t1);
} }
}; };
......
...@@ -12,30 +12,30 @@ ...@@ -12,30 +12,30 @@
#ifndef _include_LINETRAJECTORY_H #ifndef _include_LINETRAJECTORY_H
#define _include_LINETRAJECTORY_H #define _include_LINETRAJECTORY_H
#include <corsika/geometry/BaseTrajectory.h>
#include <corsika/geometry/Point.h> #include <corsika/geometry/Point.h>
#include <corsika/geometry/Vector.h> #include <corsika/geometry/Vector.h>
#include <corsika/units/PhysicalUnits.h> #include <corsika/units/PhysicalUnits.h>
namespace corsika::geometry { namespace corsika::geometry {
class LineTrajectory : public BaseTrajectory { class Line {
using VelocityVec = Vector<corsika::units::si::SpeedType::dimension_type>; using VelocityVec = Vector<corsika::units::si::SpeedType::dimension_type>;
Point const r0; Point const r0;
VelocityVec const v0; VelocityVec const v0;
public: public:
LineTrajectory(Point const& pR0, VelocityVec const& pV0) Line(Point const& pR0, VelocityVec const& pV0)
: r0(pR0) : r0(pR0)
, v0(pV0) {} , v0(pV0) {}
Point GetPosition(corsika::units::si::TimeType t) const override { Point GetPosition(corsika::units::si::TimeType t) const {
return r0 + v0 * t; return r0 + v0 * t;
} }
LengthType DistanceBetween(corsika::units::si::TimeType t1, LengthType DistanceBetween(corsika::units::si::TimeType t1,
corsika::units::si::TimeType t2) const override { corsika::units::si::TimeType t2) const {
assert(t2 >= t1); assert(t2 >= t1);
return v0.norm() * (t2 - t1); return v0.norm() * (t2 - t1);
} }
......
...@@ -17,28 +17,32 @@ ...@@ -17,28 +17,32 @@
namespace corsika::geometry { namespace corsika::geometry {
class Trajectory { template <typename T>
corsika::units::si::TimeType const fTStart, fTEnd; class Trajectory : public BaseTrajectory {
BaseTrajectory const& fTrajectory;
T fTraj;
public: public:
Trajectory(corsika::units::si::TimeType pTStart, corsika::units::si::TimeType pTEnd, Trajectory(T const& theT, corsika::units::si::TimeType pTStart,
BaseTrajectory const& pTrajectory) corsika::units::si::TimeType pTEnd)
: fTStart(pTStart) //: T(theT), fTStart(pTStart), fTEnd(pTEnd) {}
, fTEnd(pTEnd) : BaseTrajectory(pTStart, pTEnd)
, fTrajectory(pTrajectory) {} , fTraj(theT) {}
Point GetPosition(corsika::units::si::TimeType t) const { Point GetPosition(corsika::units::si::TimeType t) const {
return fTrajectory.GetPosition(t + fTStart); return fTraj.GetPosition(t + fTStart);
} }
Point GetPosition(double u) const { Point GetPosition(double u) const {
return GetPosition(fTEnd * u + fTStart * (1 - u)); return GetPosition(fTEnd * u + fTStart * (1 - u));
} }
LengthType GetDistance(corsika::units::si::TimeType t1,
corsika::units::si::TimeType t2) const {
return fTraj.DistanceBetween(t1, t2);
}
auto GetEndpoint() const { return GetPosition(fTEnd); }
auto GetStartpoint() const { return GetPosition(fTStart); }
}; };
} // namespace corsika::geometry } // namespace corsika::geometry
......
...@@ -15,7 +15,7 @@ ...@@ -15,7 +15,7 @@
#include <corsika/geometry/CoordinateSystem.h> #include <corsika/geometry/CoordinateSystem.h>
#include <corsika/geometry/Helix.h> #include <corsika/geometry/Helix.h>
#include <corsika/geometry/LineTrajectory.h> #include <corsika/geometry/Line.h>
#include <corsika/geometry/Point.h> #include <corsika/geometry/Point.h>
#include <corsika/geometry/Sphere.h> #include <corsika/geometry/Sphere.h>
#include <corsika/geometry/Trajectory.h> #include <corsika/geometry/Trajectory.h>
...@@ -144,17 +144,17 @@ TEST_CASE("Trajectories") { ...@@ -144,17 +144,17 @@ TEST_CASE("Trajectories") {
Vector<SpeedType::dimension_type> v0(rootCS, Vector<SpeedType::dimension_type> v0(rootCS,
{1_m / second, 0_m / second, 0_m / second}); {1_m / second, 0_m / second, 0_m / second});
LineTrajectory const lineTrajectory(r0, v0); Line const line(r0, v0);
CHECK((lineTrajectory.GetPosition(2_s).GetCoordinates() - CHECK((line.GetPosition(2_s).GetCoordinates() -
QuantityVector<length_d>(2_m, 0_m, 0_m)) QuantityVector<length_d>(2_m, 0_m, 0_m))
.norm() .norm()
.magnitude() == Approx(0).margin(absMargin)); .magnitude() == Approx(0).margin(absMargin));
BaseTrajectory const* base = &lineTrajectory; Trajectory<Line> base(line, 0_s, 1_s);
CHECK(lineTrajectory.GetPosition(2_s).GetCoordinates() == CHECK(line.GetPosition(2_s).GetCoordinates() ==
base->GetPosition(2_s).GetCoordinates()); base.GetPosition(2_s).GetCoordinates());
CHECK(base->DistanceBetween(1_s, 2_s) / 1_m == Approx(1)); CHECK(base.GetDistance(1_s, 2_s) / 1_m == Approx(1));
} }
SECTION("Helix") { SECTION("Helix") {
...@@ -178,10 +178,10 @@ TEST_CASE("Trajectories") { ...@@ -178,10 +178,10 @@ TEST_CASE("Trajectories") {
.norm() .norm()
.magnitude() == Approx(0).margin(absMargin)); .magnitude() == Approx(0).margin(absMargin));
BaseTrajectory const* base = &helix; Trajectory<Helix> const base(helix, 0_s, 1_s);
CHECK(helix.GetPosition(1234_s).GetCoordinates() == CHECK(helix.GetPosition(1234_s).GetCoordinates() ==
base->GetPosition(1234_s).GetCoordinates()); base.GetPosition(1234_s).GetCoordinates());
CHECK(base->DistanceBetween(1_s, 2_s) / 1_m == Approx(5)); CHECK(base.GetDistance(1_s, 2_s) / 1_m == Approx(5));
} }
} }
...@@ -12,11 +12,12 @@ ...@@ -12,11 +12,12 @@
#ifndef _corsika_setup_setuptrajectory_h_ #ifndef _corsika_setup_setuptrajectory_h_
#define _corsika_setup_setuptrajectory_h_ #define _corsika_setup_setuptrajectory_h_
#include <corsika/geometry/LineTrajectory.h> #include <corsika/geometry/Line.h>
#include <corsika/geometry/Trajectory.h>
namespace corsika::setup { namespace corsika::setup {
typedef corsika::geometry::LineTrajectory Trajectory; typedef corsika::geometry::Trajectory<corsika::geometry::Line> Trajectory;
} }
#endif #endif
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