Determine if a point is in a rotated hyperellipsoid

162 Views Asked by At

I asked this question a few days ago. I am implementing an algorithm to test whether a set of points is within an hyper-ellipsoid (HE). That post provided the correct answer for a non-rotated HE centered at the origin. However, I need a solution for an HE that is not at the origin and may be rotated in any dimension.

So, either I rotate all my data points (I assume, relative to the center of the HE?) or I account for the position and rotation of the HE. In either case, I am unsure of how to proceed.

From the answer to my previous question I get that there is a pattern that generalizes to n-dimensions:

$$ \left( {x \over a} \right) ^2 + \left({y \over b}\right) ^2 + \left({z \over c}\right) ^2 + \ldots \left( {n \over d_n }\right) ^2 \le 1 $$

And, after a bit more research, I see for the 2-dimensional case that rotations and translations can be accounted for by the following:

$$ \cfrac{ \left( \cos (\alpha)(x_p - x_0) + \sin (\alpha)(y_p - y_0) \right)^2 }{a^2} + \cfrac{ \left( \sin (\alpha)(x_p - x_0) - \cos (\alpha)(y_p - y_0) \right)^2 }{b^2} \le 1 $$

where:

  • $x_p$, $y_p$ are the data point coordinates
  • $x_0$, $y_0$ are the coordinates of the HE center
  • $\alpha$ is the rotation angle in x, y of the HE
  • $a$, $b$ are the half-length axes of the HE

So, how does this generalize to the n-dimensional case? Is there a similar predictable pattern that emerges in the trigonometry?

1

There are 1 best solutions below

5
On BEST ANSWER

First, rewrite your original test slightly as $${x_1^2\over a_1^2}+\cdots+{x_n^2\over a_n^2}-1\le0.$$ Any orientation-preserving affine transformation will preserve the direction of the inequality, so you can either transform all of your points into this “standard” coordinate system and test them against this relatively simple inequality or transform the ellipsoid into the coordinate system of the data and test the points against a potentially more complex expression. I suspect that doing the former will be more efficient.

As for how to transform the equation of the ellipsoid, it’s easiest for me to describe in matrix form. Any quadric can be expressed by an equation of the form $\mathbf x^TQ\mathbf x=0$, where $Q$ is a symmetric matrix and $\mathbf x^T=(x_1,\dots,x_n,1)$. Under the point transformation $\mathbf x'=M\mathbf x$, substituting for $\mathbf x$ gives $(M^{-1}\mathbf x')^TQ(M^{-1}\mathbf x')={\mathbf x'}^T(M^{-T}QM^{-1})\mathbf x'$, so the matrix of the transformed quadric is $Q'=M^{-T}QM^{-1}$. For your application, the standardized ellipsoid that you’d be starting with has the matrix $Q = \operatorname{diag}(1/a_1^2,\dots,1/a_n^2,-1)$.