I have code for circle-circle intersection. But I need to expand it to 3-D. How do I calculate:
- Radius and center of the intersection circle of two spheres
- Points of the intersection of a sphere and a circle?
Given two spheres (sc0,sr0)
and (sc1,sr1)
, I need to calculate a circle of intersection whose center is ci
and whose radius is ri
.
Moreover, given a sphere (sc0,sr0)
and a circle (cc0, cr0)
, I need to calulate the two intersection points (pi0, pi1)
I have checked this link and this link, but I could not understand the logic behind them and how to code them.
I tried ProGAL library for sphere-sphere-sphere intersection, but the resulting coordinates are rounded. I need precise results.
Answer
Sphere-Sphere Intersection
Let's start with the more obvious one - sphere-sphere. It's almost identical to the circle-circle case in 2D. We can project down on any plane containing the line between the spheres' centers to get an identical 2D picture:
Here the first sphere has center c_1
and radius r_1
, the second c_2
and r_2
, and their intersection has center c_i
and radius r_i
. Let d = ||c_2 - c_1||
, the distance between the spheres.
Edge Cases:
As was pointed out previously, if r_1 + r_2 < d
then there is no intersection at all. The spheres have a separation between them.
If r_1 + r_2 == d
, then the intersection is a single point, located a distance of r_1
on the line from c_1
to c_2
, or: c_i = c_1 + (c_2 - c_1) * r_1/d
This can also happen from the other side, if d + min(r_1, r_2) = max(r_1, r_2)
. In this case, the single point of intersection lies on the larger sphere, at c_i = c_larger + (c_smaller - c_larger) * r_larger/d
Finally, if one circle is entirely encompassed within the other, such that d + min(r_1, r_2) < max(r_1, r_2)
, then once again there is no intersection.
Typical intersections
Okay, now the meaty bit.
Let's define c_i = c_1 + h * (c_2 - c_1)
, for some value of h
. Note that h
can be negative or greater than one in the case where one circle encompasses the center of the other. Fortunately, the signs will "just work," so we don't need to add any special case logic for this condition.
We can express r_1
and r_2
in terms of r_i
, h
, and d
using Pythagorean Theorem, then solve the system of equations to find h:
$$\require{cancel} \begin{array}{r|l} (hd)^2 + r_i^2 = r_1^2 & \left((1-h)d\right)^2 + r_i^2 = r_2^2\\ h^2d^2 - r_1^2 = -r_i^2 & \left(1 - 2h + h^2 \right)d^2 - r_2^2 = -r_i^2 \end{array}\\ \begin{align} \cancel{h^2d^2} - r_1^2 &= d^2 - 2hd^2 + \cancel{h^2d^2} - r_2^2\\ 2hd^2 &= d^2 + r_1^2 - r_2^2\\ h &= \frac 1 2 + \frac {r_1^2 - r_2^2} {2 d^2} \end{align}$$
h = 1/2 + (r_1 * r_1 - r_2 * r_2)/(2 * d*d)
We can sub this into our formula for c_i above to find the center of the circle of intersections. Then, reversing one of our earlier Pythagorean relations:
r_i = sqrt(r_1*r_1 - h*h*d*d)
So, now we have the center and radius of our intersection. Now we can revolve this around the separating axis to get our full circle of solutions. The circle lies in a plane perpendicular to the separating axis, so we can take n_i = (c_2 - c_1)/d
as the normal of this plane.
After choosing a tangent and bitangent t_i
and b_i
perpendicular to this normal and each other, you can write any point on this circle as:
p_i(theta) = c_i + r_i * (t_i * cos(theta) + b_i sin(theta));
Because of the Hairy Ball Theorem, there's no one universal way to choose the tangent/bitangent to use. My recommendation would be to pick one of the coordinate axes not parallel to n_i
, and set t_i = normalize(cross(axis, n_i))
, and b_i = cross(t_i, n_i)
or somesuch.
Sphere-Circle Intersection
A circle, being flat/2-dimensional, can only intersect with things within its own plane. So, the first thing is: find out where that plane intersects our sphere.
Given a circle with center & radius c_c
and r_c
, respectively, in the plane with (unit) normal n
, and a sphere with center & radius c_s
and r_s
...
The plane of the circle cuts the sphere d = dot(n, c_c - c_s)
units from the sphere's center. Again, this can be negative if the normal is pointing toward the sphere's center, but the signs will work out.
If abs(d) > r_s
then there is no intersection. Our plane passes above/below the sphere entirely.
If there is an intersection, it will be between our original circle and a new one formed where this plane meets the sphere, with center c_p = c_s + d*n
.
If d == r_s
then this is the sole point of intersection with the plane; otherwise we have a circle whose radius r_p
we can find with Pythagorean Theorem:
r_p = sqrt(r_s*r_s - d*d)
And now we've reduced the problem to circle-circle, which we already know how to solve. Checking the distance/radii involved according to the edge cases above will let us take care of the non-intersecting and single point of intersection cases.
If the two circles do overlap by more than one point, we can use the formulae from the Sphere-Sphere check to find c_i
and r_i
as before.
Instead of a full circle of intersection points though, we'll have only the two points lying in the plane of the circle.
To express these points, we can define a tangent vector in the plane,
t = normalize(cross(c_p - c_c, n))
Then our intersection points are:
p_0 = c_i - t * r_i
p_1 = c_i + t * r_i
And we're done. :D
I hope that helps!
No comments:
Post a Comment