这总是很有趣:)

首先：Mohsen Nosratinia的回答是好的，只要

您不需要知道实际距离

你可以绝对肯定地保证你永远不会去极地附近

并且永远不会接近±180°子午线

对于给定的纬度，-180°和180°经度实际上是相同的点，因此仅仅观察角度之间的差异是不够的 . 这在极地地区将是一个更大的问题，因为那里的大的经度差异对实际距离的影响较小 .

球面坐标对于导航，绘图和类似的东西非常有用和实用 . 然而，对于空间计算，就像您尝试计算的表面距离一样，球面坐标实际上非常麻烦 .

虽然可以直接使用角度进行这样的计算，但我个人并不认为它非常实用：你经常需要具备球形三角学的强大背景，并且有很多经验可以知道它的许多缺陷 - 通常存在不稳定性或"special points"你需要解决(例如，极点)，因为你引入的触发函数等需要考虑的象限模糊度等 .

我已经学会了在大学里完成所有这些工作，但我也了解到球形触发方法通常会引入复杂性，从数学角度讲并不是严格要求的，换句话说，球形触发是最基本问题的最简单表示 .

例如，如果将纬度和经度转换为3D笛卡尔 X,Y,Z 坐标，然后通过简单公式找到距离，则距离问题非常微不足道

距离(a，b)= R·arccos(a / | a |·b / | b |)

其中 a 和 b 是球体上的两个这样的笛卡尔向量 . 请注意| a | = | b | = R ，地球半径为 R = 6371 .

在MATLAB代码中：

% Some example coordinates (degrees are assumed)

lon = 360*rand(2030, 1354);

lat = 180*rand(2030, 1354) - 90;

% Your point of interest

P = [4, 54];

% Radius of Earth

RE = 6371;

% Convert the array of lat/lon coordinates to Cartesian vectors

% NOTE: sph2cart expects radians

% NOTE: use radius 1, so we don't have to normalize the vectors

[X,Y,Z] = sph2cart( lon*pi/180, lat*pi/180, 1);

% Same for your point of interest

[xP,yP,zP] = sph2cart(P(1)*pi/180, P(2)*pi/180, 1);

% The minimum distance, and the linear index where that distance was found

% NOTE: force the dot product into the interval [-1 +1]. This prevents

% slight overshoots due to numerical artifacts

dotProd = xP*X(:) + yP*Y(:) + zP*Z(:);

[minDist, index] = min( RE*acos( min(max(-1,dotProd),1) ) );

% Convert that linear index to 2D subscripts

[ii,jj] = ind2sub(size(lon), index)

如果你坚持将转换转换为笛卡尔并直接使用lat / lon，则必须使用Haversine公式，例如on this website，这也是 distance() 从映射工具箱中使用的方法 .

现在，所有这一切对整个地球都有效，你发现光滑的球形地球足够精确近似 . 如果你想要包括地球的扁率或一些更高阶的形状模型(或上帝禁止，包括地形在内的距离)，你需要做更复杂的事情 . 但我认为这不是你的目标:)

PS - 我不会感到惊讶，如果你能写出我所做的一切，你可能会重新发现Haversine公式 . 我只是喜欢能够计算一些简单的东西，就像单独从第一原理开始沿着球体的距离，而不是从很久以前植入脑中的一些黑盒子公式:)