Solved: what is wrong with this LatLong->Distance code? A Trigonometry "feature".

[StephenG1CMZ]

Thank you for taking the time to look at this.
When I first saw your code I was thinking the lower-case Sin was exact and avoiding a rounding error - but it doesn't fix this example, which is still 1 km rather than 54 km.

For spherical calculations I usually use an average of polar and equatorial radius, (actually. (2a+b)/3), or repeat the calculations for both - but that detail isn't significant.

Example data:
http://www.ga.gov.au/geodesy/datums/vinc...nverse.jsp

Update: a missing atan2 call is added - the reported distance is still 1 km. Actually, that is probably just equivalent to the previous line.

Code:





EXPORT Haversine(LAT1,LON1,LAT2,LON2)
 BEGIN
  LOCAL RR:=6378.137;
  LOCAL DLON:=(LON2-LON1)/180*Pi;
  LOCAL DLAT:=(LAT2-LAT1)/180*Pi;
  LOCAL AA:=(sin(DLAT/2))^2 + 
   COS((LAT1)/180*Pi)*
   COS((LAT2)/180*Pi)*(sin(DLON/2))^2;
  LOCAL CC:=2*ASIN(MIN(1,√(AA)));//MIN 1 AVOIDS ROUNDING OUTRANGE
  //INSERT LINES
  LOCAL DD;
  CC:=2*ATAN2YX(√(AA),√(1-AA^2));
  //END INSERT
  DD:=RR*CC;//REMOVE LOCAL
  
  RETURN DD;
 END;

 EXPORT HAVERSINES()
 BEGIN
  LOCAL DIST;

  PRINT();
  DIST:=Haversine(−37°57′3.72030″,144°25′29.52440″,−37°39′10.15610″,143°55′35​.38390″);//EXPECTS 54972m
  PRINT(DIST+" km");
  PRINT("EXPECTED: 54972m");//VINCENTY AUSTRALIA
 
 END;