Trigonometría esférica y coordenadas celestes

1. El triángulo esférico

Sobre la esfera unidad, un círculo máximo es la intersección con un plano que pasa por el centro $O$; es la geodésica de la esfera. Tres puntos $A, B, C$ no situados en un mismo círculo máximo, unidos por arcos de círculo máximo menores que $\pi$, forman un triángulo esférico.

• Los lados $a, b, c$ son ángulos centrales: $a = \angle BOC$, $b = \angle COA$, $c = \angle AOB$.
• Los ángulos $A, B, C$ son los ángulos diedros entre los planos de los lados que concurren en cada vértice (equivalentemente, el ángulo entre las tangentes a los dos lados en el vértice).

A diferencia del triángulo plano, $A + B + C > \pi$; el exceso esférico $E = A + B + C - \pi$ es igual al área del triángulo (teorema de Girard).

Escribimos $\hat{\mathbf a}, \hat{\mathbf b}, \hat{\mathbf c}$ para los vectores unitarios $\overrightarrow{OA}, \overrightarrow{OB}, \overrightarrow{OC}$. Entonces $\hat{\mathbf b}\cdot\hat{\mathbf c} = \cos a$, $\hat{\mathbf c}\cdot\hat{\mathbf a} = \cos b$ y $\hat{\mathbf a}\cdot\hat{\mathbf b} = \cos c$.

2. Fórmula de los cosenos

Teorema.

Demostración. La tangente en $A$ al arco $AB$ es la componente de $\hat{\mathbf b}$ perpendicular a $\hat{\mathbf a}$, normalizada:

donde los denominadores son las normas: $\lVert\hat{\mathbf b} - \cos c\,\hat{\mathbf a}\rVert^2 = 1 - 2\cos^2 c + \cos^2 c = \sin^2 c$. El ángulo $A$ es el que forman $\mathbf u$ y $\mathbf v$:

Para lados pequeños, $\cos x \approx 1 - x^2/2$ y $\sin x \approx x$ dan $a^2 \approx b^2 + c^2 - 2bc\cos A$: la ley del coseno plana (ejercicio 1).

3. Fórmula de los senos

Teorema.

Demostración. De la fórmula de los cosenos,

Desarrollando el numerador con $\sin^2 x = 1 - \cos^2 x$:

que es simétrico en $a, b, c$. Por tanto

no depende de qué vértice elijamos. Como todos los ángulos y lados están en $(0, \pi)$, los senos son positivos y se puede extraer la raíz. $\blacksquare$

(Geométricamente, $N = [\hat{\mathbf a}\,\hat{\mathbf b}\,\hat{\mathbf c}]^2$, el cuadrado del producto mixto: el volumen del paralelepípedo, que no privilegia ningún vértice.)

Fórmula de los cinco elementos. Aplicando la de los cosenos a $b$ y sustituyendo $\cos c$ de nuevo se obtiene

que usaremos para el acimut.

3.1 El triángulo polar y la fórmula de los cosenos para los ángulos

Dado un triángulo $ABC$, su triángulo polar $A'B'C'$ tiene por vértices los polos de los lados: $A'$ es el polo del círculo máximo $BC$ situado del mismo lado que $A$, etc. Se demuestra (Todhunter, cap. II) que sus lados y ángulos son suplementarios de los ángulos y lados del original: $a' = \pi - A$, $A' = \pi - a$. Aplicando la fórmula de los cosenos al triángulo polar y sustituyendo:

Es la fórmula de los cosenos para los ángulos: permite obtener un lado conocidos los tres ángulos, algo imposible en el plano (donde los tres ángulos solo fijan la forma, no el tamaño). En la esfera el tamaño está determinado porque el exceso esférico es el área.

3.2 Condicionamiento numérico: cuándo falla la fórmula de los cosenos

La separación angular $d$ entre dos astros se calcula a menudo con $\cos d = \sin\delta_1\sin\delta_2 + \cos\delta_1\cos\delta_2\cos\Delta\alpha$ (ejercicio 2). Para $d$ pequeño es una fórmula mal condicionada: $\cos d \approx 1 - d^2/2$ se calcula con un error absoluto del orden del épsilon de máquina, $\epsilon \approx 1.1\times10^{-16}$ en doble precisión, y al invertir,

