Python：查找点是否位于多边形的边界上

Question

我有一个点-i，我希望创建一个函数来了解该点是否位于多边形的边界上。

使用：

def point_inside_polygon(x, y, poly):
    """Deciding if a point is inside (True, False otherwise) a polygon,
    where poly is a list of pairs (x,y) containing the coordinates
    of the polygon's vertices. The algorithm is called the 'Ray Casting Method'"""
    n = len(poly)
    inside = False
    p1x, p1y = poly[0]
    for i in range(n):
        p2x, p2y = poly[i % n]
        if y > min(p1y, p2y):
            if y <= max(p1y, p2y):
                if x <= max(p1x, p2x):
                    if p1y != p2y:
                        xinters = (y-p1y) * (p2x-p1x) / (p2y-p1y) + p1x
                    if p1x == p2x or x <= xinters:
                        inside = not inside
        p1x, p1y = p2x, p2y
    return inside

我只能知道这些点是否位于多边形内。

poly = [(0,0), (2,0), (2,2), (0,2)]
point_inside_polygon(1,1, poly)
True
point_inside_polygon(0,0, poly)
false
point_inside_polygon(2,0, poly)
False
point_inside_polygon(2,2, poly)
True
point_inside_polygon(0,2, poly)
True

如何编写一个函数来查找一个点是否位于多边形的边界上？

Answer 1

每个人都把事情复杂化了。这是多边形上的一个短点，假设您有一个距离函数和一个小的 EPSILON。

def pointOnPolygon(point, poly):
    for i in range(len(poly)):
        a, b = poly[i - 1], poly[i]
        if abs(dist(a, point) + dist(b, point) - dist(a, b)) < EPSILON:
            return true
    return false

Answer 2

在以下位置找到了解决方案：http://geospatialpython.com/2011/08/point-in-polygon-2-on-line.html

代码如下：

# Improved point in polygon test which includes edge
# and vertex points

def point_in_poly(x,y,poly):

   # check if point is a vertex
   if (x,y) in poly: return "IN"

   # check if point is on a boundary
   for i in range(len(poly)):
      p1 = None
      p2 = None
      if i==0:
         p1 = poly[0]
         p2 = poly[1]
      else:
         p1 = poly[i-1]
         p2 = poly[i]
      if p1[1] == p2[1] and p1[1] == y and x > min(p1[0], p2[0]) and x < max(p1[0], p2[0]):
         return "IN"

   n = len(poly)
   inside = False

   p1x,p1y = poly[0]
   for i in range(n+1):
      p2x,p2y = poly[i % n]
      if y > min(p1y,p2y):
         if y <= max(p1y,p2y):
            if x <= max(p1x,p2x):
               if p1y != p2y:
                  xints = (y-p1y)*(p2x-p1x)/(p2y-p1y)+p1x
               if p1x == p2x or x <= xints:
                  inside = not inside
      p1x,p1y = p2x,p2y

   if inside: return "IN"
   else: return "OUT"

# Test a vertex for inclusion
polygon = [(-33.416032,-70.593016), (-33.415370,-70.589604),
(-33.417340,-70.589046), (-33.417949,-70.592351),
(-33.416032,-70.593016)]
lat= -33.416032
lon= -70.593016

print point_in_poly(lat, lon, polygon)

# test a boundary point for inclusion
poly2 = [(1,1), (5,1), (5,5), (1,5), (1,1)]
x = 3
y = 1
print point_in_poly(x, y, poly2)

Answer 3

对于每一对相邻的顶点A,B：

1. 构造一个从A到B的向量，称之为p

2. 现在构造一个从 A 到测试点 X 的向量，称之为 q

3. 一对向量的点积公式是 p.q = |p||q|cosC 其中 C 是向量之间的角度。

4. 所以如果 p.q/|p||q| == 1 那么点 AX 和 AB 是共线的。在计算机上工作，您需要 1 - p.q/|p||q|

5. 还需要检查 |q|

如果 4&5 为真，则您的点在边界上。

编辑

我认为我已经看到的另一种方法是获取测试点 X，并通过 X 构建一条垂直于 A 和 B 之间的线的线。找到这条线和线 A->B 交叉的位置。算出从 X 到此交叉点的距离，如果足够小，则认为该点在直线上。

编辑——有趣的小练习！

由于我记错了一些数学，之前发布了一些错误的代码。 在回家的火车上玩 Pythonista 并想出了这个似乎基本上可行的方法。因为在 iPad 上编辑帖子很痛苦，所以留下了数学证明！

