Longitud de un grado: ¿de dónde vienen los términos en esta fórmula?

Question

Las calculadoras en línea, tales como http://www.csgnetwork.com/degreelenllavcalc.html (ver código fuente de la página) utilice el siguiente fórmulas para obtener metros por cada grado. Entiendo que, en general, cómo la distancia por grado varía dependiendo de la latitud de la ubicación, pero no entiendo cómo eso se traduce en que los de abajo. Más específicamente, donde se hacen las constantes, las 3 "cos" en cada fórmula, y los coeficientes de (2, 4, 6; 3, y 5) para "lat"?

    // Set up "Constants"
    m1 = 111132.92;     // latitude calculation term 1
    m2 = -559.82;       // latitude calculation term 2
    m3 = 1.175;         // latitude calculation term 3
    m4 = -0.0023;       // latitude calculation term 4
    p1 = 111412.84;     // longitude calculation term 1
    p2 = -93.5;         // longitude calculation term 2
    p3 = 0.118;         // longitude calculation term 3

    // Calculate the length of a degree of latitude and longitude in meters
    latlen = m1 + (m2 * Math.cos(2 * lat)) + (m3 * Math.cos(4 * lat)) +
            (m4 * Math.cos(6 * lat));
    longlen = (p1 * Math.cos(lat)) + (p2 * Math.cos(3 * lat)) +
                (p3 * Math.cos(5 * lat));

Answer 1

El director de radio del esferoide WGS84 es un = 6378137 metros y su inversa de acoplamiento es f = 298.257223563, donde el cuadrado de la excentricidad es

e2 = (2 - 1/f)/f = 0.0066943799901413165.

La parte meridional de radio de curvatura en la latitud phi es

M = a(1 - e2) / (1 - e2 sin(phi)^2)^(3/2)

y el radio de curvatura a lo largo del paralelo es

N = a / (1 - e2 sin(phi)^2)^(1/2)

Además, el radio del paralelo es

r = N cos(phi)

Estos son multiplicativos correcciones a la forma esférica de los valores de M y N, ambos de los cuales la igualdad de la esférica de radio de una, que es lo que se reduce a cuando e2 = 0.

En el punto amarillo en 45 grados de latitud norte, el azul del disco de radio de M es el osculating ("besos") círculo en el sentido de los meridianos y el disco rojo de la radio N es el osculating círculo en el sentido de los paralelos: ambos discos contienen el "abajo" de la dirección en este punto. Esta cifra se exagera el aplanamiento de la tierra por dos órdenes de magnitud.

Los radios de curvatura de determinar las longitudes de los grados: cuando un círculo tiene un radio de R, su perímetro de longitud 2 pi R cubre 360 grados, donde la longitud de un grado es pi * R / 180. Sustituyendo M y r para R -- es decir, la multiplicación de M y r por pi/180-da a simple exacto de las fórmulas para obtener el grado de longitud.

Estas fórmulas-que se basan exclusivamente en los valores dados de una y f (que se pueden encontrar en muchos lugares) y la descripción del esferoide como un elipsoide de rotación--de acuerdo con los cálculos en la pregunta dentro de 0.6 partes por millón (unos pocos centímetros), que es aproximadamente del mismo orden de magnitud de los más pequeños coeficientes en la pregunta, lo que indica que está de acuerdo. (La aproximación es siempre un poco baja.) En la parcela que el error relativo de la longitud de un grado de latitud es de color negro y de longitud es discontinua roja:

En consecuencia, podemos entender los cálculos en la pregunta son aproximaciones (a través de la serie truncada trigonométricas) para las fórmulas indicadas anteriormente.

---

Los coeficientes puede ser calculado a partir del coseno de Fourier de la serie de M y r en función de la latitud. Ellos se dan en términos de funciones elípticas de e2, que sería demasiado complicado de reproducir aquí. Para el esferoide WGS84, mis cálculos dar

  m1 = 111132.95255
  m2 = -559.84957
  m3 = 1.17514
  m4 = -0.00230
  p1 = 111412.87733
  p2 = -93.50412
  p3 = 0.11774
  p4 = -0.000165

(Usted puede adivinar cómo p4 entra en la fórmula. :) La cercanía de estos valores para los parámetros en el código da fe de la exactitud de esta interpretación. Esta mejora de la aproximación es exacta y mucho mejor que la de una parte por mil millones en todas partes.

---

Para probar esta respuesta me ejecutado R código para llevar a cabo los cálculos:

#
# Radii of meridians and parallels on a spheroid.  Defaults to WGS84 meters.
# Input is latitude (in degrees).
#
radii <- function(phi, a=6378137, e2=0.0066943799901413165) {
  u <- 1 - e2 * sin(phi)^2
  return(cbind(M=(1-e2)/u, r=cos(phi)) * (a / sqrt(u))) 
}
#
# Approximate calculation.  Same interface (but no options).
#
m.per.deg <- function(lat) {
  m1 = 111132.92;     # latitude calculation term 1
  m2 = -559.82;       # latitude calculation term 2
  m3 = 1.175;         # latitude calculation term 3
  m4 = -0.0023;       # latitude calculation term 4
  p1 = 111412.84;     # longitude calculation term 1
  p2 = -93.5;         # longitude calculation term 2
  p3 = 0.118;         # longitude calculation term 3

  latlen = m1 + m2 * cos(2 * lat) + m3 * cos(4 * lat) + m4 * cos(6 * lat);
  longlen = p1 * cos(lat) + p2 * cos(3 * lat) + p3 * cos(5 * lat);
  return(cbind(M.approx=latlen, r.approx=longlen))
}
#
# Compute the error of the approximation `m.per.deg` compared to the 
# correct formula and plot it as a function of latitude.
#
phi <- pi / 180 * seq(0, 90, 10)
names(phi) <- phi * 180 / pi
matplot(phi * 180 / pi, 10^6 * ((m.per.deg(phi) - radii(phi) * pi / 180) / 
       (radii(phi) * pi / 180)),
        xlab="Latitude (degrees)", ylab="Relative error * 10^6",lwd=2, type="l")

El cálculo exacto con radii puede ser utilizado para imprimir las tablas de las longitudes de los grados, como en

zapsmall(radii(phi) * pi / 180)

La salida es en metros y se parece a esto (con algunas líneas eliminado):

          M         r
0  110574.3 111319.49
10 110607.8 109639.36
20 110704.3 104647.09
...
80 111659.9  19393.49
90 111694.0      0.00

---

Referencias

LM Bugayevskiy y JP Snyder, Proyecciones de Mapa--Un Manual de Referencia. Taylor & Francis, 1995. (Anexo 2 y Anexo 4)

JP Snyder, Proyecciones de Mapa--Un Manual de Trabajo. USGS Profesional de Papel 1395, 1987. (Capítulo 3)

Answer 2

Esa es la fórmula del Haversine, aunque expresado de forma extraña.