Voxel Raycast Algorithm

[Elveskevtar]

If someone could help me figure out this problem that would be much appreciated!

I have a 3D Voxel grid being represented by a 1D array (for performance reasons). Each Voxel has dimensions Voxel.size * Voxel.size * Voxel.size. I am trying to figure out this algorithm: given a Ray figure out which Voxels it passes through and store the results in something like a HashSet or something of that sort.

My current best solution that I can think of is to make Box Colliders for all of the voxels nearby and see which ones intersect the Ray but I feel like there should be a more efficient way.

Thank you!

[Hyblademin]

We could always iterate along the ray in small increments and check the current box. If you will only ever be checking small lengths, it might be possible that this would work fine.

If you want to check long vectors or do a lot of them in the same update, something more efficient might be needed. You can probably find some good algorithms that do this with Google, but I’m procrastinating at work, so I’m going to try to work myself through this anyway.

I’m imagining a system where we choose an axis, find all of the unit box intervals that the ray moves through on the chosen axis, and select all of the boxes that are adjacent to those intersections based on the values on the other two axes.

I’ll try a test run in 2D. Here’s a plane divided into squares, and a ray that starts and ends within it:

If I choose the x-axis to search along and find the interval edges of the cells that are crossed while moving in x, I will find the x values represented by the vertical grid lines that intersect the ray. We can do this by finding the first interval that is included in the ray by moving forward from the x start value, then adding the length of a unit square as many times as we can without passing the x end value.

We then get the y-values on the ray at each of these intersections, giving the following:

Then select the boxes that are adjacent to these intersections ( (xn, yn) and (xn-1, yn) ), discarding overlaps:

That works fine for this example, but if we selected the y-axis for our search instead, we would have selected these:

In this case we obviously missed some. To avoid this, we should always select the axis with the highest rate of change relative to the other axes. In other words, choose the axis that has the biggest difference in the start and end values. In our example, we can confirm visually that this is the x-axis.

This is the case because after each time the other axis crosses a side, the axis with the highest rate of change is guaranteed to cross a side at least once before the other axis does again. When an axis other than our search axis crosses more than one side before our search axis crosses at least one, we will miss a square.

An important note: If the vector crosses exactly over a corner (a point where a vertical and horizontal grid line meet), we need to select all four squares in the intersection. In the algorithm, this can be done when the non-search axis values are calculated-- if they happen to be divisible by the length of our unit squares, then the point is at a node, and the surrounding cells should be selected.

Another note: If the change in both axes is the same, it doesn’t matter which is chosen for the search. If no grid lines are crossed, the result is just the cell that the start/end point is in.

To generalize this to 3 dimensions, we still need to start with the axis with the most change. We check both of the other coordinates at each interval on the search axis, and again check if they are divisible by the unit cube length. If they are, we may need to select four (intersects a cube edge) or eight (intersects a corner) cubes. If not, we select the two cubes adjacent to the found side, just like in 2 dimensions ( (xn, yn, zn) and (xn-1, yn, zn), assuming x is the search axis).

The biggest difference with 3D is that we now need to search another axis. Otherwise, it’s still possible to miss some cubes:

Here is my best drawing of a 3D ray moving through a segmented space, and all of the unit cubes it intersects. The ray has the greatest change along the y-axis. If we only search along the y-axis, we will find all of the sides filled with blue, and we will select the cubes with a blue dot. However, we will miss the cube with a red dot, unless we also search along the second-steepest axis, the z-axis, where we’ll find the side filled with red and select the remaining cube.

I think this is all pretty reasonable, since we know where all the sides, edges, and corners of our 3D grid are and can make a list of the ones intersected by the ray quickly. After that it’s just a matter of putting the adjacent cubes into a list and discarding duplicates.

[Elveskevtar]

@Hyblademin

We could always iterate along the ray in small increments and check the current box. If you will only ever be checking small lengths, it might be possible that this would work fine.

If you want to check long vectors or do a lot of them in the same update, something more efficient might be needed. You can probably find some good algorithms that do this with Google, but I’m procrastinating at work, so I’m going to try to work myself through this anyway.

I’m imagining a system where we choose an axis, find all of the unit box intervals that the ray moves through on the chosen axis, and select all of the boxes that are adjacent to that intersection based on the values on the other two axes.

I’ll try a test run in 2D. Here’s a plane divided into squares, and a ray that starts and ends within it:

