根据两个点经纬度计算距离 地理空间距离计算及优化( 二 )


=Math.sin((lat1-lat2)*0.5);
=hsinY*hsinY+
(Math.cos(lat1)*Math.cos(lat2)*hsinX*hsinX);
2*Math.atan2(Math.sqrt(h),Math.sqrt(1-h))*;
}
2)公式性能
根据2个经纬度点 , 计算这2个经纬度点之间的距离(通过经度纬度得到距离)
球面上任意两点之间的距离计算公式可以参考维基百科上的下述文章 。
值得一提的是 , 维基百科推荐使用公式 , 理由是Great- 公式用到了大量余弦函数 ,  而两点间距离很短时(比如地球表面上相距几百米的两点) , 余弦函数会得出0.999...的结果 ,  会导致较大的舍入误差 。而公式采用了正弦函数 , 即使距离很小 , 也能保持足够的有效数字 。以前采用三角函数表计算时的确会有这个问题 , 但经过实际验证 , 采用计算机来计算时 , 两个公式的区别不大 。稳妥起见 , 这里还是采用公式 。
其中
根据2个经纬度坐标 , 距离计算函数
下面就是计算球面间两点(lat1, lon1) - (lat2, lon2)之间距离的函数 。
using System;using System.Collections.Generic;using System.Linq;using System.Text;namespace HarvenSin{class Program{/// /// 根据经纬度 , 计算2个点之间的距离 。/// /// static void Main(string[] args){//39.94607,116.3279331.24063,121.42575Console.WriteLine(Distance(39.94607, 116.32793, 31.24063, 121.42575));}public static double HaverSin(double theta){var v = Math.Sin(theta / 2);return v * v;}static double EARTH_RADIUS = 6371.0;//km 地球半径 平均值 , 千米/// /// 给定的经度1 , 纬度1;经度2 , 纬度2. 计算2个经纬度之间的距离 。/// /// 经度1/// 纬度1/// 经度2/// 纬度2/// 距离(公里、千米)public static double Distance(double lat1,double lon1, double lat2,double lon2){//用haversine公式计算球面两点间的距离 。//经纬度转换成弧度lat1 = ConvertDegreesToRadians(lat1);lon1 = ConvertDegreesToRadians(lon1);lat2 = ConvertDegreesToRadians(lat2);lon2 = ConvertDegreesToRadians(lon2);//差值var vLon = Math.Abs(lon1 - lon2);var vLat = Math.Abs(lat1 - lat2);//h is the great circle distance in radians, great circle就是一个球体上的切面 , 它的圆心即是球心的一个周长最大的圆 。var h = HaverSin(vLat) + Math.Cos(lat1) * Math.Cos(lat2) * HaverSin(vLon);var distance = 2 * EARTH_RADIUS * Math.Asin(Math.Sqrt(h));return distance;}/// /// 将角度换算为弧度 。/// /// 角度/// 弧度public static double ConvertDegreesToRadians(double degrees){return degrees * Math.PI / 180;}public static double ConvertRadiansToDegrees(double radian){return radian * 180.0 / Math.PI;}}}
公式来历:
(F)=1-cos(F)
名字来历是Ha- , 即Half-  , 表示sin的一半的意思 。
hav(A) = (1-cos(A))/2 = sin(A/2)* sin(A/2)
推倒过程:
如下一个半径为1 的圆 , O是圆心 , A、B是弦(chord) 。角度AOB=theta 。则角度AOC=theta/2 。OC是垂直于AB的垂线() 。AC长度是sin(theta/2) , AB长度是2*sin(theta/2) 。
(图1)
如下地球图所示 , 假设半径R为1 , O是球心 , A (lat1,lon1) 和 B (lat2,lon2) 是我们感兴趣的2个点 。2跟经度线 lon1 , lon2相交于北极(north pole)N 。EF所在的线是赤道() 。ACBD是平面上的等腰梯形的四个顶点() 。AC和DB的弦(直线)在图上没有画出 。CD的位置是:C (lat2,lon1) and D (lat1,lon2) 。角度AOC是A点与C点的纬度差 dlat 。角度EOF是经度E点和经度F点的差dlon 。
(图2)
弦AC的长度 , 参照图1的方式 , 那么是AC=2*sin(dlat/2) , 弦BD也是一样的长度 。