Which one of these azimuth values is the right one? And why?

Question

First of all, let's make clear that I want the Azimuth on the surface of the Earth, i.e. the angle between two locations, for example New York and Moscow.

I am testing some azimuth calculations with my JS functions (shown below). For the points A(-170, -89) to B(10, 89), I get ~90º.

JS function for Azimuth on sphere (from Wikipedia)

var dLon = lon2 - lon1;
var y = Math.sin(dLon) * Math.cos(lat2);
var x = Math.cos(lat1) * Math.sin(lat2) - Math.sin(lat1) * Math.cos(lat2) * Math.cos(dLon);
var angle = Math.atan2(y, x) * 180 / Math.PI;

JS function for Azimuth on oblate spheroid (from Wikipedia)

var dLon = lon2 - lon1;
var f = 1 / 298.257223563;    /* Flattening for WGS 84 */
var b = (1 - f) * (1 - f);
var tanLat2 = Math.tan(lat2);
var y = Math.sin(dLon);
var x;
if (lat1 === 0) {
    var x = b * tanLat2;
} else {
    var a = f * (2 - f);
    var tanLat1 = Math.tan(lat1);
    var c = 1 + b * tanLat2 * tanLat2;
    var d = 1 + b * tanLat1 * tanLat1;
    var t = b * Math.tan(lat2) / Math.tan(lat1) + a * Math.sqrt(c / d);
    var x = (t - Math.cos(dLon)) * Math.sin(lat1);
}
var angle = Math.atan2(y, x) * 180 / Math.PI;

In Calculator 2, I get 90º.

In PostGIS, I get 270º

In Calculator 1, I get 180º.

I know the Azimuth gets more and more distorted near the Poles, but that's exactly why I am testing at these spots. This variety of different solutions are confusing me. Could you please help me getting the right answer for this?

Answer 1

It depends on the reference used for azimuth, e.g. map-types use 0° for North and positive is clockwise, while math-types uses 0° for East and positive is anticlockwise.

The pair of coordinates A(-170, -89) and B(10, 89) are antipodes which are a special case for finding minimum distances and azimuths. Your question can be answered with a thought exercise.

First note that the half-circumferences of the earth are:

• Equatorial: 20037.5085 km
• Meridional (north-to-south): 20003.93 km

For a pair of antipodes on the north and south poles, there are an infinite number of azimuths, since the distance is the same along each longitude. (What direction would you go from the south pole to the north pole?)

For a pair of antipodes on the equator, the shortest distances are north or south, since it is slightly shorter along the meridional direction.

For any other pair of antipodes, it is the same answer as from the equator: north or south.

---

Update

To investigate the problem a bit more with the PostGIS SQL query:

SELECT ST_Distance(A, B), degrees(ST_Azimuth(A, B))
FROM (
  SELECT 'POINT(-170 -89)'::geography A, 'POINT(10 89)'::geography B
) f;

With PostGIS 2.0 and 2.1, the incorrect results are:

   st_distance   |     degrees
-----------------+------------------
 20003900.583699 | 270.005278779849

But with PostGIS 2.2 (and PROJ 4.9.1), the correct results are now:

   st_distance    | degrees
------------------+---------
 20003931.4586255 |     180