（译）计算距离方位和更多经纬度之间的点

Posted 2021-02-12 surimn

tags:

篇首语：本文由小常识网(cha138.com)小编为大家整理，主要介绍了（译）计算距离方位和更多经纬度之间的点相关的知识，希望对你有一定的参考价值。

计算距离、方位和更多经纬度之间的点。最近在研究预测未来坐标和速度、时间之间的关系，希望这篇文章对地图应用有所帮助。

作者：狐狸家的鱼

本文链接：计算距离、方位和更多经纬度之间的点

原文链接：Calculate distance, bearing and more between Latitude/Longitude points

GitHub：sueRimn

该页面介绍了对经纬度点的各种计算，以及实现它们的公式和代码片段。

所有这些公式都是基于球形地球的计算（忽略椭球效应） - 大多数情况下都是准确的（事实上，地球是非常略微椭圆形的，使用球形模型的误差通常高达0.3%¹ - 有关详细信息，请参阅注释）。

1、距离

这使用‘ hasrsine ‘公式¹来计算两点之间的大圆弧距离 - 即地球表面上的最短距离 - 给出点之间的最短距离。

hasrsine: a = sin²(Δφ/2) + cos φ₁ ⋅ cos φ₂ ⋅ sin²(Δλ/2)

公式: c = 2 ⋅ atan2( √a, √(1−a) )

d = R ⋅ c

其中φ为纬度，λ为经度，R为地球半径(平均半径= 6,371km)。注意角度必须是弧度才能传递给三角函数!

javascript:    
var R = 6371e3; // metres
var φ1 = lat1.toRadians();
var φ2 = lat2.toRadians();
var Δφ = (lat2-lat1).toRadians();
var Δλ = (lon2-lon1).toRadians();

var a = Math.sin(Δφ/2) * Math.sin(Δφ/2) +
        Math.cos(φ1) * Math.cos(φ2) *
        Math.sin(Δλ/2) * Math.sin(Δλ/2);
var c = 2 * Math.atan2(Math.sqrt(a), Math.sqrt(1-a));

var d = R * c;

在这些脚本中，我通常用 lat/lon 表示度，以 φ/λ 表示的距离/长度。发现混合度和半径是经常让人抓耳挠腮的问题。

即使在很小的距离条件下，“haversine”公式也可以很好的进行不同于基于余弦球面定律的数值计算。

以上使用的"(重新)正矢"是 1−cosθ ，以及 “半正矢”是 (1−cosθ)/2 或者 sin²(θ/2)。曾被航海家广泛使用，罗杰·辛诺特于1984年在《天空与望远镜》杂志上对其进行了描述：辛诺特解释说，大熊星座中Mizar（大熊星座中的一个恒星）和Alcor（北斗七星中的第五等星）之间的角间距(0°11’49.69”)可以通过使用haversine在TRS-80上精确计算出来。

奇怪的是，c是弧度的角距离，a是两点弦长一半的平方。

如果 atan2 不可用, c 可以从 2 ⋅ asin( min(1, √a) ) 计算得出(包括防止四舍五入误差)。

在中等核的i5个人电脑上使用Chrome，距离计算大约需要2 - 5微秒（因此大约每秒200,000 - 500,000次）。通过公因子分解几乎没有用，因为JIT编译器可能会对它们进行优化。

2、余弦球面定律

事实上，JavaScript(以及大多数现代计算机和语言)使用提供了15个重要的精度数字的 “IEEE 754” 64位浮点数。据我估计，在简单余弦公式球面定律(cos c = cos a cos b + sin a sin b cos C)的精度下给出了良好的结果，距离小到地球表面几米。(注意，余弦定律的大地测量形式是规范重新排列的，这样就可以直接使用纬度，而不是用余纬)。

这使得简单余弦定理可以在许多大地测量目的（如果不是天文学）中作为 haversine 公式的另一个选择。这种选择可能是由编程语言、处理器、编码上下文、可用地三角函数（在不同语言中）等驱动的，并且，在非常小地距离中。等矩形的近似值可能更合适。

余弦定理：d = acos( sin φ₁ ⋅ sin φ₂ + cos φ₁ ⋅ cos φ₂ ⋅ cos Δλ ) ⋅ R

