之前遇到过这样的问题:手里有一批患者的具体居住地点以及对应的经纬度,我想分析下患者的具体分布情况,不同患者之间的距离关系,这时就需要根据某两点的经纬度来求该两点的实际距离了。下面就是计算公式和代码:
公式
球面上任意两点的距离计算公式可以参考维基百科上的下述文章。
Great-circle distance
Haversine formula
下面采用的是维基百科推荐的Haversine公式,原因:
- Great-circle distance公式用到了大量余弦函数,而两点间距离很短时(比如地球表面上相距几百米的两点),余弦函数会得出0.999…的结果,会导致较大的舍入误差。
- 而Haversine公式采用了正弦函数,即使距离很小,也能保持足够的有效数字。
其中:
- R为地球半径,可取平均值 6371km;
- φ1, φ2 表示两点的纬度;
- Δλ 表示两点经度的差值。
程序代码:
# 输入
# 计算地图上两点经纬度间的距离
from math import radians, cos, sin, asin, sqrt
# Haversine(lon1, lat1, lon2, lat2)的参数代表:经度1,纬度1,经度2,纬度2(十进制度数)
def Haversine(lon1, lat1, lon2, lat2):
# 将十进制度数转化为弧度
lon1, lat1, lon2, lat2 = map(radians, [lon1, lat1, lon2, lat2])
# Haversine公式
dlon = lon2 - lon1
dlat = lat2 - lat1
a = sin(dlat/2)**2 + cos(lat1) * cos(lat2) * sin(dlon/2)**2
c = 2 * asin(sqrt(a))
r = 6371 # 地球平均半径,单位为公里
d = c * r
print("该两点间距离={0:0.3f} km".format(d))
# 广州市人民政府 113.270714,23.13552
# 深圳市人民政府 114.064803,22.549054
Haversine(113.270714,23.13552,114.064803,22.549054)
# 输出
该两点间距离=104.280 km