点在折线的哪一侧

Question

我正在尝试编写一个代码，该代码将允许识别随机多段线的哪一侧（不是闭合的一侧）随机点是。折线被认为是无限的，第一段和最后一段延伸到无限。如果它们交叉折线应该被视为一个多边形（我还没有详细说明这种情况的代码）。逻辑如下：

1. 定义折线的顶点，到所考虑的点的距离最小。变量minimum2p 是到那个顶点的距离，I2p 是顶点的索引。

2. 定义多段线的线段，到所考虑的点的距离最小。只有那些与所考虑的点垂直相交的线段才算数（变量count_s）。变量minimum2s 是到段的最小距离； I2p 是该段的第一个顶点的索引； flag 是布尔变量，它存储关于哪个线段与所提到的垂线相交的信息。

3. 接下来只需选择合适的细分进行比较，例如，使用来自链接link-1、link-2 或link-3 的想法。我尝试了here 的方法，但它不适用于许多特殊情况。我在那里使用the best answer 作为折线的内部点。所以，我的方法如下：

4. 首先，检查它是折线的第一个顶点还是最后一个顶点。如果是这种情况，则选择的段相应地是第一个或最后一个，但前提是没有比第一个或最后一个更接近的其他段。如果有另一个片段，那么我选择了那个片段。

5. 接下来，如果步骤 4 不是这种情况，那么我检查折线的内部顶点。如果附近还有一个段关闭，那么我比较索引I2p 和I2s，如果最后一个存在。如果它们重合，则在选择要比较的适当段时没有歧义。如果它们不同，则优先选择最近的段而不是最近的顶点。

6. 最后，如果附近没有线段（在与线段相交的点垂直的意义上），那么对于内部顶点我应用the best answerhere的想法.

Here are some results为不同的折线，由其顶点的X和Y坐标定义，分别存储在'polylineX'和'polylineY'中（红色代表'left'位置，灰色代表'右'的位置，黑色是折线上的位置，蓝线代表折线）。

如您所见，对于相对平滑的折线，代码可以正常工作。然而，对于更清晰的，或者在某些方面复杂的，代码不能正常工作。我的代码中缺少什么？考虑到某些情况应该添加什么条件？

代码如下：

clear all
close all 
clc
clf
polylineX = [0 1 2 3 4 5 6 7 8];
polylineY = [-10 20 -13 18 -17 16 -21 23 -25];
hold on
title(['polylineX=[',num2str(polylineX),'], polylineY=[',num2str(polylineY),']'])
chosen = 0;
span = 60;

for ii = 10:70
for jj = 30:60
    ii
    jj
    position = -2;
    point = [(jj-round(span/2))/1 (ii-round(span/2))/1];

    axis equal
    plot(polylineX,polylineY,'.-','MarkerSize',1,'LineWidth',1);

    distance2p = zeros(1,length(polylineX)); % distances from the point to the points (2p) of the polyline
    distance2s = zeros(1,length(polylineX)-1); % distances from the point to the segments (2s) of the polyline
    flag = zeros(1,length(polylineX)-1);
    count_s = 0; % counter of segments, which are intersected by the normal pointing from the 'point'
    k = 0;

    for i = 1:length(polylineX)-1
        pos = sign((polylineX(i+1) - polylineX(i)) * (point(2) - polylineY(i)) -...
                   (polylineY(i+1) - polylineY(i)) * (point(1) - polylineX(i)));

        % computing the distances from the 'point' to all segments and mark if
        % the distance vectors intersect the segments
        [flag(i),distance2s(i)] = distanceToLine([polylineX(i) polylineX(i+1)],[polylineY(i) polylineY(i+1)],[point(1) point(2)]);

        if flag(i)
            if k == 0
                minimum2s = distance2s(i);
                I2s = i;
            end;
            k = 1;
            count_s = count_s + 1; % count segments, which are intersected by the normal pointing from the 'point'
            if distance2s(i) < minimum2s
                I2s = i;
                minimum2s = distance2s(i);
            end;
        end;
    end;


    % first compute the distances between the 'point' under consideration and the
    % points of the given polyline
    for i  = 1:length(polylineX)
        distance2p(i) = sqrt((point(1)-polylineX(i))^2+(point(2)-polylineY(i))^2);
    end;
    [minimum2p,I2p] = min(distance2p);
    clear k pos i

    % now we need to choose which segment of the polyline to compare our 'point' with. These
    % segments are either adjacent to that point of the polyline, which is the closest
    % to the 'point' of interest, or the closest to the 'point' segment, which
    % has an intersection with the normale pointing from the 'point'.

    if I2p == 1 % if the 'point' is near the start of polyline
        if exist('minimum2s','var')
            if I2p == I2s
                chosen = I2p;
            else
                chosen = I2s;
            end;
        else
            chosen = I2p;
        end;

    elseif I2p == length(polylineX) % if the 'point' is near the end of polyline
        if exist('minimum2s','var')
            if I2s == I2p-1
                chosen = I2p - 1;
            else
                chosen = I2s;
            end;
        else
            chosen = I2p - 1;
        end;
    else
        if exist('minimum2s','var')
            if I2p == I2s
                chosen = I2p;

            else
                chosen = I2s;
            end;
        else
                pos1 =  sign((polylineX(I2p) - polylineX(I2p-1)) * (point(2) - polylineY(I2p-1)) -...
                 (polylineY(I2p) - polylineY(I2p-1)) * (point(1) - polylineX(I2p-1)));
                % position of the second segment relative to the first segment
                pos2 =  sign((polylineX(I2p) - polylineX(I2p-1)) * (polylineY(I2p+1) - polylineY(I2p-1)) -...
                             (polylineY(I2p) - polylineY(I2p-1)) * (polylineX(I2p+1) - polylineX(I2p-1)));
                if (pos1 == 1 && pos2 == 1) || (pos1 == -1 && pos2 == -1)
                    chosen = I2p;
                elseif pos1 == 0 || pos2 == 0
                    chosen = I2p;
                else
                    chosen = I2p - 1;
                end;
        end;
    end;

    position = sign((polylineX(chosen+1) - polylineX(chosen)) * (point(2) - polylineY(chosen)) -...
                    (polylineY(chosen+1) - polylineY(chosen)) * (point(1) - polylineX(chosen)));

    if position == 1
        plot(point(1),point(2),'r.','MarkerSize',5)
    elseif position == -1;
        plot(point(1),point(2),'.','Color',[0.9 0.9 0.9],'MarkerSize',5) % gray color
    elseif position == 0
        plot(point(1),point(2),'k.','MarkerSize',5)
    elseif position == -2
        plot(point(1),point(2),'g.','MarkerSize',5)
    end;

    pause(0.00000001)
    clear chosen  count_s distance2p distance 2s flag I2p I2s minimum2p minimum2s point pos1 pos2 position

