I have an ellipse, possibly rotated and shifted from the origin, which is given by a parametrization similar to this one: $$ \begin{aligned} x &= x_0 + a\cos\theta\cos\alpha - b\sin\theta\sin\alpha \\ y &= y_0 + a\cos\theta\sin\alpha + b\sin\theta\cos\alpha \\ \end{aligned} $$
where $(x_0,y_0)$ is the center of the ellipse, $(a,b)$ are its semi-axes, $\alpha$ is the rotation angle and $\theta\in [0,2\pi]$.
The reference point is given by $(p_0,q_0)$. Now I want to determine the parameters of the lines $y = mx + h$ that are tangent to the ellipse and go through that reference point. To do so, I compute $h$ from the reference point, $h = q_0-mp_0$, and then plug the ellipse equations into the line equation: $$ \underbrace{y_0-q_0}_{\equiv q} + a\cos\theta\sin\alpha + b\sin\theta\cos\alpha = m\left(\underbrace{x_0-p_0}_{\equiv p} + a\cos\theta\cos\alpha - b\sin\theta\sin\alpha \right) $$ From this I can get an equation for the slope $m$ in dependence on $\theta$. I also know that the slope of a tangent to the ellipse is given by $\frac{y'(\theta)}{x'(\theta)}$, i.e.: $$ m = \frac{a\sin\theta\sin\alpha - b\cos\theta\cos\alpha}{a\sin\theta\cos\alpha + b\cos\theta\sin\alpha} $$
Now, if I equate the two above equations for $m$, I get an expression which can be solved for $\theta$. I would expect that there are two possible solutions to this equations, as there are two tangent lines. I used WolframAlpha to solve this equation (see here), but it gave me four different solutions (which differ in the combination of two signs); I'm including a screenshot of one of the solutions below since it's quite long ($x\equiv\theta, t\equiv\alpha$):
The four solutions emerge as $x = \pm\cos^{-1}((\pm\sqrt{\dots}$. There are no additional constraints mentioned, so something must have gone wrong, but I can't figure out what. I verified that two of those four solutions are valid with the help of some examples (testing various reference points p0,q0
):
import matplotlib.pyplot as plt
from matplotlib.patches import Ellipse
import numpy as np
from numpy import arccos, cos, pi, sign, sin, sqrt
a, b = 3, 1 # ellipse parameters
center = 1, 2 # center of ellipse
angle = 60 # rotation of ellipse (counter-clockwise) [deg]
p0, q0 = -2, 3 # the reference point (x-,y-coordinates)
p = center[0] - p0
q = center[1] - q0
t = angle/180*pi
thetas = [
s1*arccos((s2*sqrt(a**2*b**4*(2*p*cos(t)**3 + 2*p*sin(t)**2*cos(t) + 2*q*sin(t)**3 + 2*q*sin(t)*cos(t)**2)**2 + 4*a**2*(-b**2*sin(t)**4 - b**2*cos(t)**4 - 2*b**2*sin(t)**2*cos(t)**2 + p**2*sin(t)**2 - 2*p*q*sin(t)*cos(t) + q**2*cos(t)**2)*(a**2*p**2*sin(t)**2 - 2*a**2*p*q*sin(t)*cos(t) + a**2*q**2*cos(t)**2 + b**2*p**2*cos(t)**2 + 2*b**2*p*q*sin(t)*cos(t) + b**2*q**2*sin(t)**2)) - a*b**2*(2*p*cos(t)**3 + 2*p*sin(t)**2*cos(t) + 2*q*sin(t)**3 + 2*q*sin(t)*cos(t)**2))/(2*(a**2*p**2*sin(t)**2 - 2*a**2*p*q*sin(t)*cos(t) + a**2*q**2*cos(t)**2 + b**2*p**2*cos(t)**2 + 2*b**2*p*q*sin(t)*cos(t) + b**2*q**2*sin(t)**2)))
for s1 in [-1, 1] for s2 in [-1, 1]
]
fig, ax = plt.subplots()
ax.add_patch(Ellipse(center, 2*a, 2*b, angle, color='tab:blue', alpha=0.4))
ax.scatter([p0], [q0], color='tab:blue', s=50)
for theta in thetas:
ex = center[0] + a*cos(theta)*cos(t) - b*sin(theta)*sin(t)
ey = center[1] + a*cos(theta)*sin(t) + b*sin(theta)*cos(t)
ax.axline((p0,q0), (ex,ey), color='tab:orange', lw=1)
ax.scatter([ex], [ey], color='tab:orange', s=25)
plt.show()
So it seems that all that's left is to find a way how to pick the right combinations of signs for the two tangent points. However, I don't know what condition I should apply. I tried to discriminate based on the quadrant of $(p_0-x_0,q_0-y_0)$, but that didn't work. Also, since there are 6 different ways to pick 2 out of 4 possible sign combinations, it seems that the condition should include 6 different "sections/intervals".