matlab - Lat/Lng 点到 Minor Arc 段的距离

Question

我需要计算从 lat/lng GPS 点 P 到由另外 2 个 lat/lng GPS 点 A 和 B 描述的线段的最短距离。

'交叉航迹距离'帮助我计算 P 与 A 和 B 描述的大圆之间的最短距离。

然而，这不是我想要的。我需要 P 和 A-B 的 segment 线之间的距离，而不是整个大圆。

我使用了来自 http://www.movable-type.co.uk/scripts/latlong.html 的以下实现

Formula:    dxt = asin( sin(δ13) ⋅ sin(θ13−θ12) ) ⋅ R
where:
δ13 is (angular) distance from start point to third point
θ13 is (initial) bearing from start point to third point
θ12 is (initial) bearing from start point to end point
R is the earth’s radius

以下图片有望展示我正在尝试解决的问题:

在第一个图像中，由绿色线表示的交叉轨道距离是正确的，并且确实是到线段 AB 的最短距离。

在第二张图片中显示了交叉轨道距离的问题，在这种情况下，我希望最短距离是简单距离 AP，但是交叉轨道距离给出了红色<指示的距离/strong> 行。

如何更改我的算法以考虑到这一点，或检查点 X 是否在 AB 内。有可能通过计算来做到这一点吗？还是迭代是唯一可能的(昂贵的)解决方案？ (沿AB取N个点，计算P点到所有这些点的最小距离)

为了简单起见，图像中的所有线条都是直的。实际上，这些是大圆上的小弧

Answer 1

首先，一些术语:
我们的弧线是从 p1 到 p2 绘制的。
我们的第三点是p3。
与大圆相交的虚点是 p4。
p1 由 lat1,lon1 定义； p2 通过 lat2,lon2；等
dis12 是 p1 到 p2 的距离；等
bear12 是从 p1 到 p2 的方位角；等
dxt 是交叉轨道距离。
dxa是跨弧距离，我们的目标!

请注意，交叉轨道公式依赖于相对方位，bear13-bear12

我们有 3 个案例要处理。

案例1:相对方位角为钝角。所以，dxa=dis13。

案例 2.1:相对方位是锐角，并且 p4 落在我们的弧上。所以，dxa=dxt。

案例 2.2: 相对方位是锐角的，并且 p4 超出了我们的弧度。所以，dxa=dis23

算法:

第一步:如果相对方位角是钝角，dxa=dis13
完毕!
第 2 步:如果相对方位是锐角:
2.1:寻找dxt。
2.3:寻找dis12。
2.4:查找dis14。
2.4:如果 dis14>dis12，dxa=dis23。
完毕!
2.5:如果我们到达这里，dxa=abs(dxt)

MATLAB代码:

function [ dxa ] = crossarc( lat1,lon1,lat2,lon2,lat3,lon3 )
%// CROSSARC Calculates the shortest distance in meters 
%// between an arc (defined by p1 and p2) and a third point, p3.
%// Input lat1,lon1,lat2,lon2,lat3,lon3 in degrees.
    lat1=deg2rad(lat1); lat2=deg2rad(lat2); lat3=deg2rad(lat3);
    lon1=deg2rad(lon1); lon2=deg2rad(lon2); lon3=deg2rad(lon3);

    R=6371000; %// Earth's radius in meters
    %// Prerequisites for the formulas
    bear12 = bear(lat1,lon1,lat2,lon2);
    bear13 = bear(lat1,lon1,lat3,lon3);
    dis13 = dis(lat1,lon1,lat3,lon3);

    diff = abs(bear13-bear12);
    if diff > pi
        diff = 2 * pi - diff;
    end
    %// Is relative bearing obtuse?
    if diff>(pi/2)
        dxa=dis13;
    else
        %// Find the cross-track distance.
        dxt = asin( sin(dis13/R)* sin(bear13 - bear12) ) * R;

        %// Is p4 beyond the arc?
        dis12 = dis(lat1,lon1,lat2,lon2);
        dis14 = acos( cos(dis13/R) / cos(dxt/R) ) * R;
        if dis14>dis12
            dxa=dis(lat2,lon2,lat3,lon3);
        else
            dxa=abs(dxt);
        end   
    end
end

function [ d ] = dis( latA, lonA, latB, lonB )
%DIS Finds the distance between two lat/lon points.
R=6371000;
d = acos( sin(latA)*sin(latB) + cos(latA)*cos(latB)*cos(lonB-lonA) ) * R;
end

function [ b ] = bear( latA,lonA,latB,lonB )
%BEAR Finds the bearing from one lat/lon point to another.
b=atan2( sin(lonB-lonA)*cos(latB) , ...
    cos(latA)*sin(latB) - sin(latA)*cos(latB)*cos(lonB-lonA) );
end

示例输出:演示所有案例。请参阅下面的 map 。

>> crossarc(-10.1,-55.5,-15.2,-45.1,-10.5,-62.5)
ans =
   7.6709e+05
>> crossarc(40.5,60.5,50.5,80.5,51,69)
ans =
   4.7961e+05
>> crossarc(21.72,35.61,23.65,40.7,25,42)
ans =
   1.9971e+05

map 上那些相同的输出!

演示案例 1:

演示案例 2.1:

演示案例 2.2:

归功于:http://www.movable-type.co.uk/scripts/latlong.html
对于公式
和:http://www.darrinward.com/lat-long/?id=1788764
用于生成 map 图像。