Aunque las geodésicas se parecen un poco a las ondas sinusoidales en algunas proyecciones, la fórmula es incorrecta.
Aquí está una geodésica en una proyección Equirectangular. Está claro que no es una onda sinusoidal:
(La imagen de fondo está tomada de http://upload.wikimedia.org/wikipedia/commons/thumb/e/ea/Equirectangular-projection.jpg/800px-Equirectangular-projection.jpg .)
Como todas las proyecciones equirectangulares son transformaciones afines de ésta (donde la coordenada x es la longitud y la coordenada y es la latitud), y las transformaciones afines de las ondas sinusoidales siguen siendo ondas sinusoidales, no podemos esperar que ninguna geodésica en ninguna forma de la proyección equirectangular sea una onda sinusoidal (excepto el Ecuador, que se traza como una línea horizontal). Así pues, empecemos por el principio y elaboremos la fórmula correcta.
Sea la ecuación de dicha geodésica de la forma
latitude = f(longitude)
para una función f para ser encontrado. (Este enfoque ya ha renunciado a los meridianos, que no pueden escribirse de esta forma, pero por lo demás es totalmente general). Convirtiendo a coordenadas cartesianas 3D (x,y,z) se obtiene
x = cos(l) cos(f(l))
y = sin(l) cos(f(l))
z = sin(f(l))
donde l es la longitud y se supone un radio unitario (sin pérdida de generalidad). Como las geodésicas en la esfera son intersecciones con planos (que pasan por su centro), debe existir un constante vector (a,b,c)--que se dirige entre los polos de la geodésica--para la cual
a x + b y + c z = 0
sin importar el valor de l podría ser. Resolviendo para f(l) se obtiene
f(l) = ArcTan(-(a cos(l) + b sin(l)) / c)
proporcionado c es distinto de cero. Evidentemente, cuando c se acerca a 0, obtenemos en el límite a par de meridianos que difieren en 180 grados, precisamente las geodésicas que abandonamos al principio. Así que todo está bien. Por cierto, a pesar de las apariencias, esto utiliza sólo dos parámetros iguales a a/c y b/c.
Obsérvese que todas las geodésicas pueden girar hasta cruzar el ecuador a cero grados de longitud. Esto indica que f(l) puede escribirse en términos de f0(l-l0) donde l0 es la longitud del cruce ecuatorial y f0 es la expresión de una geodésica que cruza en el meridiano primario. De aquí se obtiene la fórmula equivalente
f(l) = ArcTan(gamma * sin(l - l0))
donde -180 <= l0 < 180 grados es la longitud del cruce ecuatorial (cuando la geodésica entra en el hemisferio norte al viajar hacia el este) y gamma es un número real positivo. Esto no incluye los pares de meridianos. Cuando gamma \= 0 designa el Ecuador con un punto de partida en la longitud l0; podemos tomar siempre l0 = 0 en ese caso si deseamos una parametrización única. Todavía hay sólo dos parámetros, dados por l0 y gamma esta vez.
Mathematica 8.0 se utilizó para crear la imagen. De hecho, creó una "manipulación dinámica" en la que se puede controlar el vector (a,b,c) y se muestra al instante la geodésica correspondiente. (Eso es muy chulo.) Primero obtenemos la imagen de fondo:
i = Import[
"http://upload.wikimedia.org/wikipedia/commons/thumb/e/ea/\
Equirectangular-projection.jpg/800px-Equirectangular-projection.jpg"]
Aquí está el código en su totalidad:
Manipulate[
{a, b, c} = {Cos[u] Cos[v], Sin[u] Cos[v], Sin[v]};
Show[Graphics[{Texture[i],
Polygon[{{-\[Pi], -\[Pi]/2}, {\[Pi], -\[Pi]/2}, {\[Pi], \[Pi]/2}, {-\[Pi], \[Pi]/2}},
VertexTextureCoordinates -> {{0, 0}, {1, 0}, {1, 1}, {0, 1}}]}],
Plot[ArcTan[(a Cos[\[Lambda]] + b Sin[\[Lambda]]) / (-c)], {\[Lambda], -\[Pi], \[Pi]},
PlotRange -> {Automatic, {-\[Pi]/2, \[Pi]/2}}, PlotStyle -> {Thick, Red}]],
{u, 0, 2 \[Pi]}, {v, -\[Pi]/2, \[Pi]}]