[ACCEPTED]-How do I calculate the Azimuth (angle to north) between two WGS84 coordinates-geospatial

The formulas that you refer to in the text 17 are to calculate the great circle distance 16 between 2 points. Here's how I calculate 15 the angle between points:

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;

• Remember to handle 14 the situation where 2 points are equal (check 13 if the result equals cNO_ANGLE, or modify 12 the function to throw an exception);

• This 11 function assumes that you're on a flat surface. With 10 the small distances that you've mentioned 9 this is all fine, but if you're going to 8 be calculating the heading between cities 7 all over the world you might want to look 6 into something that takes the shape of the 5 earth in count;

• It's best to provide this 4 function with coordinates that are already 3 mapped to a flat surface. You could feed 2 WGS84 Latitude directly into Y (and lon 1 into X) to get a rough approximation though.

Here is the C# solution. Tested for 0, 45, 90, 135, 180, 225, 270 3 and 315 angles.

Edit I replaced my previous ugly 2 solution, by the C# translation of Wouter's 1 solution:

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

I found this link

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

given in the answer to

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

This 2 looks promising, especially the flat earth 1 approximation given near the end.

This would work only for small differences. Otherwise 1 you can't just "latitudinalDifference / longitudinalDifference".

I would recommend implementing a correction 7 factor based on the longitude. I implemented 6 a simular routine once to return all geocoded 5 records within x miles of a specific spot 4 and ran into simular issues. Unfortunately 3 I don't have the code anymore, and can't 2 seem to recall how I got to the correction 1 number but you are on the right track.

Has anyone tested this? It does not return 7 the correct answers

This function assumes that you're on a flat surface. With the small distances that you've mentioned this is all fine, but if you're going to be calculating the heading between cities all over the world you might want to look into something that takes the shape of the earth in count;

Your flat earth has little 6 to do with this. The error, as you call 5 it, is because you are calculating an initial 4 azimuth from a point. Unless you are heading 3 straight to a pole, you relationship to 2 the pole will change with distance. Regardless, the 1 above program does not return correct results.