Solving stereographic projection for central latitude $\phi_1$ and central longitude $\lambda_0$

Question

For a given longitude and latitude and their projected point I want to know the central longitude and latitude of the stereographic projection. I used the formulars from wolfram alpha: http://mathworld.wolfram.com/StereographicProjection.html

$$\phi = sin^{-1}\left( \cos{c} \sin{\phi_1}+ \frac{y \sin{c} \cos{\phi_1} }{ \rho}\right)$$ $$\lambda = \lambda_0 + tan^{-1}\left(\frac{x \sin{c}}{ \rho \cos{\phi_1} \cos{c} - y \sin{\phi_1} \sin{c}}\right)$$

For now I found a recursive formular and as an approximation my starting point is $\phi_{pred} = \phi$. This approximation gives quite wrong results for certain inputs.

Recursive: $$\phi_1 = \sin^{-1}\left(-\frac{y \cos{\phi_{pred}} \sin{c} - \rho \sin{\phi}}{\rho \cos{c}}\right)$$ $$\lambda_0 = \lambda - tan^{-1}\left(\frac{x \sin{c}}{ \rho \cos{\phi_1} \cos{c} - y \sin{\phi_1} \sin{c}}\right)$$

I'm using SageMath.

Thanks for any hint on how to proceed.

Answer 1

Rather than an iterative method, we could try a fixed algorithm. To find one, let's start with a (probably) wrong answer and try to find exactly what we need to do to make it correct.

First, suppose a point at latitude $\phi,$ longitude $\lambda$ is projected onto coordinates $(x,y)$ by a stereographic projection centered at $\phi_1,\lambda_0,$ where the $y$ axis is in the plane of the $\lambda$ meridian and pointing "north." (I believe this is the assumption made on the page you linked; without such an assumption we would also need a parameter representing the orientation of the axes in the projection plane.) The coordinates $\phi,$ $\lambda,$ $x,$ and $y$ are known and given to you but you are not given $\phi_1$ or $\lambda_0.$

Let's start by assigning Cartesian coordinates $(X,Y,Z)$ to the space in which the sphere is embedded, with the origin at the center of the sphere, the $Z$ axis pointing upward through the north pole, and the $X$ axis pointing through the point at latitude zero and longitude zero. (The usual convention is to use lower-case letters for these axes, but you're already using $x$ and $y$ for a different set of coordinates.) So if the sphere has radius $R,$ the coordinates of the point at latitude $\phi$ and longitude $\lambda$ are \begin{align} X_o &= R \cos\phi \cos\lambda, \\ Y_o &= R \cos\phi \sin\lambda, \\ Z_o &= R \sin\phi. \end{align}

Suppose a stereographic projection centered at $(R,0,0)$ (latitude zero and longitude zero) projected a point $(X_c,Y_c,Z_c)$ onto the projected coordinates $(x,y).$ Then \begin{align} X_c &= (4R^2 - x^2 - y^2)RS, \\ Y_c &= 4R^2Sx, \\ Z_c &= 4R^2Sy \end{align} where $S = 1/(4R^2 + x^2 + y^2).$ (Compare this with Understanding the formula for stereographic projection of a point. Take into account that the axes are labeled differently there and that the center of that projection is on the $z$ axis rather than the $X$ axis.)

Now let's rotate the point $(X_c,Y_c,Z_c)$ around the $Y$ axis by an angle $\alpha$ so that its image is a point on the plane $Z = Z_o,$ that is, so that its image has latitude $\phi,$ the same as the point you were given. If the original data you were given were from an actual stereographic projection as described, such rotations will exist and you can choose one with angle $\phi_1$ so that $-\frac\pi2 \leq \phi_1 \leq \frac\pi2.$ That is, the rotation takes the point $(R,0,0)$ to the latitude $\phi_1$ of the center of the desired stereographic projection.

Let the image of $(X_c,Y_c,Z_c)$ after rotation be $(X_m,Y_m,Z_o).$ Rotate this point around the $Z$ axis so that its image is $(X_o,Y_o,Z_o).$ That is, this rotation takes the center of the projection to the correct longitude to project the given point at latitude $\phi,$ longitude $\lambda$ onto the coordinates $(x,y).$ Let $\lambda_0$ be the angle of this rotation.

You then have found $\phi_1$ and $\lambda_0,$ the coordinates of the center of the desired stereographic projection. (Note that you need to define the positive and negative directions of rotation around each axis so that the signs of $\phi_1$ and $\lambda_0$ will be correct.)

Note that the angle $\phi_1$ is not always uniquely determined by this method. There may be two choices. In those cases either choice will produce a stereographic projection with the desired results. You cannot know which one is the "real" projection unless you have other information to help you decide. The coordinates of a second projected point might be sufficient.