# Optimizing Ray-Sphere Intersections

This topic is 2135 days old which is more than the 365 day threshold we allow for new replies. Please post a new topic.

## Recommended Posts

I've made a few ray tracers in the past, but they where normally focussed on small aspects, taking shortcuts for many things, however recently I needed to make one for a uni project and decided to go all out. So after setting up my math libs, getting all the ground work into place (which most importantly includes a 4DOF FPS camera) it came time to trace a few spheres to test both the ray casting and get ready to implement a BVH/BSP scheme.

First thing I tried was this snippet of code (take from this wonderful set of articles):

Vector dst = vCenter - r.GetOrigin();
float b = dst.dot(r.GetDirection());
float c = b * b - dst.dot(dst) + fRadiusSquared;

if(c < 0)
return INTERSECT_NONE;

float d = std::sqrtf(c), r1 = b - d, r2 = b + d;
if(r2 < 0)
return INTERSECT_NONE;

if(r1 < 0)
{
fDistance = r2;
return INTERSECT_INSIDE;
}

fDistance = r1;
return INTERSECT_OUTSIDE;


it works, but far off objects come out very "fuzzy" and jagged, clearing up as you get closer:

After trying quite a few variants (such as this) which all suffered the same problem, I suspected a problem in my maths code; however, using the following code (adapted from here, based on an unoptimized quadratic formula, massively slows down rendering):

	//Compute A, B and C coefficients
Vector o = r.GetOrigin() - vCenter;
float c = o.length_squared_fast() - fRadiusSquared;
if(c < 0)
return INTERSECT_INSIDE;

float a = r.GetDirection().dot(r.GetDirection());
float b = 2 * r.GetDirection().dot(o);

//Find discriminant
float disc = b * b - 4 * a * c;

// if discriminant is negative there are no real roots, so return
// false as ray misses sphere
if (disc < 0)
return INTERSECT_NONE;

// compute q as described above
float disc_sqrt = sqrtf(disc);
float q = (b < 0) ? (-b - disc_sqrt) / 2.0f : (-b + disc_sqrt) / 2.0f;

// compute t0 and t1
float t0 = q / a;
float t1 = c / q;

// make sure t0 is smaller than t1
if (t0 > t1)
std::swap(t0,t1);

// if t1 is less than zero, the object is in the ray's negative direction
// and consequently the ray misses the sphere
if (t1 < 0)
return INTERSECT_NONE;

// if t0 is less than zero, the intersection point is at t1
if (t0 < 0)
{
fDistance = t1;
return INTERSECT_OUTSIDE;
}
// else the intersection point is at t0
else
{
fDistance = t0;
return INTERSECT_OUTSIDE;
}


it comes out correctly:

and for the life of me a can't figure out why the optimized variants seem to degrade so badly...

Anybody got any hints/advice as to whats going wrong here? is this just a side-effect from approximating it geometrically?

##### Share on other sites

The first version doesn't account for so-called cancellation: One root will be inaccurate if sqrt(D) is close to b (subtracting almost equal numbers will "cancel" lots of significant bits). This is even worse when working with single precision floats, of course. This wiki entry explains it quite well.

You trade accuracy for speed here.

Edit: Sorry, you were right, it's (also) the different approach, not just a numerical issue. Nonetheless: Speed vs accuracy still holds.

Edited by unbird

##### Share on other sites

I may be somewhat debilitated at the moment, but seeing square root calls seems rather ... in need of optimization?

I could be wrong, but why not perform a 'closest point on a line' test between each ray and sphere, checking to see if the distance squared is less/equal to the sphere's radius squared ?

• 18
• 18
• 11
• 21
• 16