If I choose the x-axis to search along and find the interval edges of the cells that are crossed while moving in x, I will find the x values represented by the vertical grid lines that intersect the ray. We can do this by finding the first interval that is included in the ray by moving forward from the x start value, then adding the length of a unit square as many times as we can without passing the x end value.

We then get the y-values on the ray at each of these intersections, giving the following:

Then select the boxes that are adjacent to these intersections ( (xn, yn) and (xn-1, yn) ), discarding overlaps:

That works fine for this example, but if we selected the y-axis for our search instead, we would have selected these:

In this case we obviously missed some. To avoid this, we should always select the axis with the highest rate of change relative to the other axes. In other words, choose the axis that has the biggest difference in the start and end values. In our example, we can confirm visually that this is the x-axis.

This is the case because after each time the other axis crosses a square, the axis with the highest rate of change is guaranteed to cross a square at least once before the other axis does again. When an axis other than our search axis crosses more than one square before our search axis crosses at least one, we will miss a square.

An important note: If the vector crosses exactly over a corner (a point where a vertical and horizontal grid line meet), we need to select all four squares in the intersection. In the algorithm, this can be done when the non-search axis values are calculated-- if they happen to be divisible by the length of our unit squares, then the point is at a node, and the surrounding cells should be selected.

Another note: If the change in both axes is the same, it doesn’t matter which is chosen for the search. If no grid lines are crossed, the result is just the cell that the start/end point is in.

To generalize this to 3 dimensions, we still need to start with the axis with the most change. We check both of the other coordinates at each interval on the search axis, and again check if they are divisible by the unit cube length. If they are, we may need to select four (intersects a cube edge) or eight (intersects a corner) cubes. If not, we select the two cubes adjacent to the found side, just like in 2 dimensions ( (xn, yn, zn) and (xn-1, yn, zn), assuming x is the search axis).

The biggest difference with 3D is that we now need to search another axis otherwise it’s still possible to miss some cubes:

Here is my best drawing of a 3D ray moving through a segmented space, and all of the unit cubes it intersects. The ray has the greatest change along the y-axis. If we only search along the y-axis, we will find all of the sides filled with blue, and we will select the cubes with a blue dot. However, we will miss the cube with a red dot, unless we also search along the second-steepest axis, the z-axis, where we’ll find the side filled with red and select the remaining cube.

I think this is all pretty reasonable, since we know where all the sides, edges, and corners of our 3D grid are and can make a list of the ones intersected by the ray quickly. After that it’s just a matter of putting the adjacent cubes into a list and discarding duplicates.

This is exactly what I was looking for. I found Bresenham’s Algorithm online and made a 3D implementation of it but it is just an approximation algorithm that describes what you mention in the beginning about taking samples at every few points. I am using this for a line of sight algorithm between my Third Person Camera and my player to detect camera occlusion with the terrain so the vectors are never very long.

Since Bresenham is an approximation algorithm it obviously misses a few voxels which I can clearly see as I have made gizmos to display wireframes of the cube positions that I am calculating the camera to be passing through. Removing the duplicates should not be a problem as I am using a HashSet. I will try your version of this algorithm as soon as work slows down and I solve another problem with actually capturing the right voxel - something weird is happening with my modular math.

Regardless thank you very much! This helps a ton and I’ll let you know how it works and I’ll probably even post my code in this thread once it is finished!

[Hyblademin]

I noticed something I didn’t account for before: It’s possible for the first and last cells to be missed by the algorithm above. Here’s an example:

In this example, the x-axis has the greatest change, so it’s the best choice for the primary search axis. We find the two blue-filled cells, but miss the red one because the first crossing of a side happens in the y direction.

The easiest solution to this is to just check the start and end points of the ray and add their cells to the list by truncating the coordinates off to the nearest unit length.

Another thought:

I’m wondering if this can be improved by removing the need to search a second axis in 3 or more dimensions. Except in the case that there is only one intersected cube, each cube is always adjacent to another intersected cube. Note: This follows easily from the fact that the ray is a continuous line, and there is no was for it to skip a cube. Also, in our algorithm, crossing an edge or corner selects all adjacent cubes, so this case is covered.

In the 3D example I drew before, we missed one cube after the first search:

While adding cubes during the first search, we could check that each new cube found is adjacent to another one. If we do, we’ll raise an alarm after we miss that cube. We could then solve this by doing either an incremental search, or a small-scale version of a secondary axis search between the last two cubes found.

Using this method, I think there’s still potential for missing cubes under certain specific circumstances, but I’m not a mathematician, so I don’t know how to check it for sure. But, maybe it will lead somewhere.