JavaScript:    
var φ1 = lat1.toRadians(), φ2 = lat2.toRadians(), Δλ = (lon2-lon1).toRadians(), R = 6371e3; // gives d in metres
var d = Math.acos( Math.sin(φ1)*Math.sin(φ2) + Math.cos(φ1)*Math.cos(φ2) * Math.cos(Δλ) ) * R;

Excel:  
=ACOS( SIN(lat1)*SIN(lat2) + COS(lat1)*COS(lat2)*COS(lon2-lon1) ) * 6371000

(or with lat/lon in degrees): 
=ACOS( SIN(lat1*PI()/180)*SIN(lat2*PI()/180) + COS(lat1*PI()/180)*COS(lat2*PI()/180)*COS(lon2*PI()/180-lon1*PI()/180) ) * 6371000

虽然比较简单,但在我的测试中,余弦定理比 haversine 略慢。

3、等矩形近似值

如果性能是一个问题,准确性就不那么重要了。在很小距离上，勾股定理可以在这个等矩形投影中使用：

公式： x = Δλ ⋅ cos φm

y = Δφ

d = R ⋅ √x² + y²

JavaScript:    
var x = (λ2-λ1) * Math.cos((φ1+φ2)/2);
var y = (φ2-φ1);
var d = Math.sqrt(x*x + y*y) * R;

这仅仅使用一个三角函数和一个根号函数，而cos定律中有6个三角函数，而 haversine 中有7个三角函数+ 2根号函数。精度有点复杂:沿着子午线没有误差，否则它们依赖于距离、方位和纬度，但对于许多用途来说足够小(通常与球面近似本身相比微不足道)。

此外，还可以使用极坐标平面地球公式：使用余纬 θ1 = π/2−φ1 and θ2 = π/2−φ2，然后 θ1 = π/2−φ1 and θ2 = π/2−φ2。我没有比较精度。

4、方位

一般来说，当你沿着一个大圆弧路径(直角方向)移动时，当前朝向会发生变化。根据距离和纬度的不同，最终的航向与最初的航向会有不同程度的差异（比如从35°N 45°E(≈Baghdad)到35°N 135°E(≈Osaka),你会从一个60°的航向开始，最终到达一个120°的航向)。

这个公式适用于如果初始方位(有时称为前方位角)沿着大圆弧走一条直线，你会从起点走到终点^1：

公式： θ = atan2( sin Δλ ⋅ cos φ₂ , cos φ₁ ⋅ sin φ₂ − sin φ₁ ⋅ cos φ₂ ⋅ cos Δλ )

φ₁,λ₁是起点,φ₂,λ₂是终点(Δλ是经度的差异)

JavaScript:
(all angles in radians)
var y = Math.sin(λ2-λ1) * Math.cos(φ2);
var x = Math.cos(φ1)*Math.sin(φ2) -
        Math.sin(φ1)*Math.cos(φ2)*Math.cos(λ2-λ1);
var brng = Math.atan2(y, x).toDegrees();

Excel:
(all angles in radians)
=ATAN2(COS(lat1)*SIN(lat2)-SIN(lat1)*COS(lat2)*COS(lon2-lon1), 
       SIN(lon2-lon1)*COS(lat2))

注意，Excel将参数转换为ATAN2—请参阅下面的注释

由于 atan2 返回值范围中的值-π……+π(-180°…+180°)，将结果正常化为罗盘方位(范围为0°…360°，- ve值转换为180°范围…360°),转换为度,然后使用(θ+ 360)% 360,%是(浮点数)模。

对最后的方位,只是把初始方位从终点到起点并反转它(使用θ=θ+ 180)% 360)。

5、中点

这是两点之间沿大圆路径的中点¹。

公式： B_x = cos φ₂ ⋅ cos Δλ

B_y = cos φ₂ ⋅ sin Δλ

φ_m = atan2( sin φ₁ + sin φ₂, √(cos φ₁ + B_x)² + B_y² )
λ_m = λ₁ + atan2(B_y, cos(φ₁)+B_x)

JavaScript:
(all angles in radians)
var Bx = Math.cos(φ2) * Math.cos(λ2-λ1);
var By = Math.cos(φ2) * Math.sin(λ2-λ1);
var φ3 = Math.atan2(Math.sin(φ1) + Math.sin(φ2),
              Math.sqrt( (Math.cos(φ1)+Bx)*(Math.cos(φ1)+Bx) +By*By ) );
var λ3 = λ1 + Math.atan2(By, Math.cos(φ1) + Bx);

