¿Cómo calculo el acimut (ángulo hacia el norte) entre dos coordenadas WGS84

Question

Tengo dos coordenadas WGS84, latitud y longitud en grados. Estos puntos están bastante juntos, p. solo un metro de distancia.

¿Hay una manera fácil de calcular el acimut de la línea entre estos puntos, es decir, el ángulo hacia el norte?

El enfoque ingenuo sería asumir un sistema de coordenadas cartesianas (porque estos puntos están muy juntos) y simplemente usar

sin (a) = abs (L2-L1) / sqrt (sqr (L2-L1) + sqr (B2-B1))

a = acimut L1, L2 = longitud B1, B2 = latitud

El error será mayor a medida que las coordenadas se alejen del ecuador porque allí la distancia entre dos grados longitudinales se vuelve cada vez más pequeña que la que está entre dos grados latitudinales (que permanece constante).

Encontré algunas fórmulas bastante complejas que realmente no quiero implementar porque parecen ser excesivas para los puntos que están tan juntos y no necesito una precisión muy alta (dos decimales son suficientes, uno probablemente esté bien) ya que hay otros factores que reducen la precisión de todos modos, como el que devuelve el GPS).

Tal vez podría determinar un factor de corrección longitudinal aproximado dependiendo de la latitud y usar algo como esto:

sin (a) = abs (L2 * f-L1 * f) / sqrt (sqr (L2 * f-L1 * f) + sqr (B2-B1))

donde f es el factor de corrección

¿Alguna pista?

(No quiero usar ninguna biblioteca para esto, especialmente las que requieren licencias de tiempo de ejecución. Cualquier fuente MPLed Delphi sería genial).

Answer 1

Las fórmulas a las que se refiere en el texto son para calcular la distancia del gran círculo entre 2 puntos. Así es como calculo el ángulo entre puntos:

uses Math, ...;
...

const
  cNO_ANGLE=-999;

...

function getAngleBetweenPoints(X1,Y1,X2,Y2:double):double;
var
  dx,dy:double;
begin
  dx := X2 - X1;
  dy := Y2 - Y1;

  if (dx > 0) then  result := (Pi*0.5) - ArcTan(dy/dx)   else
  if (dx < 0) then  result := (Pi*1.5) - ArcTan(dy/dx)   else
  if (dy > 0) then  result := 0                          else
  if (dy < 0) then  result := Pi                         else
                    result := cNO_ANGLE; // the 2 points are equal

  result := RadToDeg(result);
end;

• Recuerde manejar la situación donde 2 puntos son iguales (verifique si el resultado es igual a cNO_ANGLE, o modifique la función para lanzar una excepción);

• Esta función asume que estás en una superficie plana. Con las pequeñas distancias que ha mencionado, todo está bien, pero si va a calcular el rumbo entre ciudades de todo el mundo, es posible que desee buscar algo que tome la forma de la tierra en cuenta;

• Es mejor proporcionar esta función con coordenadas que ya están asignadas a una superficie plana. Sin embargo, puede alimentar WGS84 Latitude directamente en Y (y lon en X) para obtener una aproximación aproximada.

Answer 2

Aquí está la solución C #. Probado para 0, 45, 90, 135, 180, 225, 270 y 315 ángulos.

Editar Reemplacé mi solución fea anterior, por la traducción en C # de la solución de Wouter:

public double GetAzimuth(LatLng destination)
{
    var longitudinalDifference = destination.Lng - this.Lng;
    var latitudinalDifference = destination.Lat - this.Lat;
    var azimuth = (Math.PI * .5d) - Math.Atan(latitudinalDifference / longitudinalDifference);
    if (longitudinalDifference > 0) return azimuth;
    else if (longitudinalDifference < 0) return azimuth + Math.PI;
    else if (latitudinalDifference < 0) return Math.PI;
    return 0d;
}

public double GetDegreesAzimuth(LatLng destination)
{
    return RadiansToDegreesConversionFactor * GetAzimuth(destination);
}

Answer 3

Encontré este enlace

http://williams.best.vwh.net/avform.htm

dado en la respuesta a

Lat / Lon + Distancia + Rumbo - > Lat / Lon

Esto parece prometedor, especialmente la aproximación de tierra plana dada cerca del final.

Answer 4

Esto funcionaría solo para pequeñas diferencias. De lo contrario, no puede simplemente '' latitudinalDifference / longitudinalDifference ''.

Answer 5

Recomendaría implementar un factor de corrección basado en la longitud. Implementé una rutina simular una vez para devolver todos los registros geocodificados dentro de x millas de un lugar específico y me encontré con problemas similares. Desafortunadamente, ya no tengo el código, y parece que no puedo recordar cómo llegué al número de corrección, pero estás en el camino correcto.

Answer 6

¿Alguien ha probado esto? No devuelve las respuestas correctas

Esta función asume que estás en una superficie plana. Con las pequeñas distancias que ha mencionado, todo está bien, pero si va a calcular el rumbo entre ciudades de todo el mundo, es posible que desee buscar algo que tome la forma de la tierra en cuenta;

Tu tierra plana tiene poco que ver con esto. El error, como lo llamas, es porque estás calculando un acimut inicial desde un punto. A menos que se dirija directamente a un poste, su relación con el poste cambiará con la distancia. En cualquier caso, el programa anterior no devuelve resultados correctos.