A two star fix, using only the measured altitude of two stars, usually produces two possible location for
an observer. Usually one of the two points can be disregarded based on the observer's knowledge of their
approximate location, with one point being near by, and the other being thousands of miles away. It's
also possible to perform a three star fix, which is just doing the two star fix twice, but with a third
star to disambiguate which of the two points is the correct one.
\(\lambda_i, \delta_i, \zeta_i\) are the geographic position longitude, declination, and zenith distance of
the two observations.
$$
\begin{align*}
a_i &= \cos \lambda_i \cos \delta_i \\
b_i &= \sin \lambda_i \cos \delta_i \\
c_i &= \sin \delta_i \\
p_i &= \cos \zeta_i \\
\\
A &= a_1 b_2 - a_2 b_1 \\
B &= b_2 c_1 - b_1 c_2 \\
C &= b_2 p_1 - b_1 p_2 \\
D &= a_1 c_2 - a_2 c_1 \\
E &= b_1 c_2 - b_2 c_1 \\
F &= c_2 p_1 - c_1 p_2 \\
\\
0 & = x^2\left(1+\frac{A^2}{B^2} + \frac{D^2}{E^2}\right) + x\left(\frac{-2AC}{B^2} - \frac{2FD}{E^2}\right) + \frac{C^2}{B^2} + \frac{F^2}{E^2} - 1 \\
\end{align*}
$$
The above equation is a quadradic equation, and is solved using the quadradic formula, which yields two
solutions: \(x_1, x_2\).
$$
\begin{align*}
a &=1+\frac{A^2}{B^2} + \frac{D^2}{E^2} \\
\\
b &=\frac{-2AC}{B^2} - \frac{2FD}{E^2} \\
\\
c &=\frac{C^2}{B^2} + \frac{F^2}{E^2} - 1 \\
\\
x_1 &=\frac{-b+\sqrt{b^2 - 4ac}}{2a} \\
\\
x_2 &=\frac{-b-\sqrt{b^2 - 4ac}}{2a} \\
\end{align*}
$$
Now use the values of \(x_i\) to compute the y and z values.
$$
\begin{align*}
y_i &=\frac{F-Dx_i}{E} \\
\\
z_i &=\frac{C-Ax_i}{B} \\
\end{align*}
$$
\(x_i, y_i, z_i\) are the rectangular coordinates of the matching points. Convert them to polar
coordinates to get latitude and longitude.