Para $d = 1''$ ($4.8\times10^{-6}$ rad) el error es $2\times10^{-11}$ rad $= 5\times10^{-6}\,''$, tolerable; pero para $d = 1$ mas (milisegundo de arco, la precisión de Gaia), $\delta d \approx 2\times10^{-8}$ rad $= 5$ mas: el error supera al resultado. La cura es una fórmula que no reste cantidades casi iguales. Con el semiverseno $\operatorname{hav}\,x = \sin^2(x/2) = (1 - \cos x)/2$, la fórmula de los cosenos se reescribe (ejercicio 8) como

en la que todos los términos son pequeños y positivos cuando $d$ lo es, y $d = 2\arcsin\sqrt{\operatorname{hav}\,d}$ se obtiene con error relativo $O(\epsilon)$. Equivalente y aún más robusta en todo el rango es la forma vectorial: con $\hat{\mathbf a}$ y $\hat{\mathbf b}$ los vectores unitarios, $d = \operatorname{atan2}\!\left(\lVert\hat{\mathbf a}\times\hat{\mathbf b}\rVert,\ \hat{\mathbf a}\cdot\hat{\mathbf b}\right)$, que es exacta al nivel $\epsilon$ tanto para $d \to 0$ como para $d \to \pi$ (donde el haversine, a su vez, se degrada). Es la que usan SOFA y Astropy.

4. El triángulo de posición

Para un observador en latitud $\varphi$ consideramos tres puntos de la esfera celeste:

• $P$, el polo norte celeste (tomando $\varphi$ con signo, negativo al sur, las fórmulas valen también en el hemisferio sur),
• $Z$, el cenit,
• $X$, el astro, con declinación $\delta$ y ángulo horario $H$ (medido hacia el oeste desde el meridiano).

Sus lados son

y sus ángulos en $P$ y en $Z$ son $\angle ZPX = H$ y $\angle PZX = 2\pi - A$ si el astro está al oeste ($0 < H < \pi$), con el acimut $A$ medido desde el norte hacia el este. (Si está al este, $\angle ZPX = 2\pi - H$ y $\angle PZX = A$; las fórmulas finales son las mismas.) Que la altura del polo sea $\varphi$ es inmediato: el arco $PZ$ es la colatitud.

4.1 Derivación de la transformación ecuatorial → horizontal

Altura. Fórmula de los cosenos para el lado $ZX$, opuesto a $P$:

Acimut, seno. Fórmula de los senos con los pares $(H, ZX)$ y $(2\pi - A, PX)$:

Acimut, coseno. Fórmula de los cinco elementos (o la de los cosenos para el lado $PX$, sustituyendo $\sin h$):

Dividiendo, y usando la función de dos argumentos para obtener el cuadrante correcto,

Estas son exactamente las fórmulas de equatorialToHorizontal en src/lib/math/coordinates.ts (Meeus ec. 13.5–13.6, con acimut desde el norte en lugar de desde el sur).

Forma matricial. Las tres ecuaciones equivalen a una rotación de ángulo $\tfrac{\pi}{2} - \varphi$ alrededor del eje este–oeste que lleva el polo al cenit (ejercicio 3); en la práctica se encadenan así las rotaciones de precesión, nutación y tiempo sidéreo.

5. Tiempo sidéreo

El ángulo horario de un astro es $H = \theta - \alpha$, donde $\theta$ es el tiempo sidéreo local: el ángulo horario del punto vernal. Se obtiene del de Greenwich sumando la longitud $\lambda$ (positiva al este):

El tiempo sidéreo medio de Greenwich en función de la fecha juliana $\mathrm{JD}$ (UT1; las escalas de tiempo se tratan en la lección siguiente) es (Meeus ec. 12.4)

con $T$ en siglos julianos desde J2000.0. El coeficiente $360.9856…°$ por día solar medio refleja que en un día solar la Tierra gira algo más de una vuelta respecto de las estrellas: el día sidéreo dura $86\,164.09\ \text{s} = 23\ \text{h}\ 56\ \text{min}\ 4.09\ \text{s}$. Por eso cada estrella sale unos 3 min 56 s antes cada noche, y el cielo de una misma hora va cambiando con las estaciones.

6. Salida y puesta de las estrellas

Imponiendo $h = 0$ en la ecuación de la altura:

$H_0 \in [0, \pi]$ es el semiarco diurno: la estrella está sobre el horizonte durante $2H_0$ (en tiempo sidéreo, $2H_0/15°$ horas). Sale cuando $\theta = \alpha - H_0$ y se pone cuando $\theta = \alpha + H_0$.

