从Lat/Lng点到Minor Arc段的距离

Chr*_*ker 16 gis matlab gps distance latitude-longitude

我需要计算从lat/lng GPS点P到2个其他lat/lng GPS点A和B描述的线段的最短距离.

"跨轨道距离"帮助我计算P和A和B描述的大圆之间的最短距离.

但是,这不是我想要的.我需要P和AB 线段之间的距离,而不是整个大圆.

我使用了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
Run Code Online (Sandbox Code Playgroud)

以下图片有望展示我想要解决的问题: 跨轨道距离正确 跨轨距离不正确

在第一幅图像中,由绿线表示的跨轨道距离是正确的,并且实际上是到线段AB的最短距离.

在第二个图像中显示了跨轨道距离的问题,在这种情况下,我希望最短距离是简单距离AP,但是跨轨道距离给出了由红线指示的距离.

如何更改算法以考虑这一点,或检查点X是否在AB内.是否可以通过计算方式完成此操作?或者迭代是唯一可能的(昂贵的)解决方案?(沿着AB取N点并计算从P到所有这些点的最小距离)

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

小智 19

首先,一些命名法:
我们的弧线从p1到p2.
我们的第三点是p3.
与大圆相交的虚点是p4.
p1由lat1,lon1定义; p2 by lat2,lon2; 等
dis12是从P1到P2的距离; 等
bear12是从P1到P2的轴承; 等
dxt是跨轨道距离.
dxa是跨弧距离,我们的目标!

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

我们有3个案件需要处理.

案例1:相对轴承是钝的.所以,dxa = dis13.

情况1

案例2.1:相对轴承是锐角,并且p4落在我们的弧上.所以,dxa = dxt.

案例2.1

案例2.2:相对轴承是锐角,并且p4超出了我们的弧度.所以,dxa = dis23

在此输入图像描述

算法:

第1步:如果相对方位是钝的,则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);

    %// Is relative bearing obtuse?
    if abs(bear13-bear12)>(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
Run Code Online (Sandbox Code Playgroud)

样本输出:展示所有案例.见下面的地图.

>> 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
Run Code Online (Sandbox Code Playgroud)

地图上的那些相同的输出!:

展示案例1:

案例1在地图上

展示案例2.1:

案例2.1在地图上

展示案例2.2:

案例2.2在地图上

感谢:http ://www.movable-type.co.uk/scripts/latlong.html
公式
和:http ://www.darrinward.com/lat-long/?id = 1788764
用于生成地图图像.

  • 谢谢!如果我没有记错当前代码有轻微的bug,这个条件太弱了`if abs(bear13-bear12)>(pi/2)`,因为可以有镜像的情况。因此,首先最好计算差异——`diff = abs(bear13-bear12)`。然后检查差异是否大于 180 度——`if diff>PI then diff = 2*PI-diff`(这不是标准化,而是镜像翻转)。其余的遵循`if diff>pi/2...` (2认同)