经度可以使用(lon+540)%360-180归一化为-180…+180

正如初始方位可能与最终方位不同，中点也可能不在经纬度之间。35°N、45°E和35°N、135°E之间的中点在45°N、90°E附近。

6、中间点

也可以计算出两点之间沿大圆弧路径上任意处的中间点¹。

公式：a = sin((1−f)⋅δ) / sin δ

b = sin(f⋅δ) / sin δ

x = a ⋅ cos φ₁ ⋅ cos λ₁ + b ⋅ cos φ₂ ⋅ cos λ₂

y = a ⋅ cos φ₁ ⋅ sin λ₁ + b ⋅ cos φ₂ ⋅ sin λ₂

z = a ⋅ sin φ₁ + b ⋅ sin φ₂

φ_i = atan2(z, √x² + y²)

λ_i = atan2(y, x)

f是沿大圆航线的一部分(f = 0是点1,f = 1是点2),δ是两点之间的角距离d / R。

7、目标点给出从起点到终点的距离和方位方位

已知一个起始点、初始方位和距离，这将计算沿着(最短距离)大圆弧运动的终点和最终方位。

公式：φ₂ = asin( sin φ₁ ⋅ cos δ + cos φ₁ ⋅ sin δ ⋅ cos θ )

λ₂ = λ₁ + atan2( sin θ ⋅ sin δ ⋅ cos φ₁, cos δ − sin φ₁ ⋅ sin φ₂ )

φ是纬度，λ是经度，θ是轴承方位(顺时针从北),δ是角距离d / R，d是经过的距离，R是地球的半径

JavaScript:
(all angles in radians)
var φ2 = Math.asin( Math.sin(φ1)*Math.cos(d/R) +
                    Math.cos(φ1)*Math.sin(d/R)*Math.cos(brng) );
var λ2 = λ1 +Math.atan2(Math.sin(brng)*Math.sin(d/R)*Math.cos(φ1),
                         Math.cos(d/R)-Math.sin(φ1)*Math.sin(φ2));

经度可以使用 (lon+540)%360-180 规范化为 -180…+180

Excel:
(all angles in radians)
lat2: =ASIN(SIN(lat1)*COS(d/R) + COS(lat1)*SIN(d/R)*COS(brng))
lon2: =lon1 + ATAN2(COS(d/R)-SIN(lat1)*SIN(lat2), SIN(brng)*SIN(d/R)*COS(lat1)) * Remember that Excel reverses the arguments to ATAN2 – see notes below

对于最终方位，只需将初始方位从起点移到起点，然后用(brng+180)%360将其反转。

9、给定起点和方位两条路径的交点

这是一个比这一页上的大多数计算要复杂得多的计算，但是我已经被问过很多次了。这来自埃德·威廉的航空公式。参见下面的JavaScript。

公式： δ₁₂ = 2⋅asin( √(sin²(Δφ/2) + cos φ₁ ⋅ cos φ₂ ⋅ sin²(Δλ/2)) ) （角距离. p1–p2）

θ_a = acos( ( sin φ₂ − sin φ₁ ⋅ cos δ₁₂ ) / ( sin δ₁₂ ⋅ cos φ₁ ) ) （初始/ 最终方位）

θ_a = acos( ( sin φ₂ − sin φ₁ ⋅ cos δ₁₂ ) / ( sin δ₁₂ ⋅ cos φ₁ ) ) （点1和点2之间）

    if sin(λ2−λ1) > 0
       θ12 = θa
       θ21 = 2π − θb
    else
       θ12 = 2π − θa
       θ21 = θb    
       α1 = θ13 − θ12          （角度 p2–p1–p3）
       α2 = θ21 − θ23          （角度 p1–p2–p3）

　　　α₃ = acos( −cos α₁ ⋅ cos α₂ + sin α₁ ⋅ sin α₂ ⋅ cos δ₁₂ ) （角度 p1–p2–p3）

　　 δ₁₃ = atan2( sin δ₁₂ ⋅ sin α₁ ⋅ sin α₂ , cos α₂ + cos α₁ ⋅ cos α₃ ) （角距离 p1–p3）

φ₃ = asin( sin φ₁ ⋅ cos δ₁₃ + cos φ₁ ⋅ sin δ₁₃ ⋅ cos θ₁₃) （p3纬度）

