5

我有一个点 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

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

4

6 回答 6

2

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

  1. 编写一个函数来判断一个点是否在线段上
  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
于 2013-07-19T16:24:42.163 回答
1

对于每对相邻顶点 A,B:

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

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

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

  4. 所以如果 pq/|p||q| == 1 那么点 AX 和 AB 是共线的。在计算机上工作,您将需要 1 - pq/|p||q| < some_small_value 取决于您想要的准确度。

  5. 还需要检查 |q| < |p| (即X比B更接近A)

如果 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 
于 2013-07-19T15:38:41.397 回答
1

我没有对此进行测试,但总体思路是:

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
于 2013-07-19T15:42:36.340 回答
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
于 2015-08-16T12:19:44.993 回答
0

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

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

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

于 2013-07-19T18:21:52.757 回答
-1

在以下位置找到了解决方案: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)
于 2013-08-14T13:34:36.203 回答