Sign in to follow this  
taby

Double-checking ray-AABB intersection code

Recommended Posts

I've got some ray-AABB intersection code that I was hoping someone could take a look at to see if they see any major errors. As far as I know, it works, but I thought I would run it past the community. The code is fair game for anyone who would like to use it as well.
void ray_aabb_in(const vector_3 &start, const vector_3 &normalized_dir, const vector_3 &aabb_min, const vector_3 &aabb_max, vector_3 &entry_point)
{
	// Axis aligned bounding box face normals.
	static const vector_3 x_plus_normal(1, 0, 0);
	static const vector_3 x_minus_normal(-1, 0, 0);
	static const vector_3 y_plus_normal(0, 1, 0);
	static const vector_3 y_minus_normal(0, -1, 0);
	static const vector_3 z_plus_normal(0, 0, 1);
	static const vector_3 z_minus_normal(0, 0, -1);

	// If within the box, the entry point is the ray origin.
	if(within(start.x, aabb_min.x, aabb_max.x) && within(start.y, aabb_min.y, aabb_max.y) && within(start.z, aabb_min.z, aabb_max.z))
	{
		entry_point = start;
		return;
	}

	// Find which plane normal points opposite to normalized_dir.

	// Test against +x plane.
	if(start.x > aabb_max.x && 0 > normalized_dir.dot(x_plus_normal))
	{
		// Extend the unit line segment to whatever length that meets the +x plane.
		vector_3 final_dir = normalized_dir;
		final_dir *= d(start.x, aabb_max.x) / d(start.x, start.x + normalized_dir.x);

		vector_3 potential_point = start;
		potential_point += final_dir;

		// Is the end of the extended line segment within the volume defined by the infinitely extended yz plane?
		if(within(potential_point.y, aabb_min.y, aabb_max.y) && within(potential_point.z, aabb_min.z, aabb_max.z))
		{
			entry_point = potential_point;
			return;
		}
	}

	// Test again -x plane.
	else if(start.x < aabb_min.x && 0 > normalized_dir.dot(x_minus_normal))
	{
		vector_3 final_dir = normalized_dir;
		final_dir *= d(start.x, aabb_min.x) / d(start.x, start.x + normalized_dir.x);

		vector_3 potential_point = start;
		potential_point += final_dir;

		if(within(potential_point.y, aabb_min.y, aabb_max.y) && within(potential_point.z, aabb_min.z, aabb_max.z))
		{
			entry_point = potential_point;
			return;
		}
	}

	// Test against +y plane.
	if(start.y > aabb_max.y && 0 > normalized_dir.dot(y_plus_normal))
	{
		vector_3 final_dir = normalized_dir;
		final_dir *= d(start.y, aabb_max.y) / d(start.y, start.y + normalized_dir.y);

		vector_3 potential_point = start;
		potential_point += final_dir;

		if(within(potential_point.x, aabb_min.x, aabb_max.x) && within(potential_point.z, aabb_min.z, aabb_max.z))
		{
			entry_point = potential_point;
			return;
		}
	}

	// Test against -y plane.
	else if(start.y < aabb_min.y && 0 > normalized_dir.dot(y_minus_normal))
	{
		vector_3 final_dir = normalized_dir;
		final_dir *= d(start.y, aabb_min.y) / d(start.y, start.y + normalized_dir.y);

		vector_3 potential_point = start;
		potential_point += final_dir;

		if(within(potential_point.x, aabb_min.x, aabb_max.x) && within(potential_point.z, aabb_min.z, aabb_max.z))
		{
			entry_point = potential_point;
			return;
		}
	}

	// Test against +z plane.
	if(start.z > aabb_max.z && 0 > normalized_dir.dot(z_plus_normal))
	{
		vector_3 final_dir = normalized_dir;
		final_dir *= d(start.z, aabb_max.z) / d(start.z, start.z + normalized_dir.z);

		vector_3 potential_point = start;
		potential_point += final_dir;

		if(within(potential_point.x, aabb_min.x, aabb_max.x) && within(potential_point.y, aabb_min.y, aabb_max.y))
		{
			entry_point = potential_point;
			return;
		}
	}

	// Test against -z plane.
	else if(start.z < aabb_min.z && 0 > normalized_dir.dot(z_minus_normal))
	{
		vector_3 final_dir = normalized_dir;
		final_dir *= d(start.z, aabb_min.z) / d(start.z, start.z + normalized_dir.z);

		vector_3 potential_point = start;
		potential_point += final_dir;

		if(within(potential_point.x, aabb_min.x, aabb_max.x) && within(potential_point.y, aabb_min.y, aabb_max.y))
		{
			entry_point = potential_point;
			return;
		}
	}

	// If here, no match found, and entry_point is invalidated.
	entry_point.x = numeric_limits<double>::infinity();
}