Δλ₁₃ = atan2( sin θ₁3 ⋅ sin δ₁₃ ⋅ cos φ₁ , cos δ₁₃ − sin φ₁ ⋅ sin φ₃ ) （经度p1-p3）

λ₃ = λ₁ + Δλ₁₃ (p3 经度)

φ₁、λ₁、θ₁₃：1起点&(初始)轴承从1日指向交点

φ₂、λ₂、θ₂₃：2日起点&(初始)轴承2指向交点

φ₃、λ₃：交点

% =(浮点数)模

如果sin α₁= 0 以及 sin α₂ = 0：无穷多解，如果sin α₁ ⋅ sin α₂ < 0：歧义解。这个公式并不总是适用于经向线或赤道线。

用向量比用球面三角简单多了： latlong-vectors.html.

10、航迹距离

这是一个新的问题：有时我被问到关于一个点到大圆弧路径的距离(有时称为交叉轨迹误差)。

　　　　　　公式： d_xt = asin( sin(δ₁₃) ⋅ sin(θ₁₃−θ₁₂) ) ⋅ R

　　　　δ_13是(角)距离起点第三点

　　　　　　　θ₁3是(初始)轴承从起点到第三点

θ₁₂是(初始)轴承从起点到终点

R是地球的半径

JavaScript:    
var δ13 = d13 / R;
var dAt = Math.acos(Math.cos(δ13)/Math.cos(dXt/R)) * R;

11、离两极最近的点

“克莱劳公式”会用给出的方位 θ 和最大圆上的纬度 φ，计算出一个大圆弧路径的最大纬度：

公式：φ_max = acos( | sin θ ⋅ cos φ | )

JavaScript:    
var φMax = Math.acos(Math.abs(Math.sin(θ)*Math.cos(φ)));

12、恒向线

“恒向线”(或斜航线)是一条恒定方位的路径，它以相同的角度穿过所有子午线。

水手们过去常常(现在有时仍然)沿着恒向线航行，因为跟随着一个恒定的罗盘方位比跟随绕一个大圆需要不断调整方位要容易得多。等高线是墨卡托投影图上的直线(也有助于导航)。

等高线一般比大圆线长。例如，从伦敦到纽约，沿一条斜脊线比沿一个大圆要长4%——这对航空燃料很重要，但对帆船来说并不特别重要。从纽约到北京——接近可能的最极端的例子(虽然不能航行)——沿着一条斜脊线要长30%。

计算恒向线的关键是范德曼函数，它给出了已知纬度的墨卡托投影地图上的高度：ln(tanφ + secφ) 或是 ln( tan(π/4+φ/2) )。当然，这在两极是趋于无穷的(与墨卡托投影一致)。对于强迫症患者,甚至有一个椭圆形的版本,“等量纬度”：ψ= ln (tan(π/ 4 +φ/ 2)/ ((1−e⋅sinφ)/ (1 + e⋅sinφ)]^{e / 2})，或其更好的等值于 ψ= atanh (sinφ)−e⋅atanh (e⋅sinφ)。

该公式从球形纬度和经度坐标推算出墨卡托投影以东和以北值，然后：

E = R ⋅ λ

N = R ⋅ ln( tan(π/4 + φ/2) )

以下公式来源于埃德·威廉姆斯的航空公式¹。

1.距离

由于恒向线是墨卡托投影上的直线，沿恒向线的两点之间的距离就是该线的长度(毕达哥拉斯)，但是项目的失真需要得到补偿。

按固定纬度航行(东西航行)，cosφ 补偿很简单，通常情况下，它是Δφ/Δψ，Δψ = ln( tan(π/4 + φ2/2) / tan(π/4 + φ1/2) ) （“投影”纬度差）

公式：Δψ = ln( tan(π/4 + φ₂/2) / tan(π/4 + φ₁/2) ) （“投影”纬度差）

q = Δφ/Δψ (or cosφ for E-W line)

d = √(Δφ² + q²⋅Δλ²) ⋅ R

φ是纬度,λ是经度,Δλ是正在走的最短路线(< 180°),R是地球的半径,ln是自然对数

JavaScript:
(all angles in radians)
var Δψ =Math.log(Math.tan(Math.PI/4+φ2/2)/Math.tan(Math.PI/4+φ1/2));
var q = Math.abs(Δψ) > 10e-12 ? Δφ/Δψ : Math.cos(φ1); // E-W course becomes ill-conditioned with 0/0

