Approximating an Ellipse given 4 points.

Question

I am facing a problem in which I need to create masks for a certain region. This region is not a perfect ellipse, but for all intents and purpose, the ellipse needs to encapsulate these 4 coordinates.

I am finding the center of these 4 coordinates by using two triangles created from the 4 points. The resulting centre point I use to create a set of coordinates for the mask.

The following is matlab code as to how I do this:

triangleOneCenterX = (firstX + fourthX + thirdX)/3;
triangleOneCenterY = (firstY + fourthY + thirdY)/3;

triangleTwoCenterX = (fourthX + secondX + thirdX)/3;
triangleTwoCenterY = (fourthY + secondY + thirdY)/3;

finalCenterX = (triangleOneCenterX + triangleTwoCenterX)/2;
finalCenterY = (triangleOneCenterY + triangleTwoCenterY)/2;

xRadius = sqrt((fourthX-thirdX)^2 + (fourthY-thirdY)^2)/2;
yRadius = sqrt((firstX-secondX)^2 + (firstY-secondY)^2)/2;

theta = 0 : 0.01 : 2*pi;
x = xRadius * cos(theta) + finalCenterX;
y = yRadius * sin(theta) + finalCenterY;

The resulting mask, and image is as follows: (Disregard the inner points)

As you can see two of the points do not lie on the ellipse, This issue happens either on the "pseudo-y" or "x" axis. Two points would line up perfectly, but not the rest.

I am struggling to fix this mistake, my theory as to why it is wrong is how I am computing the centre point. Either, the triangle method is wrong, as I computing the centre of the 4 points, not necessarily the ellipse that satisfies the 4 points.

Any insight would be extremely helpful! I thank you regardless for taking time out of your day and reading this post this far :)

Kind regards,

Answer 1

You have $4$ points $P_k = (x_k, y_k ) , k = 1, 2,3,4 $. And you want to generate the equation of an ellipse that passes through these four points. Such an ellipse would have to satisfy the general conic equation:

$ A x^2 + B xy + C y^2 + D x + E y + F = 0 $

For an ellipse $A$ cannot be zero, thus for each of the four points, the equation becomes

$ x_k^2 + B x_k y_k + C y_k^2 + D x_k + E y_k + F = 0, k = 1, 2, 3, 4 $

And this is a $4 \times 5$ linear system in the five unknowns $B, C, D, E, F$,

$ Q a = b $

where

$ Q = \begin{bmatrix} x_1 y_1 && y_1^2 && x_1 && y_1 && 1 \\ x_2 y_2 && y_2^2 && x_2 && y_2 && 1 \\ x_3 y_3 && y_3^2 && x_3 && y_3 && 1 \\ x_4 y_4 && y_4^2 && x_4 && y_4 && 1 \end{bmatrix} $

$ a = [ B, C, D, E, F]^T $

$ b = [ -x_1^2 , - x_2^2 , - x_3^2 , -x_4^2 ]^T $

Assuming $Q$ has full row rank (i.e. assuming it has $4$ linearly independent rows), then, we can write

$ a = Q^T u + w$

where $w$ is in the null space of $Q$, and $u \in \mathbb{R}^4$. It follows that

$ Q (Q^T u + w) = b $

which reduces to

$ {Q Q}^T u = b $

Therefore,

$ u = ({Q Q}^T )^{-1} b $

And

$ a = Q^T u + w = V + t W \hspace{35pt} (*)$

where $V = Q^T ({Q Q}^T )^{-1} b$, and $W$ is a fixed vector in the one-dimensional null space of $Q$.

To ensure that we get an ellipse, we have to impose the condition:

$ C - \dfrac{1}{4} B^2 > 0 $

So that

$ ( V_2 + t W_2 ) - \dfrac{1}{4} ( V_1 + t W_1 )^2 \gt 0 $

And this re-arranges to

$ W_1^2 t^2 + (2 V_1 W_1 - 4 W_2 ) t + (V_1^2 - 4 V_2) \lt 0 $

And the solution of this inequality is the interval $[t_1, t_2]$ where $t_1, t_2 $ are the two roots of the quadratic on the left, $t_1 \lt t_2 $.

Finally, once you select $ t \in [t_1, t_2] $ and obtain $a$, then you have the ellipse equation in algebraic form

$ x^2 + B x y + C y^2 + D x + E y + F = 0 $

To draw this ellipse, you need to find the vector equation of the ellipse, i.e. you want to express the points $(x,y)$ in terms the center and semi-minor and semi-major axes vectors.

Define

$ r = [x, y]^T $

$ Q = \begin{bmatrix} 1 && \dfrac{B}{2} \\ \dfrac{B}{2} && C \end{bmatrix} $

This is a new $Q$ not to be confused with the $Q$ described in the previous analysis.

Also, define the vector

$G = \begin{bmatrix} D \\ E \end{bmatrix} $

Now the algebraic equation of the ellipse is

$ r^T Q r + r^T G + F = 0 $

Start by finding the center of the ellipse. It is given by

$ r_0 = - \dfrac{1}{2} Q^{-1} G $

Substituting this, the above algebraic equation of the ellipse becomes

$ ( r - r_0)^T Q (r - r_0) = r_0^T Q r_0 - F $

Dividing through by the right hand side, yields

$ (r - r_0)^T P (r - r_0) = 1 $

where $ P = \dfrac{1}{ r_0^T Q r_0 - F} Q $

Now diagonalize matrix $P$ into $ P = R D R^T $, with $R$ orthogonal and $D$ diagonal.

The formulas are

$ R = \begin{bmatrix} \cos(\theta) && - \sin(\theta) \\ \sin(\theta) && \cos(\theta) \end{bmatrix} $

where

$ \theta = \dfrac{1}{2} \tan^{-1} \left( \dfrac{ 2 P_{12}}{P_{11} - P_{22}} \right) $

The diagonal entries of $D$ are

$D_{11} = \dfrac{1}{2} ( P_{11} + P_{22} ) + \dfrac{1}{2} ( P_{11} - P_{22} ) \cos(2 \theta) + P_{12} \sin(2 \theta) $

$D_{22} = \dfrac{1}{2} ( P_{11} + P_{22} ) - \dfrac{1}{2} ( P_{11} - P_{22} ) \cos(2 \theta) - P_{12} \sin(2 \theta) $

finally, define the vector $z$ by $ r = r_0 + R z $ , then it follows

$ z^T D z = 1 $

i.e.

$ D_{11} z_1^2 + D_{22} z_2^2 = 1$

Its solution, leads to the parameterization of $r$. The solution is

$ z = \begin{bmatrix} \dfrac{1}{\sqrt{D_{11}}} \cos(\phi) \\ \dfrac{1}{\sqrt{D_{22}}} \sin(\phi) \end{bmatrix} $

Then from above

$ r = r_0 + R z $

Inspecting this equation, one realizes that the semi-axes of the ellipse are the columns of $R$ but scaled by $\dfrac{1}{\sqrt{D_{11}}}$ and $\dfrac{1}{\sqrt{D_{22}}}$.

You may choose to find the minimum area ellipse among all possible ellipses, and this is easily achievable because the ellipses are a function of a single variable $t$. Using a function minimizer, the critical value of $t$ that minimizes the ellipse area can be found.

As an example, let the four points be $(1,0), (9,3), (7,10), (3, 6) $

The image below shows the minimum area ellipse that passes through these four points.