与直角参考线的距离

Question

几周前，我收到了一个脚本，我认为该脚本计算了一个点与参考线（直角）的距离，并返回沿参考线的距离以及从参考线到该点的距离观点。老实说，我认为它是这样做的，但是我并不完全理解这个脚本，即使它只有 26 行长。

由于我认为脚本是另一个依赖它的脚本变得不稳定的问题，我认为真正了解正在发生的事情实际上很方便。我曾希望有人能帮我解决这个问题。

脚本（用 Matlab 编写，但如果需要，我也可以在 python 中获得）：

function [sp,np]=locate(s,x,y,xp,yp)
ipp=0;
ns=length(s);
for ip=1:length(xp);
    for i=1:ns-2
        cosa=(x(i+2)-x(i))/(s(i+2)-s(i));
        sina=(y(i+2)-y(i))/(s(i+2)-s(i));
        sproj=s(i)+(xp(ip)-x(i))*cosa+(yp(ip)-y(i))*sina;
        if sproj<s(1)
            ipp=ipp+1;
            sp(ipp)=sproj;
            np(ipp)=-(xp(ip)-x(1))*sina+(yp(ip)-y(1))*cosa;
            break
        elseif sproj>=s(i)&sproj<=s(i+2)
            ipp=ipp+1;
            sp(ipp)=sproj;
            np(ipp)=-(xp(ip)-x(i))*sina+(yp(ip)-y(i))*cosa;
            break
        elseif sproj>s(ns)
            ipp=ipp+1;
            sp(ipp)=sproj;
            np(ipp)=-(xp(ip)-x(ns))*sina+(yp(ip)-y(ns))*cosa;
            break
        end
    end
end

这里 s 是沿参考线的距离，x 和 y 是参考线上的 x 和 y 点，xp 和 yp 也是与线的距离（直角）的点的坐标因为需要计算沿参考线的距离。脚本调用如下：

dist(1)=0;
for i=2:length(xref);
    dist(i)=dist(i-1)+sqrt((xref(i)-xref(i-1))^2+(yref(i)-yref(i-1))^2);
end

%% Create computational grid
ds=(dist(end)-dist(1))/(ns-1);  % stepsize
s=0:ds:dist(end);               % distance
xr=spline(dist,xref,s);            % x of line grid points
yr=spline(dist,yref,s);            % y of line grid points

%% Compute locations of initial line
[si,ni]=locate(s,xr,yr,xi,yi);
n=interp1(si,ni,s,'linear','extrap');    % distance change at right angle

有没有人可以帮助我解释第一个脚本中究竟发生了什么，因为我真的无法理解它。

Answer 1

我对第一个函数的看法取决于数组 x,y 描述的点是否位于一条线上（如你所说）（如第二个函数所暗示的那样；在第二个脚本中，它们是从样条插值获得的。

如果 x,y 是直线上的连续点

那么循环内的sina、cosa的计算是多余的，因为角度总是相同的。对于数组 xp,yp 中的每个点，该函数计算其在直线上的正交投影。它返回沿线的位置，从第一个点 x(1),y(1) 开始测量，以及有符号的法线距离。对于这项任务，所有这些都以一种相当低效的方式发生。

如果 x,y 不在一条直线上

那么sina、cosa的计算是不正确的。它们是从 (x(i),y(i)) 到 (x(i+2),y(i+2)) 方向的单位向量的分量。但是分母是s(i+2)-s(i)，[看第二个函数]一般会大于(x(i),y(i))到(x( i+2),y(i+2))。

以下会更有意义：

cosa = (x(i+2)-x(i))/sqrt((x(i+2)-x(i))^2 + (y(i+2)-y(i))^2);
sina = (y(i+2)-y(i))/sqrt((x(i+2)-x(i))^2 + (y(i+2)-y(i))^2);

直线距离和沿曲线距离的混淆一直存在于函数的其余部分：sproj 与 s(i+2) 的比较，尤其是与 s(ns) 的比较，并非出于几何动机。

无论如何，该函数试图定位 xp,yp 在曲线上的正交投影附近的一个点（这不是唯一的开始）；根据曲线的摆动程度，它可能会成功或失败。

Answer 2

如果我明白你在寻找什么，那就是它。

这个函数接收三个点，前两个作为 x 值向量和 y 值向量（我将它们称为 p1 和 p2）和第三个点 (p3)。该函数返回 p3 到与 p1&p2 相交的直线的最小距离以及从 (p1 到交点) 或 (p2 到交点)的距离的最小值。

例子：

定义线的两个点是 (1,2) 和 (5,10)。定义它们的线是 y = 2x。如果第三点是 (1,0)。垂线是 y = -(x-1)/2。两条线的交点将是 (2/5, 1/5)。 定义参考线的最近点是 (1,2)，从 (2/5, 1/5) 到 (1,2) 的距离是 1.7，从 (2/5, 1/5) 到 ( 1,0) 为 0.89。该函数将返回 [1.7, 0.89]。

如果这是你想要的，请告诉我。

function [ dRef, dInt ] = StackExchange( x, y, xP, yP )
%StackExchange: Calculate line intersecting (x(1), y(1)) & (x(2), y(2))
%finds perpendicular line which intersects (xP, yP) which occurs as
%(IP(1), IP(2)) then calculated the minimum distance from the reference points
%(x(1), y(1)) or (x(2), y(2)) to (IP(1), IP(2))
%and from (IP(1), IP(2)) to (xP, yP)

%Calculate  slope of reference line, the slope of the perpendicular line
%will be -1/m
m = diff(y) ./ diff(x) ;

%IP(1) is the x-coordinate of the intersection point
%This formula is solving the equation
% m * (x - x(1)) + y(1) = -(1/m) * (x - xP) + yP
%for x
IP(1) = (m^2 * x(1) + xP + m * (yP - y(1))) ./ (m^2 + 1);
%IP(2) is the y-coordinate of the intersection point
IP(2) = m * (IP(1) - x(1)) + y(1);

%Minimum distance from the reference points to the intersection point
dRef = sqrt(min((IP(1) - x).^2 + (IP(2) - y).^2));

%Distance from the intersection point to  point of interest
dInt = sqrt(sum(([xP, yP] - IP).^2));

end