最近工作需要，网上搜索了下根据经纬度计算两地距离的方法，发现要么是几何法，画图、作一堆辅助线，然后证明推理，要么二话不说直接套公式。这篇文章介绍一种容易理解的方式来求这个距离。
0b00 思路
地球是个不规则的椭球体、为了简便我们当作球体来计算。球体上两地的最短距离就是经过两点的大圆的劣弧长度。
思路如下:
弧长 ← 弦长 (两点距离) ← 两点坐标 (直角坐标) ← 经纬度
0b01 计算
1. 坐标转换
设
地球半径为 $R$
地心到 E 0° N 0° 的连线为 x 轴
地心到 E 90° N 0° 的连线为 y 轴
地心到 E 0° N 90° 的连线为 z 轴
地球表面有一点 $A$, 经度为 $e$, 纬度为 $n$, 单位为弧度
则 $A$ 的坐标可表示为:
$$
x = R cdot cos(n) cdot cos(e)\
y = R cdot cos(n) cdot sin(e)\
z = R cdot sin(n)
$$
代码
const R = 6371
const {cos, sin, PI} = Math
let getPoint = (e, n) => {
// 首先将角度转为弧度
e *= PI/180
n *= PI/180
reutrn {
x: R*cos(n)*cos(e),
y: R*cos(n)*sin(e),
z: R*sin(n)
}
}
2. 根据坐标计算两点距离
这个太简单，跳过
3. 根据弦长求弧长
这个可以画个图，帮助理解：
现在已知弦长 $c$, 半径 $R$, 要求弧 $r$ 的长度这很简单, 只需先求出 $∠alpha$ 的大小：
$$
alpha = arcsin(c/2/R)\
r = 2alpha cdot R
$$
代码
const {asin} = Math
const R = 6371
r = asin(c/2/R)*2*R
0b10 最终代码
/**
* 获取两经纬度之间的距离
* @param {number} e1 点 1 的东经, 单位: 角度, 如果是西经则为负
* @param {number} n1 点 1 的北纬, 单位: 角度, 如果是南纬则为负
* @param {number} e2
* @param {number} n2
*/
function getDistance(e1, n1, e2, n2){
const R = 6371
const {sin, cos, asin, PI, hypot} = Math
/** 根据经纬度获取点的坐标 */
let getPoint = (e, n) => {
e *= PI/180
n *= PI/180
// 这里 R* 被去掉, 相当于先求单位圆上两点的距, 最后会再将这个距离放大 R 倍
return {x: cos(n)*cos(e), y: cos(n)*sin(e), z: sin(n)}
}
let a = getPoint(e1, n1)
let b = getPoint(e2, n2)
let c = hypot(a.x – b.x, a.y – b.y, a.z – b.z)
let r = asin(c/2)*2*R
return r
}