I've seen lots of spins on the ray-box intersection test, and a lot of them seem to boil down to this code: (omitting some details here)

  invRay.x = 1.0 / ray.x;

  tx1 = (box.min.x - rayOrigin.x) * invRay.x;
  tx2 = (box.max.x - rayOrigin.x) * invRay.x;

  tmin = min(tx1, tx2);
  tmax = max(tx1, tx2);

  // Repeat for Y and Z
  // Then check to ensure tmin < tmax

(I got mine specifically from this discussion and this source. Also note that if you use std::min and std::max, it's far from branchless)

However, this raises a significant problem.

What if the ray's origin in one or more dimensions is equal to the min or max points?

What I found is that I get this situation on an example ray and box:

  • Ray origin: (3.225, 0.243, 2.221)
  • Ray direction: (0, 0, -1)
  • Box min: (3.223, 0.241, 2.119)
  • Box max: (3.225, 0.243, 2.121)

Is that when tx1 is calculated, it correctly gets -inf. However, when tx2 is calculated, it gets (3.225 - 3.225) * inf which is coming out to -nan. This means that later, when we do:

float tmax = std::max(tx1, tx2);

We get max(-inf, -nan) which comes out to -inf and from that point on in the tests, tmax will always be -inf, making the overall ray-box test fail.

This is confirmed by some comments on the page where I got the code:

Something to consider here is that 0 * inf =nan which occurs when the ray starts exactly on the edge of a box


Having read both this and part 2, I'm unable to work out what changes I need to make in order to reliably detect intersections when a ray is exactly on the edge of an AABB.

Does anyone know a fast method for ray-box intersection with edge and corner hits?

I should note this is not being used with a BVH. I need an exact intersection test, I cannot increase the size of the box at all.

  • $\begingroup$ Nuts, I have this bug in some of my code... I think, in this case, you want 0 * inf to be 0, since it happens at the origin of the Ray. That should result in a ray lying on one face of the box to be a hit, the opposite face a miss. But you need to check for the zero. $\endgroup$
    – user2500
    Feb 1, 2020 at 15:27
  • 1
    $\begingroup$ I think the pragmatic fix for BVH use would be to inflate the boxes by 1 ulps. I suspect that would beat any fix I can think of as all of them involve several extra instructions. I will keep thinking about it though, it is an interesting problem. $\endgroup$
    – Olivier
    Feb 2, 2020 at 0:30
  • $\begingroup$ Making the boxes physically larger in my case would make the results totally worthless. I can check if the result is NAN but that doesn't guarantee a hit, only that the origin started on one of those dimensions.Also, making them larger would just mean it would be the same problem, just at a different point in space, no? Or do you mean, if I detect a nan, extend the box and then try again? $\endgroup$ Feb 3, 2020 at 2:50
  • $\begingroup$ @TylerShellberg I meant always make them larger. This will give a tiny % of false positives, which for a BVH is of no importance. If you require an exact intersection test, you might want to edit your question to specify that and also if you care about fixing the problem on both the min/max side of each axis or only on one side. $\endgroup$
    – Olivier
    Feb 4, 2020 at 1:13
  • $\begingroup$ @Oliver I am not using this for a BVH. $\endgroup$ Feb 4, 2020 at 18:25

2 Answers 2


Here is how I found this can be solved.

(The below also happens at the start, just fyi)

  float tmin = std::numeric_limits<float>::lowest();
  float tmax = std::numeric_limits<float>::max();

Add an additional check:

if(rayDir.x == 0.0f)
  if (rayOrigin.x() < box.min.x() || rayOrigin.x() > box.max.x())
    return false;
  tx1 = (box.min.x - rayOrigin.x) * invRay.x;
  tx2 = (box.max.x - rayOrigin.x) * invRay.x;

  tmin = std::max(tmin, std::min(tx1, tx2));
  tmax = std::min(tmax, std::max(tx1, tx2));

This way the tx1/tx2 values never get calculated if the ray's direction vector is zero in that direction, aka if it's parallel to that plane.

It's not particularly compact or pretty, but it does neatly solve the issue.

  • $\begingroup$ You can specialize the intersection function based on the zero components (with 7 different versions), and avoid the runtime cost of the checks as you traverse for a given ray. $\endgroup$
    – user2500
    Feb 3, 2020 at 21:00
  • $\begingroup$ @DanielMGessel Could you clarify what you mean by that? $\endgroup$ Feb 3, 2020 at 21:09
  • $\begingroup$ Yeah, guess that was a bit sloppy. If you’re traversing a BVH for a single ray (on a scalar processor, not a gpu) consider having a version of the traversal code for the general case where no component of the ray direction is zero (most rays), eliminating your added checks, and a version with the checks. Select traversal code based on ray direction, Might speed things up on a cpu. $\endgroup$
    – user2500
    Feb 4, 2020 at 0:25
  • $\begingroup$ Just grokked the pbr book’s solution. I’m slow today. Check it out: pbr-book.org/3ed-2018/Shapes/Basic_Shape_Interface.html $\endgroup$
    – user2500
    Feb 4, 2020 at 0:44
  • $\begingroup$ @DanielMGessel How is checking what type of ray, then choosing a function, rather than the function determining it any faster? Either way, the exact same checks need to happen on the ray. It doesn't seem like less work is being done. $\endgroup$ Feb 4, 2020 at 15:48

This is a C++ implementation. Slab is the space between a pair of parallel planes.Test if ray R(t) = p + t * d intersect against bounding box bb.

bool IntersectRayBB( const double* p, const double* d, const BoundingBox& bb, double& tmin, double* q )
    double tolerance = 1e-9;
    tmin = 0.0;
    double tmax = FLOAT32_MAX; // Set to max distance ray can travel
    for( int i = 0; i < 3; i++ )
        if( Abs( d[i] ) < tolerance )
            // Ray is parallel to slab. No hit if origin not within slab.
            if( p[i] < bb.min[i] || p[i] > bb.max[i] )
                return false;
            // Compute intersection t value of ray within near and far plane of slab
            double denom = 1.0 / d[i];
            double t1 = ( bb.min[i] - p[i] ) * denom;
            double t2 = ( bb.max[i] - p[i] ) * denom;
            // Make t1 be intersection with near plane, t2 with far plane
            if( t1 > t2 )
                std::Swap( t1, t2 );
            // Compute the intersection of slab intersection intervals
            if( t1 > tmin )
                tmin = t1;
            if( t2 > tmax )
                tmax = t2;
            // Exit with no collision as soon as slab intersection becomes empty
            if( tmin > tmax )
                return false;

    //Ray intersects all 3 slabs. Return point (q) and intersesction t value (tmin)
    q = p + d * tmin;
    return true;

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.