Take the 2-minute tour ×
Geographic Information Systems Stack Exchange is a question and answer site for cartographers, geographers and GIS professionals. It's 100% free, no registration required.

I am attempting to calculate the distance between two latitude/longitude points. I have a piece of code that mostly works that I yanked from this post but I do not really understand how it works.

Here is the code:

<?php
// POINT 1
$thisLat = deg2rad(44.638);
$thisLong = deg2rad(-63.587);

// POINT 2
$otherLat = deg2rad(44.644);
$otherLong = deg2rad(-63.911);

$MeanRadius = 6378 - 21 * sin($lat1);

$xa = (Cos($thisLat)) * (Cos($thisLong));
$ya = (Cos($thisLat)) * (Sin($thisLong));
$za = (Sin($thisLat));

$xb = (Cos($otherLat)) * (Cos($otherLong));
$yb = (Cos($otherLat)) * (Sin($otherLong));
$zb = (Sin($otherLat));

$distance = $MeanRadius * Acos($xa * $xb + $ya * $yb + $za * $zb);

echo $distance;
?>

I have a couple questions:

  1. what are xa, ya, za? I understand that they are points on a 3D cartesian plane but where are they relative to? The center of the earth?
  2. How does this cos($xa * $xb + $ya * $yb + $za * $zb) calculate the distance between the points? I know that in 2D I would do this:

alt text

Pythagorean Theorem 
distance^2 = b^2 + a^2
distance = sqr((y2-y1)^2 + (x2 - x1)^2)
  1. How accurate will this be? There was some discussion about that on the other page. But I specifically want to use the distance to tell if users are within the something like 10m, 20m or 50m of each other. Will I be able to do this with good accuracy?
  2. What should I use for $MeanRadius? Is that a reasonable value? I think that that value assumes that the earth is a elipse.
share|improve this question
add comment

2 Answers 2

up vote 14 down vote accepted

This is terrible code for general-purpose use because it can give erroneous results or even fail altogether for short distances. Use the Haversine Formula instead.

(The formula on which your code is based converts two points on the sphere (not an ellipse) into their 3D Cartesian coordinates (xa,ya,za) and (xb,yb,zb) on the unit sphere and forms their dot product, which will equal the cosine of the angle between them. The ACos function returns that angle, which when scaled by the earth's radius will estimate the distance. The problem is that the cosine of a small angle, say of size 'e' in radians, differs from 1 by an amount close to e^2/2. This disappears into the floating point error cloud when e is smaller than the square root of twice the floating point precision. If you're computing in single precision this means values of e less than 0.001 -- about one kilometer -- will be confused with zero! In double precision the cutoff is around e = 10^-8, but by the time e = 10^-4 or so (about 10 meters) you potentially can lose so much precision that you need to worry, unless the implementation of the ACos function is unusually accurate (e.g., has some high-precision internal computations built in)).

share|improve this answer
add comment

Take a look at this site, http://www.movable-type.co.uk/scripts/latlong.html

If you implement it and you get the wrong answers you probably have the wrong units. I think for that site most operations are done in radians instead of decimal degrees.

share|improve this answer
    
+1 Nice reference. –  whuber Nov 30 '10 at 15:55
add comment

Your Answer

 
discard

By posting your answer, you agree to the privacy policy and terms of service.

Not the answer you're looking for? Browse other questions tagged or ask your own question.