// if dLon over 180° take shorter rhumb line across the anti-meridian:
if (Math.abs(Δλ) > Math.PI) Δλ = Δλ>0 ? -(2*Math.PI-Δλ) : (2*Math.PI+Δλ);

var dist = Math.sqrt(Δφ*Δφ + q*q*Δλ*Δλ) * R;

2.方位

恒向线是墨卡托投影上的一条直线，其投影角度与罗盘方位相等。

公式：Δψ = ln( tan(π/4 + φ2/2) / tan(π/4 + φ1/2) ) （“投影”纬度差）

θ = atan2(Δλ, Δψ)

φ是纬度,λ是经度,Δλ是正在走的最短路线(< 180°),R是地球的半径,ln是自然对数

JavaScript:
(all angles in radians)
var Δψ =Math.log(Math.tan(Math.PI/4+φ2/2)/Math.tan(Math.PI/4+φ1/2));

// if dLon over 180° take shorter rhumb line across the anti-meridian:
if (Math.abs(Δλ) > Math.PI) Δλ = Δλ>0 ? -(2*Math.PI-Δλ) : (2*Math.PI+Δλ);

var brng = Math.atan2(Δλ, Δψ).toDegrees();

3.终点

给定一个起点和一个距离d和不变的方位θ,这将计算出终点。如果你沿着一条等高线保持恒定的方位，你就会逐渐地向其中一个极点螺旋靠近。

公式：δ = d/R （角距离）

φ₂ = φ₁ + δ ⋅ cos θ

　　　　　　　　　　　　　　　　　　　　 Δψ = ln( tan(π/4 + φ₂/2) / tan(π/4 + φ₁/2) ) （“投影”纬度差）

　　　　　　　　 q = Δφ/Δψ (or cos φ for E-W line)

Δλ = δ ⋅ sin θ / q

λ2 = λ1 + Δλ

φ是纬度,λ是经度,Δλ是正在走的最短路线(< 180°),ln是自然对数,R是地球的半径

JavaScript:
(all angles in radians)
var δ = d/R;
var Δφ = δ * Math.cos(θ);
var φ2 = φ1 + Δφ;

var Δψ =Math.log(Math.tan(φ2/2+Math.PI/4)/Math.tan(φ1/2+Math.PI/4));
var q = Math.abs(Δψ) > 10e-12 ? Δφ / Δψ : Math.cos(φ1); // E-W course becomes ill-conditioned with 0/0

var Δλ = δ*Math.sin(θ)/q;
var λ2 = λ1 + Δλ;

// check for some daft bugger going past the pole, normalise latitude if so
if (Math.abs(φ2) > Math.PI/2) φ2 = φ2>0 ? Math.PI-φ2 : -Math.PI-φ2;

4.中间点

这个计算“斜向中点”的公式是由罗伯特·希尔(Robert Hill)和克莱夫·图特(Clive Tooth1) (thx Axel!)推导出来的。

公式：φm = (φ1+φ2) / 2

　　　　f1 = tan(π/4 + φ1/2)

　　　　f2 = tan(π/4 + φ2/2)

　　　　fm = tan(π/4+φm/2)

　　　　　　　　　　　　　　　　　　　　λm = [ (λ2−λ1) ⋅ ln(fm) + λ1 ⋅ ln(f2) − λ2 ⋅ ln(f1) ] / ln(f2/f1)

　　　　　　　　φ是纬度,经度λ,ln是自然对数

JavaScript:
(all angles in radians)
if (Math.abs(λ2-λ1) > Math.PI) λ1 += 2*Math.PI; // crossing anti-meridian

var φ3 = (φ1+φ2)/2;
var f1 = Math.tan(Math.PI/4 + φ1/2);
var f2 = Math.tan(Math.PI/4 + φ2/2);
var f3 = Math.tan(Math.PI/4 + φ3/2);
var λ3 = ( (λ2-λ1)*Math.log(f3) + λ1*Math.log(f2) - λ2*Math.log(f1) ) / Math.log(f2/f1);

if (!isFinite(λ3)) λ3 = (λ1+λ2)/2; // parallel of latitude

经度可以使用(lon+540)%360-180规范化为-180…+180。

5.在web页面中使用脚本

在web页面中使用这些脚本大致如下:

