Clever Ray Sphere Intersection

This notebook presents a much faster and simpler means of computing the intersection between a sphere and a ray. There are several 'key insights' about the relative geometry of a ray and sphere that must be understood before the mathematics of this approach makes geometric sense. To try to help visualize these connections this notebook uses explicit values for an example.

Ross Beveridge, October 9, 2018

var('r');
r = 3;
Cv = vector(SR, 3, (5,5,5));
Lv = vector(SR, 3, (3,0,0));
Uv = vector(SR, 3, (0,4,4)); Uv = Uv / Uv.norm();
Tv = vector(SR, 3, (Cv - Lv));
pretty_print("Sphere center: C = ", Cv);
pretty_print("Ray start:     L = ", Lv);
pretty_print("Ray direction: U = ", Uv);
pretty_print("Base to Center: T = ", Cv - Lv);

$$\newcommand{\Bold}[1]{\mathbf{#1}}\verb|Sphere|\phantom{\verb!x!}\verb|center:|\phantom{\verb!x!}\verb|C|\phantom{\verb!x!}\verb|=| \left(5,\,5,\,5\right)$$

$$\newcommand{\Bold}[1]{\mathbf{#1}}\verb|Ray|\phantom{\verb!x!}\verb|start:|\phantom{\verb!xxxxx!}\verb|L|\phantom{\verb!x!}\verb|=| \left(3,\,0,\,0\right)$$

$$\newcommand{\Bold}[1]{\mathbf{#1}}\verb|Ray|\phantom{\verb!x!}\verb|direction:|\phantom{\verb!x!}\verb|U|\phantom{\verb!x!}\verb|=| \left(0,\,\frac{1}{2} \, \sqrt{2},\,\frac{1}{2} \, \sqrt{2}\right)$$

$$\newcommand{\Bold}[1]{\mathbf{#1}}\verb|Base|\phantom{\verb!x!}\verb|to|\phantom{\verb!x!}\verb|Center:|\phantom{\verb!x!}\verb|T|\phantom{\verb!x!}\verb|=| \left(2,\,5,\,5\right)$$

In this develoopment the base of the ray plays a critical role. Although it may not at first be obvious, in order to project the center of the sphere onto the ray with position L playing the role of the origin we need only take the dot product of T onto U.

The term 'disc' is short hand for comparing the radius of the circle made by the sphere intersection the ray-sphere-center-plane and the diameter of the circle centered on the sphere center and reaching out to the point where the ray comese closest to the sphere center. While the above sentence's language is a tad complex, it is corret and worth disecting while looking at the picture in PowerPoint.

Do recognize that when 'disc' is greater than zero its square root is the key measure 'd' shown on the PowerPoint slide.

v    = Tv.dot_product(Uv); 
csq  = Tv.dot_product(Tv); 
pretty_print("v = ", v, ", c*c = ", csq);
disc = r^2 - (csq - v^2);
pretty_print("disc = ", LatexExpr("{\;\;}"), disc);

$$\newcommand{\Bold}[1]{\mathbf{#1}}\verb|v|\phantom{\verb!x!}\verb|=| 5 \, \sqrt{2} \verb|,|\phantom{\verb!x!}\verb|c*c|\phantom{\verb!x!}\verb|=| 54$$

$$\newcommand{\Bold}[1]{\mathbf{#1}}\verb|disc|\phantom{\verb!x!}\verb|=| {\;\;} 5$$

If the intermediate value disc is greater than zero the ray intersects the sphere.

bs     = 10.0;
d      = 0.0 if (disc < 0) else sqrt(disc);
Pt     = Lv + (v-d) * Uv;
Ptfar  = Lv + sqrt(1.5 * (bs^2)) * Uv;
Ptprj  = Lv + v * Uv;
contents = [
  arrow3d(Lv, Lv + Uv, width=3),
  point(Pt, color='green', size=16),
  sphere(Cv, size=r, color = Color('#F1C40F'), opacity=0.33),
  polygon3d([Lv, Cv, Ptprj], color=Color('#444444'), alpha=0.33),
  point(Cv, color='gray', size=12),
  point(Ptprj, color='red', size=17),
  line3d([Lv, Cv], thickness=5, color=Color('#A04000')),
  line3d([Cv, Ptprj], thickness=5, color=Color('#0E6655')),
  line3d([Lv, Ptfar], thickness=5, color=Color('#6C3483')),
  line3d([(0,0,0),(bs, 0,  0)],thickness=5,color='red'),
  line3d([(0,0,0),(0, bs,  0)],thickness=5,color='green'),
  line3d([(0,0,0),(0,  0, bs)],thickness=5,color='blue')
  ];
show(sum(contents), aspect_ratio=(1,1,1));