没有做太多测试，没有测试除以零等，警告用户。

    # we determine the point of intersection X between
    # the line between A and B and a line through T
    # that is perpendicular to the line AB (can't draw perpendicular
    # in ascii, you'll have to imagine that angle between AB and XT is 90
    # degrees.
    #
    #       B
    #      /
    #.    X  
    #    / \
    #   /   T
    #  A
    # once we know X we can work out the closest the line AB
    # comes to T, if that distance is 0 (or small enough)
    # we can consider T to be on the line
    import math


    # work out where the line through test point t
    # that is perpendicular to ab crosses ab
    #
    # inputs must be 2-tuples or 2-element lists of floats (x,y)
    # returns (x,y) of point of intersection
    def intersection_of_perpendicular(a,b,t):

    if a[0] == b[0]:
            return (a[0],t[1])

    if a[1] == b[1]:
            return (t[0],a[1])

    m = (a[1] - b[1])/(a[0] - b[0]) #slope of ab

    x_inter = (t[1] - a[1] + m*a[0] + (1/m)*t[0])*m/(m**2 + 1)
    y_inter = m*(x_inter - a[0]) + a[1]
    y_inter2 = -(1/m)*(x_inter - t[0]) + t[1]

    #print '...computed ',m,(x_inter, y_inter), y_inter2
    return (x_inter, y_inter)

    # basic Pythagorean formula for distance between two points
    def distance(a,b):
        return math.sqrt( (a[0]-b[0])**2 + (a[1]-b[1])**2 )

    # check if a point is within the box defined by a,b at
    # diagonally opposite corners
    def point_in_box(a,b,t):
        xmin = min(a[0],b[0])
        xmax = max(a[0],b[0])
        ymin = min(a[1],b[1])
        ymax = max(a[1],b[1])

        x_in_bounds = True
        if xmax != xmin:
            x_in_bounds = xmin <= t[0] <= xmax
        y_in_bounds = True
        if ymax != ymin:
            y_in_bounds = ymin <= t[1] <= ymax
        return x_in_bounds and y_in_bounds

    # determine if point t is within 'tolerance' distance
    # of the line between a and b
    # returns Boolean
    def is_on_line_between(a,b,t,tolerance=0.01):
        intersect = intersection_of_perpendicular(a,b,t)
        dist = distance(intersect, t)
        in_bounds = point_in_box(a,b,t)
        return in_bounds and (dist < tolerance)


a = (0,0)
b = (2,2)
t = (0,2)

p = intersection_of_perpendicular(a,b,t)
bounded = point_in_box(a,b,t)
print 'd ',distance(p,t), ' p ',p, bounded

a = (0,2)
b = (2,2)
t = (1,3)

p = intersection_of_perpendicular(a,b,t)
bounded = point_in_box(a,b,t)
print 'd ',distance(p,t),' p ',p, bounded

a = (0.0,2.0)
b = (2.0,7.0)
t = (1.7,6.5)

p = intersection_of_perpendicular(a,b,t)
bounded = point_in_box(a,b,t)
on = is_on_line_between(a,b,t,0.2)
print 'd ',distance(p,t),' p ',p, bounded,on 

Answer 4

对于凸多边形，您可以在 O(log n) 时间内解决问题，方法是按顺时针顺序排列顶点并为每个顶点存储顶点与多边形内部的点 c 之间的角度。然后对于查询点 x，您获得从 c 到 x 的角度，并通过二分搜索找到唯一的相邻顶点对 (v1,v2)，使得到 x 的角度介于到 v1 和到 v2 的角度之间。那么 x 要么在边 (v1,v2) 上，要么 x 不在边界上。

如果您有一个更复杂的多边形，那么您可以尝试通过添加一些内部边将多边形分解为凸多边形的并集（例如，首先三角剖分，然后删除边以获得更大的凸多边形）；如果你得到的凸多边形的数量很小（比如 k），那么你可以测试每个凸多边形以查看一个点是否在边缘上，并且总运行时间为 O(k lg n)，其中 n 是总数多边形中的顶点数。

或者，如果您不担心使用额外的空间，并且您确实想快速确定您是否在一条边上，那么您可以通过沿每条边添加额外的“顶点”将每条边分成等间距的段;很难说多少就足够了（听起来像一个有趣的数学问题），但很明显，如果你在每条边上添加足够多的额外顶点，那么你可以通过从你的集合中找到最近的邻居来判断一个点必须位于哪条边上顶点（原始顶点和您添加的顶点），然后只需测试最近邻顶点所在的一或两条边。如果您使用二维 kd-tree（您将构建树作为预处理步骤，然后该树支持快速 k-最近邻查询），您可以非常快速地在 2-d 中找到 k-最近邻，并且kd-tree 树只使用线性空间。

Answer 5

将问题分解为三个步骤可能会有所帮助：

1. 写一个可以判断if a point is on a line segment的函数。
2. 计算构成多边形边界的所有线段。
3. 如果点位于任何线段上，则该点位于边框上。

这里有一些 python 代码，假设你已经为isPointOnLineSegmentBetweenPoints 编写或找到了合适的候选人：

def pointOnPolygon(point, polygonVertices):
    n = len(polygonVertices)
    for i in range(n):
        p1 = polygonVertices[i]
        p2 = polygonVertices[-n+i+1]
        if isPointOnLineSegmentBetweenPoints(point, p1, p2):
            return true
    return false

Answer 6

我没有测试过这个，但总体思路是：

def pointOnBorder(x, y, poly):
    n = len(poly)
    for(i in range(n)):
        p1x, p1y = poly[i]
        p2x, p2y = poly[(i + 1) % n]
        v1x = p2x - p1x
        v1y = p2y - p1y #vector for the edge between p1 and p2
        v2x = x - p1x
        v2y = y - p1y #vector from p1 to the point in question
        if(v1x * v2y - v1y * v2x == 0): #if vectors are parallel 
            if(v2x / v1x > 0): #if vectors are pointing in the same direction
                if(v1x * v1x + v1y * v1y >= v2x * v2x + v2y * v2y): #if v2 is shorter than v1
                    return true
    return false