Line / cube intersection=

1 
pub struct Vec3 {
  x: f32;
  y: f32;
  z: f32;
}

pub struct Ray {
  starting_point: Vec3;
  direction: Vec3
}

Let i, j, k be integers with 0 <= i, j, k <= 7 (forming a 8x8x8 cube).

Given a ray R, we want to find all (i, j, k) such that R intersects the cube defined by [i, i+1], [j, j+1], [k, k + 1].

====

There seems to be all types of edge cases: ray starting inside/outside the cube; ray parallel to an edge / plane, etc ...

Is there a simple + numerically stable way to calculate the indices of all cubes we intersect in the 8x8x8 cube ? There seems to be lots of room for division by 0 when the ray is parallel or close to parallel to one of the edges / planes.

Thanks!

2

I would look at a classical line-drawing algorithm to see if it can be adapted to voxel space. In the worst case, you can project the problem onto two orthogonal planes, and then use a more precise test for the voxels in the intersection of those results.

3

Easy peasy: "Intersection of a Line and a Box" https://www.geometrictools.com/Documentation/IntersectionLineBox.pdf

:slight_smile:

4

The go-to algorithm for this problem is “A Fast Voxel Traversal Algorithm for Ray Tracing”, by John Amanatides and Andrew Woo, 1987. Compared to a sequence of box intersection tests, it has the advantage of being much cheaper per step, and it cannot skip any cubes due to numerical error or the ray perfectly hitting a corner, because it always steps from one cube to an adjacent cube (so, for example, raytracing an opaque shell will never have “light leaks” through the corners).

The paper is rather terse about how to initialize the algorithm, and there are some edge cases in doing that correctly, but I've sorted them out in my Rust implementation, all-is-cubes/src/raycast.rs. It is thoroughly tested including by fuzz testing. It does have some extra features you may not need (reporting faces, reporting distances, fast-forward to where the ray enters the bounds of the volume of interest), but if you strip out the code for them then they should not have any performance impact.

(The trickiest part of the initialization is captured in the function scale_to_integer_step(). Also note that when the ray direction is zero on one axis, part of t_delta becomes infinite, but the rest of the algorithm handles that just fine.)

4 Likes