• Si $-\tan\varphi\tan\delta \le -1$, nunca se pone: circumpolar. Esto equivale a $|\varphi + \delta| \ge \pi/2$ (ejercicio 4).
• Si $-\tan\varphi\tan\delta \ge 1$, nunca sale.

Con $h = 0$ en la ecuación del coseno para el lado $PX$, $\sin\delta = \cos\varphi\cos A$: el acimut del orto es

y el ocaso es simétrico, en $2\pi - A_{\text{orto}}$. Para Sirio ($\delta = -16.72°$) desde Bogotá ($\varphi = 4.71°$): $\cos A = -0.2886$, $A = 106.8°$, entre el este y el sureste.

Para el Sol y la Luna hay que usar $h_0 = -0.833°$ (refracción de $34'$ más semidiámetro de $16'$) o $h_0 \approx +0.125°$ para la Luna (que además tiene paralaje de $\sim 57'$). La fórmula general es $\cos H_0 = (\sin h_0 - \sin\varphi\sin\delta)/(\cos\varphi\cos\delta)$ (Meeus cap. 15).

7. La proyección estereográfica del planisferio

El simulador dibuja el hemisferio visible proyectándolo desde el nadir $N$ sobre el plano del horizonte (que pasa por el centro $O$ de la esfera). Un punto de distancia cenital $z$ va a un punto a distancia

del centro, en la dirección de su acimut. Prueba: el ángulo central entre $Z$ y $X$ es $z$; el ángulo inscrito $\angle ZNX$ que subtiende el mismo arco vale $z/2$, y en el triángulo rectángulo $NOX'$ (con $|NO| = 1$) se tiene $r = \tan(z/2)$. El cenit va a $r = 0$ y el horizonte a $r = 1$. Como el cielo se mira desde abajo, el plano se dibuja con el norte arriba y el este a la izquierda: $x = -r\sin A$, $y = r\cos A$ (src/lib/math/proyeccion.ts).

7.1 Conformidad

En coordenadas $(z, A)$ la métrica de la esfera es $ds^2 = dz^2 + \sin^2 z\,dA^2$, y en el plano, en polares, $d\sigma^2 = dr^2 + r^2 dA^2$. Con $r = \tan(z/2)$:

Los factores de escala en la dirección radial y en la azimutal son iguales, $k(z) = \tfrac{1}{2}\sec^2(z/2)$, luego $d\sigma^2 = k^2\,ds^2$: la proyección es conforme (preserva ángulos, y por tanto la forma de figuras pequeñas: las constelaciones no se deforman, solo crecen hacia el borde, hasta $k = 1$ en el horizonte frente a $k = \tfrac12$ en el cenit).

7.2 Las circunferencias van a circunferencias

Teorema. La imagen de toda circunferencia de la esfera que no pasa por $N$ es una circunferencia del plano; si pasa por $N$, es una recta.

Demostración. Sea $\iota$ la inversión del espacio de centro $N$ y radio $\sqrt2$: $X \mapsto X'$ en la semirrecta $NX$ con $|NX|\,|NX'| = 2$. Una inversión de radio $k$ transforma una esfera que pasa por su centro, de diámetro $2\rho$, en el plano perpendicular a ese diámetro situado a distancia $k^2/(2\rho)$ del centro. Aquí $k^2 = 2$ y $2\rho = |NZ| = 2$, así que la esfera celeste va al plano a distancia 1 de $N$: el plano del horizonte por $O$. Como cada punto se mueve a lo largo de la recta $NX$, $\iota$ restringida a la esfera es la proyección estereográfica (comprobación: $Z \mapsto O$ porque $|NZ| = 2$, y el horizonte, a distancia $\sqrt2$ de $N$, queda fijo). Ahora bien:

1. una circunferencia de la esfera es la intersección de la esfera con un plano $\Pi$;
2. la inversión lleva $\Pi$ a una esfera que pasa por $N$ (o a sí mismo si $N \in \Pi$);
3. por tanto la imagen de la circunferencia es la intersección del plano del horizonte con una esfera (circunferencia) o con un plano (recta, cuando $N \in \Pi$). $\blacksquare$

Consecuencias para el planisferio: los círculos de altura constante son circunferencias centradas; el ecuador celeste y los paralelos de declinación (circunferencias de la esfera) son arcos de circunferencia, que es lo que se ve al animar el simulador: cada estrella recorre un arco de circunferencia alrededor de la imagen del polo.

Problemas

1. Demuestra que para un triángulo esférico pequeño la fórmula de los cosenos se reduce a la ley del coseno plana, y estima el error relativo para lados de $10°$.

2. Calcula la separación angular entre Ácrux ($\alpha = 12^{\mathrm h}26^{\mathrm m}36^{\mathrm s}$, $\delta = -63°05'57''$) y Gácrux ($\alpha = 12^{\mathrm h}31^{\mathrm m}10^{\mathrm s}$, $\delta = -57°06'48''$).

3. Escribe la rotación que lleva el sistema ecuatorial local (ejes: intersección meridiano–ecuador hacia el sur, oeste, polo norte) al horizontal (sur, oeste, cenit) y recupera de ella las tres fórmulas de §4.1.

4. Demuestra que $-\tan\varphi\tan\delta \le -1 \iff |\varphi + \delta| \ge \pi/2$ para $\varphi, \delta \in (-\pi/2, \pi/2)$, y concluye qué estrellas son circumpolares en Bogotá ($\varphi = 4.71°$).

5. Calcula la duración del día en Bogotá en el solsticio de junio ($\delta_\odot = +23.44°$), primero con el horizonte matemático y luego con $h_0 = -0.833°$.

6. Demuestra que en la proyección estereográfica el área se amplifica en un factor $k^2 = \tfrac14\sec^4(z/2)$ y calcula cuánto más grande se ve una constelación en el horizonte que en el cenit.

7. (Simulador) En Bogotá, pon la hora en que Sirio sale y comprueba con el puntero que su acimut está cerca de $107°$. Luego cambia a Buenos Aires ($\varphi = -34.6°$) y predice primero con la fórmula el nuevo acimut del orto.

8. (Computacional) Demuestra la identidad del semiverseno de §3.2 a partir de la fórmula de los cosenos y comprueba numéricamente el condicionamiento: calcula la separación entre $(\alpha, \delta) = (0, 0)$ y $(\alpha, \delta) = (x, 0)$ para $x = 1°, 1', 1'', 1\ \text{mas}, 1\ \mu\text{as}$ con las tres fórmulas (cosenos, haversine, vectorial) y tabula el error relativo frente al valor exacto $x$.

import numpy as np
def cosenos(a1, d1, a2, d2):
    return np.arccos(np.sin(d1)*np.sin(d2) + np.cos(d1)*np.cos(d2)*np.cos(a2 - a1))
def haversine(a1, d1, a2, d2):
    h = np.sin((d2-d1)/2)**2 + np.cos(d1)*np.cos(d2)*np.sin((a2-a1)/2)**2
    return 2*np.arcsin(np.sqrt(h))
def vectorial(a1, d1, a2, d2):
    u = np.array([np.cos(d1)*np.cos(a1), np.cos(d1)*np.sin(a1), np.sin(d1)])
    v = np.array([np.cos(d2)*np.cos(a2), np.cos(d2)*np.sin(a2), np.sin(d2)])
    return np.arctan2(np.linalg.norm(np.cross(u, v)), np.dot(u, v))
for nombre, x in [("1°", np.radians(1)), ("1'", np.radians(1/60)), ('1"', np.radians(1/3600)),
                  ("1 mas", np.radians(1e-3/3600)), ("1 µas", np.radians(1e-6/3600))]:
    print(nombre, [abs(f(0, 0, x, 0) - x)/x for f in (cosenos, haversine, vectorial)])

Lecturas

• Smart (1977), Textbook on Spherical Astronomy, caps. I–II: el desarrollo clásico de las fórmulas y del triángulo de posición. Todhunter (1886): triángulo polar y todas las identidades, con demostraciones completas.
• Green (1985), Spherical Astronomy, caps. 1–3: tratamiento moderno con vectores y matrices de rotación.
• Meeus (1998), caps. 12, 13 y 15: las fórmulas tal como se programan, con ejemplos numéricos.
• Snyder (1987), §21, y Needham (1997), cap. 3: la proyección estereográfica desde la cartografía y desde el análisis complejo.
• Sinnott (1984), Sky & Telescope 68, 159: el artículo que popularizó el haversine en astronomía computacional.