Querying the closest point on a transformed sphere

Question

Given the 4x4 transform matrix of a sphere and a point in space, I want to find the closest point on the sphere's surface.

Normally I would trace a line between the point and the sphere's center, and use the sphere's radius to get my solution, but here i am dealing with a non-uniformly scaled sphere. Here's a quick example in Python:

import numpy as np
from numpy.core.umath_tests import inner1d

# 4x4 transform matrix of a sphere with the following components:
# Scale XYZ = 1,5,1 (Scaled only in Y axis to keep this example simple)
# Rotation XYZ = 0,0,45 (Simple tilt for this example)
# Position XYZ = -1,3,0 (Position along XY plane, again for simplicity)
M = np.array([[ 0.70710678,  0.70710678,  0.        ,  0.        ],
              [-3.53553391,  3.53553391,  0.        ,  0.        ],
              [ 0.        ,  0.        ,  1.        ,  0.        ],
              [-1.        ,  3.        ,  0.        ,  1.        ]])

# Query point p0
p0 = np.array([-2,6,0])

# Transform the point into a unit sphere
I = np.linalg.inv(M)
p1 = np.array(p)-M[3,:3]
p1 = np.dot(p1,I)

# Normalize the point so it is on the surface of the unit sphere
mag = np.sqrt(inner1d(p1,p1)) # magnitude
p1 /= mag

# Transform back into 3D space
p1 = np.dot(p1,M[:3,:3]) + M[3,:3] #result [-1.65653216, 4.96959649, 0.]

This solution is fast and works well when the query point is already close to the sphere, but not so much when it is distant. See in the image above: point p2 which would be the desired result.

see my quest for [ray and ellipsoid intersection accuracy improvement](http://stackoverflow.com/q/25470493/2521214) — Spektre, May 03 '17 at 07:21

score 4 · Answer 1 · answered May 05 '17 at 09:33

I don't know if the iterative procedure below will work, but it comes to mind intuitively. Maybe is it worth giving it a try.

From the given point, draw a line to the center of the elliposid and get the intersection point with the surface.
Construct the normal plane at the intersection point and project the given point onto it.
Join the projection to the center, find the intersection with the surface and repeat from 2., until convergence.

I can't tell if convergence to the closest point can be guaranteed. All I can say is that when there is convergence, you have found an orthogonal projection of the given point onto the surface.

The picture shows the first iterates (intersections in green, projections in orange).

This is probably equivalent to Newton's iterations.

score 3 · Answer 2 · answered May 02 '17 at 23:14

You want to look at David Eberley's "Distance from a Point to an Ellipse, an Ellipsoid, or a Hyperellipsoid" (PDF download.) Ultimately you are finding roots of a quartic polynomial for a 2D ellipse, and roots of a 6-th degree polynomial for a 3D ellipsoid, so this is by no means a simple problem.

Given this complexity, you might want to instead seek an approximate result, e.g., by meshing the ellipsoid and seeking the closest mesh vertex.

Querying the closest point on a transformed sphere

2 Answers2