Pretty Ray-to-AABB casting algorithm

1 files changed, 37 insertions, 10 deletions

diff --git a/src/Math/Ray.hpp b/src/Math/Ray.hpp
index f83a2b5..c2c2274 100644
--- a/src/Math/Ray.hpp
+++ b/src/Math/Ray.hpp
@@ -11,21 +11,48 @@ struct Ray {
         Vec3 normal;
     };
 
-    // https://gdbooks.gitbooks.io/3dcollisions/content/Chapter3/raycast_aabb.html
-    RaycastResult cast(AABB box) const {
+    // References:
+    // [1]: https://gdbooks.gitbooks.io/3dcollisions/content/Chapter3/raycast_aabb.html
+    // [2]: https://tavianator.com/2015/ray_box_nan.html
+    RaycastResult cast(AABB box, Real max_distance = 0) const {
+        RaycastResult hit{ true };
+
+        // `t` is the time at which the ray intersects a plane of the AABB.
+        // A single vector defines 3 planes, one for each axis.
+        // Our AABB ends up having 6 planes, taken from the `from` and `to` points,
+        // hence the 6 `t`s.
         Vec3 t_from = (box.min - origin) / direction;
         Vec3 t_to = (box.max - origin) / direction;
 
-        Real t_min = t_from.zip(t_to, LAMBDA(std::min, 2)).reduce(LAMBDA(std::max, 2));
-        Real t_max = t_from.zip(t_to, LAMBDA(std::max, 2)).reduce(LAMBDA(std::min, 2));
+        // NOTE: We have to use std::fmin and std::fmax here, since
+        // std::min and std::max are not IEEE 754 compliant. [2]
+        Vec3 smaller_t_values = t_from.zip(t_to, LAMBDA(std::fmin, 2));
+        Vec3 bigger_t_values = t_from.zip(t_to, LAMBDA(std::fmax, 2));
+
+        Real biggest_min_t = smaller_t_values.reduce(LAMBDA(std::fmax, 2));
+        Real smallest_max_t = bigger_t_values.reduce(LAMBDA(std::fmin, 2));
+
+        // The normal of the collision is a normal of the plane which is intersected
+        // by the ray at the biggest smaller `t` value.
+        // Since a plane always has two normals, we know which one to pick by looking
+        // at the direction of the ray.
+        hit.normal = smaller_t_values.zip(direction, [=](Real t, Real d) {
+            return Math::floats_equal(t, biggest_min_t) ? Math::sign(d) : 0;
+        });
+
+        // If the smallest `t` is negative, the ray is pointing away from the AABB.
+        // If the biggest `t` is smaller than the smallest `t`, the ray is missing the AABB,
+        // since that means we do not reach all the planes defined by `from` before we reach
+        // a plane defined by `to`, which only happens if there is no intersection.
+        if (smallest_max_t < 0 || biggest_min_t > smallest_max_t) return {};
 
-        if (t_max < 0 || t_min > t_max) return {};
+        // Get the intersection point from the collision time `t`
+        // with the formula `p = o + dt`.
+        hit.point = origin + direction * biggest_min_t;
+        if (max_distance > 0 && (hit.point - origin).magnitude() > max_distance)
+            return {};
 
-        return {
-            .hit = true,
-            .point = origin + direction * t_min,
-            .normal = t_from.zip(t_to, [](auto a, auto b) { return a < b ? -1.0 : 1.0; })
-        };
+        return hit;
     }
 
     Vec3 origin;