void ray_aabb_out(const vector_3 &start, const vector_3 &normalized_dir, const vector_3 &aabb_min, const vector_3 &aabb_max, vector_3 &exit_point)
{
	// Axis aligned bounding box face normals.
	static const vector_3 x_plus_normal(1, 0, 0);
	static const vector_3 x_minus_normal(-1, 0, 0);
	static const vector_3 y_plus_normal(0, 1, 0);
	static const vector_3 y_minus_normal(0, -1, 0);
	static const vector_3 z_plus_normal(0, 0, 1);
	static const vector_3 z_minus_normal(0, 0, -1);

	// Find which plane normal points in the same direction as normalized_dir.

	// Test against +x plane.
	if(start.x < aabb_max.x && 0 < normalized_dir.dot(x_plus_normal))
	{
		// Extend the unit line segment to whatever length that meets the +x plane.
		vector_3 final_dir = normalized_dir;
		final_dir *= d(start.x, aabb_max.x) / d(start.x, start.x + normalized_dir.x);

		vector_3 potential_point = start;
		potential_point += final_dir;

		// Is the end of the extended line segment within the volume defined by the infinitely extended yz plane?
		if(within(potential_point.y, aabb_min.y, aabb_max.y) && within(potential_point.z, aabb_min.z, aabb_max.z))
		{
			exit_point = potential_point;
			return;
		}
	}

	// Test against -x plane.
	else if(start.x > aabb_min.x && 0 < normalized_dir.dot(x_minus_normal))
	{
		vector_3 final_dir = normalized_dir;
		final_dir *= d(start.x, aabb_min.x) / d(start.x, start.x + normalized_dir.x);

		vector_3 potential_point = start;
		potential_point += final_dir;

		if(within(potential_point.y, aabb_min.y, aabb_max.y) && within(potential_point.z, aabb_min.z, aabb_max.z))
		{
			exit_point = potential_point;
			return;
		}
	}

	// Test against +y plane.
	if(start.y < aabb_max.y && 0 < normalized_dir.dot(y_plus_normal))
	{
		vector_3 final_dir = normalized_dir;
		final_dir *= d(start.y, aabb_max.y) / d(start.y, start.y + normalized_dir.y);

		vector_3 potential_point = start;
		potential_point += final_dir;

		if(within(potential_point.x, aabb_min.x, aabb_max.x) && within(potential_point.z, aabb_min.z, aabb_max.z))
		{
			exit_point = potential_point;
			return;
		}
	}

	// Test against -y plane.
	else if(start.y > aabb_min.y && 0 < normalized_dir.dot(y_minus_normal))
	{
		vector_3 final_dir = normalized_dir;
		final_dir *= d(start.y, aabb_min.y) / d(start.y, start.y + normalized_dir.y);

		vector_3 potential_point = start;
		potential_point += final_dir;

		if(within(potential_point.x, aabb_min.x, aabb_max.x) && within(potential_point.z, aabb_min.z, aabb_max.z))
		{
			exit_point = potential_point;
			return;
		}
	}

	// Test against +z plane.
	if(start.z < aabb_max.z && 0 < normalized_dir.dot(z_plus_normal))
	{
		vector_3 final_dir = normalized_dir;
		final_dir *= d(start.z, aabb_max.z) / d(start.z, start.z + normalized_dir.z);

		vector_3 potential_point = start;
		potential_point += final_dir;

		if(within(potential_point.x, aabb_min.x, aabb_max.x) && within(potential_point.y, aabb_min.y, aabb_max.y))
		{
			exit_point = potential_point;
			return;
		}
	}

	// Test against -z plane.
	else if(start.z > aabb_min.z && 0 < normalized_dir.dot(z_minus_normal))
	{
		vector_3 final_dir = normalized_dir;
		final_dir *= d(start.z, aabb_min.z) / d(start.z, start.z + normalized_dir.z);

		vector_3 potential_point = start;
		potential_point += final_dir;

		if(within(potential_point.x, aabb_min.x, aabb_max.x) && within(potential_point.y, aabb_min.y, aabb_max.y))
		{
			exit_point = potential_point;
			return;
		}
	}

	// If here, no match found, and exit_point is invalidated.
	exit_point.x = numeric_limits<double>::infinity();
}