end;
end;

Answer 1

我的想法是从点到线整合角度。积分从无穷远的一侧开始，然后从所有点到无穷远的另一侧。那只是一堆 atan2 函数。

曲线的边由积分符号确定。即使曲线重叠，这也应该有效。

在python中：

from math import atan2,pi
#import matplotlib.pyplot as plt

# difference of two angles should be always -pi < x < pi
def fixed_angle(a):
    if a > pi:
       return a - 2*pi
    elif a < (-1*pi):
       return a + 2*pi
    assert(-1*pi < a < pi) # just check, to be sure
    return a

# polyline
xs = [0, 1, 2, 3, 4, 5, 6, 7, 8];
ys = [-10, 20, -13, 18, -17, 16, -21, 23, -25];

# start point
x = 4
y = 0

#from first two points
angle_start = atan2(ys[0]-ys[1],xs[0]-xs[1])
#last angle is angle of last section
angle_end = atan2(ys[-1]-ys[-2],xs[-1]-xs[-2]) 

integ = 0
prev =  angle_start
for i in range(len(xs)):
    a = atan2(ys[i]-y,xs[i]-x)
    integ += fixed_angle(a-prev)
    prev = a

integ += fixed_angle(angle_end - prev)

if integ > 0:
    print("point is left")
else:
    print("point is right")

#plt.plot(xs,ys)
#plt.show()

Answer 2

基于“winding number”概念的替代解决方案怎么样？这将涉及计算从测试点查看的边界折线的所有线段所对向的总角度。

如果您想象这条折线在无穷远处被半圆延伸以形成闭合轮廓，那么对于轮廓外的点，角度和将为零，但对于轮廓内的点，角度和将为 2*Pi。显然，半圆是在左半平面还是在右半平面将决定直线的哪一侧是“内”或“外”。如果从角度和中排除半圆本身，因为它贡献了 +/-Pi，那么直线一侧的点应该有 +Pi 的角度和，而另一边应该给出 -Pi。是 +Pi 还是 -Pi 取决于折线中顶点序列定义的方向。

这种方法的一个微妙之处是计算折线的每段所对的（有符号的）角度。我相信这可以通过查看将测试点连接到给定段的两端的两个向量来完成。通过获取这些向量的点积（并除以它们的长度），您可以计算该段所对角度的余弦值。在正角和负角之间取反余弦将是不明确的。但是，通过取两个向量的cross-product（将它们的 z 分量视为零），您可以计算角度的正弦，其符号将允许您计算反余弦的正确分支。

Answer 3

一种不计算角度的方法是考虑由测试点和连续折线段构成的三角形的有符号面积。

三角形的有符号区域可以理解为（按照您采用的约定），一个顶点位于相对线的“左”或“右”。它由公式0.5 * (-x2*y1 + x3*y1 + x1*y2 - x3*y2 -x1*y3 + x2*y3) 计算三个顶点(xi, yi)。至关重要的是，您遍历顶点的顺序会影响符号。

假设你的折线是一个列表[v1, v2, ...., vn]。对于每个有序的相邻顶点对(v_i, v_(i+1))，计算三角形(v_i, your-test-point, v_(i+1)) 的有符号面积。似乎要关心“边”决策的是绝对面积最小的三角形：即该点“最接近”折线的位置。根据您的约定，通过考虑其符号区域来确定该三角形相对于您的折线是左侧还是右侧。

*edit - 面积为零的三角形表示您的测试点与折线段共线；你需要单独测试它是否真的在段上。

Answer 4

我的解决方案在图中不言自明。

1. 投影点，例如P₁ 到折线 [O₁,O₂,O₃,O_{4 sub>]，投影点为 P^'1，位于 O1O2 线处。}
2. 如果 O₁P^'₁P₁ 和 P₁P ^'₁O₂ 都是逆时针方向的，那么点 P₁ 在左边折线。
3. 如果点（例如P₂）在折线的右侧，那么O₃P^'_{2P2 和 P2P^'2O4 都是顺时针方向。}