Как рассчитать азимут (угол на север) между двумя координатами WGS84?

StackOverflow https://stackoverflow.com/questions/642555

  •  22-07-2019
  •  | 
  •  

Вопрос

У меня есть две координаты WGS84: широта и долгота в градусах.Эти точки довольно близки друг к другу, например.всего один метр друг от друга.

Есть ли простой способ вычислить азимут линии между этими точками, то есть угол на север?

Наивным подходом было бы принять декартову систему координат (поскольку эти точки расположены очень близко друг к другу) и просто использовать

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

a = азимут L1, L2 = долгота B1, B2 = широта

Ошибка будет больше по мере удаления координат от экватора, поскольку там расстояние между двумя градусами долготы становится все меньше, чем расстояние между двумя градусами широты (которое остается постоянным).

Я нашел несколько довольно сложных формул, которые на самом деле не хочу реализовывать, потому что они кажутся излишними для точек, которые расположены так близко друг к другу, и мне не нужна очень высокая точность (двух десятичных знаков достаточно, один, вероятно, тоже подойдет, поскольку есть есть и другие факторы, которые в любом случае снижают точность, например, тот, который возвращает GPS).

Возможно, я мог бы просто определить приблизительный поправочный коэффициент по долготе в зависимости от широты и использовать что-то вроде этого:

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

где f — поправочный коэффициент

Есть какие-нибудь подсказки?

(Я не хочу использовать для этого какие-либо библиотеки, особенно те, которые требуют лицензий времени выполнения.Подойдет любой исходный код Delphi с MPL.)

Это было полезно?

Решение

Формулы, на которые вы ссылаетесь в тексте, предназначены для расчета расстояния по большому кругу между двумя точками.Вот как я вычисляю угол между точками:

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;

  • Не забудьте обработать ситуацию, когда 2 точки равны (проверьте, равен ли результат cNO_ANGLE, или измените функцию, чтобы она выдавала исключение);

  • Эта функция предполагает, что вы находитесь на ровной поверхности.С небольшими расстояниями, о которых вы упомянули, все в порядке, но если вы собираетесь рассчитывать направление между городами по всему миру, возможно, вам захочется изучить что-то, что учитывает форму Земли;

  • Лучше всего предоставить этой функции координаты, которые уже сопоставлены с плоской поверхностью.Однако вы можете подать широту WGS84 непосредственно в Y (а lon в X), чтобы получить грубое приближение.

Другие советы

Вот решение C #. Проверено на углы 0, 45, 90, 135, 180, 225, 270 и 315.

Изменить . Я заменил свое предыдущее уродливое решение переводом Wouter на C #: 

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);
}

Я нашел эту ссылку

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

дано в ответе на

Lat / Lon + Distance + Heading - > Lat / Lon

Это выглядит многообещающе, особенно в приближении плоской земли, приведенном ближе к концу.

Это будет работать только для небольших различий. В противном случае вы не можете просто "latitudeDifference / longticalDifference".

Я бы рекомендовал ввести поправочный коэффициент на основе долготы. Я однажды реализовал симуляционную процедуру, чтобы вернуть все геокодированные записи в пределах x миль от конкретного места, и столкнулся с симуляцией. К сожалению, у меня больше нет кода, и я не могу вспомнить, как я добрался до номера коррекции, но вы на правильном пути.

Кто-нибудь проверял это? Он не возвращает правильные ответы

Эта функция предполагает, что вы находитесь на плоской поверхности. С небольшими расстояниями, о которых вы упомянули, это все хорошо, но если вы собираетесь рассчитывать курс между городами по всему миру, вам, возможно, захочется взглянуть на что-то, что примет форму земли;

Ваша плоская земля не имеет к этому никакого отношения. Ошибка, как вы это называете, заключается в том, что вы вычисляете начальный азимут из точки. Если вы не направляетесь прямо к полюсу, ваше отношение к полюсу будет меняться с расстоянием. Несмотря на это, вышеприведенная программа не возвращает правильные результаты.

Лицензировано под: CC-BY-SA с атрибуция
Не связан с StackOverflow
scroll top