IAP GITLAB

Skip to content
Snippets Groups Projects
Commit f79ed358 authored by Fan Hu's avatar Fan Hu Committed by ralfulrich
Browse files

fix logical bugs

parent b4d74c69
No related branches found
No related tags found
No related merge requests found
......@@ -18,7 +18,7 @@ ObservationCubic<TTracking, TOutput>::ObservationCubic(
template <typename TTracking, typename TOutput>
template <typename TParticle, typename TTrajectory>
inline ProcessReturn ObservationCubic<TTracking, TOutput>::doContinuous(
TParticle &particle, TTrajectory &, bool const stepLimit) {
TParticle &particle, TTrajectory & step, bool const stepLimit) {
/*
The current step did not yet reach the ObservationCubic, do nothing now and
wait:
......@@ -42,7 +42,7 @@ inline ProcessReturn ObservationCubic<TTracking, TOutput>::doContinuous(
}
HEPEnergyType const energy = particle.getEnergy();
Point const pointOfIntersection = particle.getPosition();
Point const pointOfIntersection = step.getPosition(1);
DirectionVector const dirction = particle.getDirection();
// add our particles to the output file stream
......@@ -69,7 +69,7 @@ inline LengthType ObservationCubic<TTracking, TOutput>::getMaxStepLength(
CORSIKA_LOG_TRACE("getMaxStepLength, particle={}, pos={}, dir={}, cubic={}",
particle.asString(), particle.getPosition(),
particle.getDirection(), cubic_.asString());
particle.getDirection(), asString());
auto const intersection =
TTracking::intersect(particle, static_cast<Cubic const>(*this));
......
......@@ -100,10 +100,14 @@ inline Intersections Tracking::intersect(TParticle const &particle,
SpeedType vx = velocity.getX(cs);
SpeedType vy = velocity.getY(cs);
SpeedType vz = velocity.getZ(cs);
CORSIKA_LOG_TRACE("particle in cubic coordinate: position: ({:.3f}, {:.3f}, "
"{:.3f}) m, veolocity: ({:.3f}, {:.3f}, {:.3f}) m/ns",
x0 / 1_m, y0 / 1_m, z0 / 1_m, vx / (1_m / 1_ns),
vy / (1_m / 1_ns), vz / (1_m / 1_ns));
auto get_intersect_min_max = [](LengthType x0, SpeedType v0, LengthType dx) {
auto t1 = (x0 - dx) / v0;
auto t2 = (x0 + dx) / v0;
auto t1 = (dx - x0) / v0;
auto t2 = (-dx - x0) / v0;
if (t1 > t2)
return std::make_pair(t1, t2);
else
......@@ -115,8 +119,10 @@ inline Intersections Tracking::intersect(TParticle const &particle,
auto [tz_max, tz_min] = get_intersect_min_max(z0, vz, cubic.getZ());
TimeType t_exit = std::min(std::min(tx_max, ty_max), tz_max);
TimeType t_enter = std::min(std::min(tx_min, ty_min), tz_min);
TimeType t_enter = std::max(std::max(tx_min, ty_min), tz_min);
CORSIKA_LOG_DEBUG("t_enter: {} ns, t_exit: {} ns", t_enter / 1_ns,
t_exit / 1_ns);
if ((t_exit > t_enter)) {
if (t_enter < 0_s && t_exit > 0_s)
CORSIKA_LOG_DEBUG("numericallyInside={}", cubic.contains(position));
......
......@@ -12,12 +12,12 @@ namespace corsika {
// a CoordinateSystemPtr to specify the orintation of coordinate
Cubic(Point const& center, CoordinateSystemPtr cs, LengthType const x, LengthType const y, LengthType const z)
: center_(center)
, cs_(make_translation(cs, center.getCoordinates(center.getCoordinateSystem())))
, cs_(make_translation(cs, center.getCoordinates(cs)))
, x_(x), y_(y), z_(z) {}
Cubic(Point const& center, CoordinateSystemPtr cs, LengthType const side)
: center_(center)
, cs_(make_translation(cs, center.getCoordinates(center.getCoordinateSystem())))
, cs_(make_translation(cs, center.getCoordinates(cs)))
, x_(side/2), y_(side/2), z_(side/2) {}
......
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