bool ray_aabb_intersection(const vector_3 &start, const vector_3 &normalized_dir, const vector_3 &aabb_min, const vector_3 &aabb_max, vector_3 &entry_point, vector_3 &exit_point)
{
	// Find potentially valid entry point.
	ray_aabb_in(start, normalized_dir, aabb_min, aabb_max, entry_point);

	// Is entry point valid?
	if(numeric_limits<double>::infinity() == entry_point.x)
		return false;

	// If so, find exit point.
	ray_aabb_out(start, normalized_dir, aabb_min, aabb_max, exit_point);

	return true;
}


Share this post


Link to post
Share on other sites
Here's a visual of many rays (top) intersecting with a series of axis-aligned bounding boxes (bottom).


Share this post


Link to post
Share on other sites
I haven't waded through the entire code so maybe you need something else, but the following is about as simple as you can get:


/*
* Ray-box intersection using IEEE numerical properties to ensure that the
* test is both robust and efficient, as described in:
*
* Amy Williams, Steve Barrus, R. Keith Morley, and Peter Shirley
* "An Efficient and Robust Ray-Box Intersection Algorithm"
* Journal of graphics tools, 10(1):49-54, 2005
*
*/
static bool intersect_bbox(const Vec3f bbox[2],const Ray& ray, float t0, float t1)
{
float tmin = (bbox[ray.sign[0]][0] - ray.org[0]) * ray.inv_dir[0];
float tmax = (bbox[1-ray.sign[0]][0] - ray.org[0]) * ray.inv_dir[0];

float tymin = (bbox[ray.sign[1]][1] - ray.org[1]) * ray.inv_dir[1];
float tymax = (bbox[1-ray.sign[1]][1] - ray.org[1]) * ray.inv_dir[1];

if ( (tmin > tymax) || (tymin > tmax) )
return false;
if (tymin > tmin)
tmin = tymin;
if (tymax < tmax)
tmax = tymax;

float tzmin = (bbox[ray.sign[2]][2] - ray.org[2]) * ray.inv_dir[2];
float tzmax = (bbox[1-ray.sign[2]][2] - ray.org[2]) * ray.inv_dir[2];

if ( (tmin > tzmax) || (tzmin > tmax) )
return false;
if (tzmin > tmin)
tmin = tzmin;
if (tzmax < tmax)
tmax = tzmax;
return ( (tmin < t1) && (tmax > t0) );

}


You can get an intersection point by grabbing tmin (if return is true).

ray.inv_dir[i] is 1.0 / ray.dir[i] and ray.sign[i] = (ray.inv_dir[i] < 0);

Share this post


Link to post
Share on other sites
Hey thank you very much. I will try this code and see how it compares (hopefully it's identical!).

Share this post


Link to post
Share on other sites

Create an account or sign in to comment

You need to be a member in order to leave a comment

Create an account

Sign up for a new account in our community. It's easy!

Register a new account

Sign in

Already have an account? Sign in here.

Sign In Now

Sign in to follow this