Finding Intersection of an ellipse with another ellipse when both are rotated

Solution 1:

J. Richter-Gebert's Perspectives on Projective Geometry section 11.4 gives a general recipe for intersecting conics. It also agrees with what Wikipedia has to say on intersecting conics. A rough outline is the following:

1. Turn your conics into symmetric matrices $A$ and $B$ in $\mathbb R^{3\times3}$ such that a point with coordinates $(x,y)$ lies on them if $$(x,y,1)\cdot A\cdot\begin{pmatrix}x\\y\\1\end{pmatrix}=0 \quad\text{and}\quad (x,y,1)\cdot B\cdot\begin{pmatrix}x\\y\\1\end{pmatrix}=0$$ You can compute this form from your equations, with a bit of work. The comment below will help you with this; for parabolas see also this answer.

2. Find a linear combination $C=\lambda A+\mu B$ such that $(\lambda,\mu)\neq(0,0)$ and $\det(C)=0$. This involves solving a cubic equation, so there will in general be three solutions. In this solution counting I've identified multiples of the same matrix, since these represent the same conic. Unless $A$ itself is already degenerate, you can assume $\mu=1$ so that you have a cubic equation in $\lambda$ only. Some solutions might not be real, which accounts for some of the different cases you mentioned.

3. These three solutions correspond to three degenerate conics through the four points of intersection of $A$ and $B$. A degenerate conic in this sense is a product of two lines. These two lines might coincide, or be complex, which again corresponds to some of the case distinctions you mentioned.

4. To obtain a pair of lines from the matrix of a degenerate conic, you have to add an antisymmetric matrix in such a way that the resulting matrix has rank $1$. This means that the determinants of all $2\times2$ minors of the resulting matrix have to be zero. This antisymmetric matrix has three parameters, for which you can formulate nine equations.

If you don't feel like solving nine (non-linear) equations, you can also do the following: compute the adjugate matrix of $C$. This will be a symmetric matrix of rank $1$. Take any nonzero row or column. Rows or columns only differ by a scalar factor. These coordinates you just obtained correspond to the homogenous coordinates of the point where the two lines intersect, which I'll call $p$. Then a suitably chosen multiple of these coordinates will correspond to the entries of your antisymmetric matrix. Written in formulas:

$$ \operatorname{adj}(C) = pp^T \qquad p = \begin{pmatrix}p_1\\p_2\\p_3\end{pmatrix} \qquad P = \begin{pmatrix} 0 & p_3 & -p_2 \\ -p_3 & 0 & p_1 \\ p_2 & -p_1 & 0 \end{pmatrix} \\ C' = C+\alpha P \qquad \operatorname{rank}(C') = 1 $$

In this case, you only need to solve a single quadratic equation to compute $\alpha$. You can use an arbitrary $2\times2$ minor to obtain this equation for $\alpha$. You will obtain two solutions which differ only by their sign, and you can choose either one.

Using either approach (3 parameters or using the adjugate matrix), you end up with a rank $1$ matrix describing the same degenerate conic. Any non-zero row of that matrix will be one of its component lines, any non-zero column the other. So now you have split the degenerate conic into two constituent lines.

5. You can now either intersect these lines with the conic (a conic-line intersection is described in section 11.3) or use two different solutions (from step 2.) to obtain two different pairs of lines and then intersect these lines, one from each pair, to obtain the four points of intersection. The approach intersecting pairs of lines might require the use of complex lines even if their intersection points might turn out to be real, so you'd best do all the computations using complex numbers, and only discard complex final results at the very end. The easiest way to intersect these lines is by computing their cross product and dehomogenize the result.

Cross references:

• Intersection of conics using matrix representation has an example with some nicely tuned numbers.
• My answer to Find intersection of hyperbola and ellipse. follows the above steps exactly, but the numbers there are rounded since that instance of the problem is not designed to have rational numbers.
• Decomposition of a degenerate conic focuses on an example for the decomposition in step 4.
• Algorithm: Intersection of two conics starts with a quadratic polynomial representation that is closer to the matrix from my step 1 than the original question asked here. Answers and contents there mention the book by Richter-Gebert which is the foundation of my answer as well.

Solution 2:

Another approach ...

Your first ellipse can be described by the parametric equations: $$ x = H_1 + a_1\cos A \cos\theta - b_1 \sin A \sin\theta \\ y = K_1 + a_1\sin A \cos\theta + b_1 \cos A \sin\theta $$ If we let $t = \tan \tfrac12\theta$, then $$ \cos\theta = \frac{1-t^2}{1+t^2} \quad ; \quad \sin\theta = \frac{2t}{1+t^2} $$ Substititing into the equations above gives rational quadratic expressions for $x$ and $y$ that describe the first ellipse.

Substitute these expressions for $x$ and $y$ into the equation of the second ellipse, and simplify. After a lot of algebra, you will get a quartic equation in $t$ whose roots describe the intersection points.

There are formulas that give the solutions of a quartic equation, but I would recommend that you use numerical methods, instead. There are many available software packages that will do the job. Ask again here if you need help finding one.

Solution 3:

In addition to what @bubba said, you can ease the computation as follows: translate the first ellipse center to the origin, counter rotate by angle $A$, and apply a rescaling with factors $\frac1{a_1}$, $\frac1{b1}$. This will turn the first ellipse in the unit circle.

Then, you can rewrite the second equation in its implicit form

$$\frac{((a_1x-\delta x)\cos(B-A)+(b_1y-\delta y)\sin(B-A))^2}{a_2^2}+\frac{((a_1x-\delta x)\sin(B-A)-(b_1y-\delta y)\cos(B-A))^2}{b_2^2} -1=0\\ =Ax^2+Bxy+Cy^2+Dx+Ey+F.$$

and plug the Weierstrass expressions for the circle, $x=\dfrac{1-t^2}{1+t^2}$ and $y=\dfrac{2t}{1+t^2}$, giving

$$(A+F-D)t^4+2(E-B)t^3+2(F+2C-A)t^2+2(B+E)t+A+D+F=0.$$

For an easy numerical solution, you can also use the parametric representation of the ellipse $2$ and discretize it using line segments or, better, circular arcs. This can be done by sampling points in triples or quads along the ellipse and bounding the deviation between the true curve and the linear or curvilinear approximation.

Then you reduce the problem to that of line/arc or arc/arc intersection, which is quite tractable.