MySQL Forums

Re: WGS84 Data & distance queries
Posted by: Pierre Herman
Date: September 01, 2007 02:33PM

Here is a function in PHP that returns the distance between two points on earth based on the WGS84 datum. This is nothing new. It is called the Vincenty's folmula and the theory can be found at http://www.ngs.noaa.gov/PUBS_LIB/inverse.pdf

This is maybe not what you're looking for but it can be easily translated in any other language and is accurate to less than a meter. You just need to send the latitude and longitude in decimal degrees.

function distInvVincenty(\$lat1, \$lon1, \$lat2, \$lon2) {
/* WGS84 stuff */
\$a = 6378137;
\$b = 6356752.3142;
\$f = 1/298.257223563;
/* end of WGS84 stuff */
\$sinU1 = sin(\$U1);
\$cosU1 = cos(\$U1);
\$sinU2 = sin(\$U2);
\$cosU2 = cos(\$U2);

\$lambda = \$L;
\$lambdaP = 2*pi();
\$iterLimit = 20;
while ((abs(\$lambda-\$lambdaP) > pow(10, -12)) && (\$iterLimit-- > 0)) {
\$sinLambda = sin(\$lambda);
\$cosLambda = cos(\$lambda);
\$sinSigma = sqrt((\$cosU2*\$sinLambda) * (\$cosU2*\$sinLambda) + (\$cosU1*\$sinU2-\$sinU1*\$cosU2*\$cosLambda) * (\$cosU1*\$sinU2-\$sinU1*\$cosU2*\$cosLambda));
if (\$sinSigma==0) {return 0;}
\$cosSigma = \$sinU1*\$sinU2 + \$cosU1*\$cosU2*\$cosLambda;
\$sigma = atan2(\$sinSigma, \$cosSigma);
\$sinAlpha = \$cosU1 * \$cosU2 * \$sinLambda / \$sinSigma;
\$cosSqAlpha = 1 - \$sinAlpha*\$sinAlpha;
\$cos2SigmaM = \$cosSigma - 2*\$sinU1*\$sinU2/\$cosSqAlpha;
if (is_nan(\$cos2SigmaM)) {\$cos2SigmaM = 0;}
\$C = \$f/16*\$cosSqAlpha*(4+\$f*(4-3*\$cosSqAlpha));
\$lambdaP = \$lambda;
\$lambda = \$L + (1-\$C) * \$f * \$sinAlpha *(\$sigma + \$C*\$sinSigma*(\$cos2SigmaM+\$C*\$cosSigma*(-1+2*\$cos2SigmaM*\$cos2SigmaM)));
}
if (\$iterLimit==0) {return NaN;} // formula failed to converge

\$uSq = \$cosSqAlpha * (\$a*\$a - \$b*\$b) / (\$b*\$b);
\$A = 1 + \$uSq/16384*(4096+\$uSq*(-768+\$uSq*(320-175*\$uSq)));
\$B = \$uSq/1024 * (256+\$uSq*(-128+\$uSq*(74-47*\$uSq)));
\$deltaSigma = \$B*\$sinSigma*(\$cos2SigmaM+\$B/4*(\$cosSigma*(-1+2*\$cos2SigmaM*\$cos2SigmaM)- \$B/6*\$cos2SigmaM*(-3+4*\$sinSigma*\$sinSigma)*(-3+4*\$cos2SigmaM*\$cos2SigmaM)));

\$s = \$b*\$A*(\$sigma-\$deltaSigma);
return "Distance: ".\$s;
}

Subject
Views
Written By
Posted
9490
April 25, 2006 07:26PM
6704
December 19, 2006 03:42AM
5893
April 18, 2007 01:07AM
4867
May 09, 2007 03:47PM
4682
May 17, 2007 12:34PM
4471
May 21, 2007 09:13AM
Re: WGS84 Data & distance queries
9407
September 01, 2007 02:33PM

Sorry, you can't reply to this topic. It has been closed.

Content reproduced on this site is the property of the respective copyright holders. It is not reviewed in advance by Oracle and does not necessarily represent the opinion of Oracle or any other party.