<!doctype html>
<html lang="en">
<head>
    <title>Using the scripts in web pages</title>
    <meta charset="utf-8">
    <script type="module">
        import LatLon from ‘https://cdn.jsdelivr.net/gh/chrisveness/[email protected]/latlon-spherical.min.js‘;
        document.addEventListener(‘DOMContentLoaded‘, function() {
            document.querySelector(‘#calc-dist‘).onclick = function() {
                calculateDistance();
            }
        });
        function calculateDistance() {
            const p1 = LatLon.parse(document.querySelector(‘#point1‘).value);
            const p2 = LatLon.parse(document.querySelector(‘#point2‘).value);
            const dist = parseFloat(p1.distanceTo(p2).toPrecision(4));
            document.querySelector(‘#result-distance‘).textContent = dist + ‘ metres‘;
        }
    </script>
</head>
<body>
<form>
    Point 1: <input type="text" name="point1" id="point1" placeholder="lat1,lon1">
    Point 2: <input type="text" name="point2" id="point2" placeholder="lat2,lon2">
    <button type="button" id="calc-dist">Calculate distance</button>
    <output id="result-distance"></output>
</form>
</body>
</html>

请参阅下面JavaScript源代码，也可以在GitHub上获得。完整的文档以及测试套件都是可用的。

6.笔记

（1）

精度：由于地球并不完全是一个球体，所以在使用球面几何时存在小的误差。地球实际上大致是椭球体(或者更精确地说，是扁圆球体)，半径在6,378千米(赤道)和6,357千米(极地)之间变化，局部曲率半径从6,336千米(赤道子午线)到6,399千米(极地)之间变化。6371千米是地球平均半径的普遍接受值。这意味着，根据纬度和旅行方向的不同，假设球面几何形状穿越赤道的误差可能高达0.55%，但一般低于0.3% (whuber在stackexchange上对此进行了详细的研究)。对于我来说，1km中精度超过3m就足够了，但是如果你想要更高的精度，你可以使用Vincenty公式来计算椭球体上的测地线距离，结果可以精确到1mm以内。(完全出于反常，我从来没有需要过这么精确的代码，我查了一下这个公式，发现JavaScript实现比我预想的要简单)。

（2）

三角函数以弧度为参数，所以纬度、经度、和方位(十进制或度/分钟/秒)需要转换为弧度，rad =π.deg / 180。当弧度转换回度时(deg = 180.rad/π),如果使用带符号的十进制度，则West为负。对于方位，其在-π到+π(-180°+ 180°)的范围需要转换为0到+ 2π(0°-360°)。这可以通过(brng+2.π)%2.π(或brng+360)%360)，%(浮点数)是模运算符(请注意,不同的语言以不同的方式实现模操作)。

（3）

所有方位均相对于正北，0°=N, 90°=E等。如果你用指南针工作，地磁北极从正北极以一种复杂的方式围绕地球发生变化，这种差异必须由本地地图上显示的差异来补偿。

（4）

这里广泛使用的atan2()函数接受两个参数，atan2(y, x)，并计算比值y/x的反正切。它比atan(y/x)更灵活,因为它处理x = 0,并且返回在所有4个象限的值-π+π(反正切函数返回在范围-π/ 2 +π/ 2之间的值)。

（5）

如果你在电子表格（Microsoft Excel, LibreOffice Calc, Google Sheets, Apple Numbers）中使用任何公式包括atan2，你将需要反转参数，如Excel等，它们与JavaScript相反，常规顺序为atan(y,x)，但是Excel使用atan2(x, y)。在（VBA）宏观中使用 atan2，你可以用 WorksheetFunction.Atan2()。

（6）

如果您正在使用谷歌地图，这些函数中的一些(computeDistanceBetween()、computeHeading()、computeOffset()、插值()等）现在在谷歌地图 API V3“球面”库中提供。注意，它们使用的默认地球半径为6,378,137米)。

（7）

如果你使用英国地形测量网格参考，我已经实现了一个脚本之间的转换Lat/Long和OS网格参考。

（8）

如果您使用UTM坐标或MGRS网格引用，我已经实现了用于在Lat/Long、UTM和MGRS之间转换的脚本。

（9）

我从美国人口普查局的GIS FAQ中学到了很多，现在已经没有了，所以我复制了一份。

（10）

感谢埃德·威廉姆斯的航空公式。

（11）

对于英里，用km除以1.609344。