If you don't mind ignoring the slight oblateness of the Earth (and your posted Haversine code does just that anyway) consider pre-converting all of your spherical (lat/long) coordinates into 3D unit-length cartesian coordinates first, per:
http://en.wikipedia.org/wiki/Spherical_coordinate_system
Then your spherical distance between cartesian coordinates p1
and p2
is simply:
r * acos(p1 . p2)
Since p1
and p2
will have unit length this reduces to four multiplications, two additions and one inverse trig operation per pair.
Also note that the calculation of dot products is an ideal candidate for optimisation, e.g. via GPU, MMX extensions, vector libraries, etc.
Furthermore, if your intent is to order the pairs by distance, potentially ignoring more distant pairs, you can defer the expensive r*acos()
part of the equation by sorting the list just on the dot product value since for all valid inputs (i.e. the range [-1, 1]
) it's guaranteed that:
acos(x) < acos(y) if x > y
You then just take the acos()
of the values you're actually interested in.
Re: the potential inaccuracies with using acos()
, those are really only significant if you're using single-precision float
variables. Using a double
with 16 significant digits should get you distances accurate to within one metre or less.
与恶龙缠斗过久,自身亦成为恶龙;凝视深渊过久,深渊将回以凝视…