我需要用 Matlab 解决一个最小化问题,我想知道哪个是最简单的解决方案。我一直在考虑的所有潜在解决方案都需要大量的编程工作。
假设我有一个纬度/经度坐标点 (A,B),我需要在纬度/经度坐标 map 中搜索离该点最近的点。
特别地,纬度和经度数组是两个 2030x1354 元素矩阵(1 公里距离),其想法是在这些矩阵中找到唯一索引,使到坐标 (A,B) 的距离最小,即找到最接近给定坐标 (A,B) 的值。
如有任何帮助,我们将不胜感激。
谢谢!
最佳答案
这总是一个有趣的:)
首先:Mohsen Nosratinia的回答是可以的,只要
- 你不需要知道实际距离
- 你可以绝对肯定地保证你永远不会靠近极地地区
- 永远不会靠近±180°子午线
对于给定的纬度,-180° 和 +180° 经度实际上是同一点,因此仅查看角度之间的差异是不够的。这在极地地区将是一个更大的问题,因为那里较大的经度差异对实际距离的影响较小。
球坐标对于导航、绘图等目的非常有用和实用。然而,对于空间计算,就像您要计算的表面距离一样,球坐标实际上使用起来非常麻烦。
虽然直接使用角度进行此类计算可能,但我个人认为这不是很实用:您通常需要有很强的球面背景三角学,以及了解它的许多陷阱的丰富经验——经常有你需要解决的不稳定性或“特殊点”(例如极点),由于你引入了三角函数,你需要考虑象限模糊,等等
我在大学里学会了做这一切,但我也了解到球面三角法经常引入复杂性,数学并不是严格要求的,换句话说,球面三角法是< strong>不是底层问题的最简单表示。
例如,如果您将纬度和经度转换为 3D 笛卡尔 X,Y,Z
坐标,然后通过简单的公式求出距离,那么您的距离问题就非常简单了
distance (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)
如果您坚持跳过向笛卡尔坐标系的转换并直接使用纬度/经度,则必须使用 Haversine 公式,如概述 on this website例如,这也是映射工具箱中 distance()
使用的方法。
现在,所有这些都适用于整个地球,前提是您发现光滑的球形地球足够精确的近似值。如果你想包括地球的扁率或一些更高阶的形状模型(或者上帝保佑,包括地形的距离),你需要做更复杂的事情。但我不认为这是你的目标:)
PS - 如果您将我所做的一切都写出来,我不会感到惊讶,您可能会重新发现 Haversine 公式。我只是更喜欢能够仅根据第一性原理来计算像沿着球体的距离这样简单的东西,而不是根据很久以前你植入脑海的一些黑匣子公式:)
关于matlab - 如何使用 MATLAB 找到与给定坐标最近的点?,我们在Stack Overflow上找到一个类似的问题: https://stackoverflow.com/questions/17411274/