numpy:在不创建掩码的情况下获取多边形内的索引
numpy: get the indices within a polygon without creating a mask
我有一个格式如下的多边形矢量 (x1,y1,x2,y2, ...., xn,yn)
。例如,考虑这个多边形数组:
polyPoints = [3,5,7,8,9,5]
如何获取从这些点生成的多边形内的所有索引(或坐标)?
到目前为止,我看到的答案要求您先创建 2D 蒙版,然后才能获得多边形内的索引。
你可以使用 scikit-image:
import numpy as np
from skimage.draw import polygon
points = [3,5,7,8,9,5]
r, c = polygon(points[1::2], points[::2])
print(r, c)
输出是:
[5 5 5 5 5 5 6 6 6 6 7 7] [3 4 5 6 7 8 5 6 7 8 6 7]
使用口罩可能是您所能达到的最有效的方法。这是一些效率相当低的算法,但可能可以优化以接近掩码方法。这本质上做了一个掩码,但在线上。
方法是:
求所有边的直线方程
查找边界框
对于边界框内的每个 y(或 x,以较小者为准),计算在该 y 处与水平线 (y=yi) 相交的边,并找出它们与 x 相交的边.
对于边界框内的每个 x,找出 x 右侧与直线 y=yi 相交的边数。如果边数为奇数,则点 (x,y) 位于多边形内部。
它确实适用于简单的正方形几何体。
import numpy as np
# taken from:
def line(p1, p2):
A = (p1[1] - p2[1])
B = (p2[0] - p1[0])
C = (p1[0]*p2[1] - p2[0]*p1[1])
return A, B, -C
def intersection(L1, L2):
D = L1[0] * L2[1] - L1[1] * L2[0]
Dx = L1[2] * L2[1] - L1[1] * L2[2]
Dy = L1[0] * L2[2] - L1[2] * L2[0]
if D != 0:
x = Dx / D
y = Dy / D
return x,y
else:
return False
# polyPoints = np.array([0, 0, 4, 0,4, 4, 0, 4])
polyPoints = np.array([[3,5,7,8,9,5]])
polyPoints = polyPoints.reshape(-1, 2)
npoints = polyPoints.shape[0]
polyEgdes = []
for i in range(npoints):
point1, point2 = polyPoints[i, :], polyPoints[(i+1) % npoints, :]
polyEgdes.append(line(point1, point2))
# bounding box
boundingBox = np.vstack((polyPoints.min(axis=0), polyPoints.max(axis=0)))
inside_points = []
for y in range(boundingBox[0, 1], boundingBox[1, 1]):
x_intersect = []
for l in polyEgdes:
# y_ins should be same as y
insect_point = intersection(l, [0, y, 0])
if insect_point:
x_intersect.append(insect_point[0])
x_intersect = np.array(x_intersect)
for x in range(boundingBox[0, 0]+1, boundingBox[1, 0]-1):
x_int_points = x_intersect[(x_intersect - x) >= 0]
if len(x_int_points) % 2 == 1:
inside_points.append((x, y))
print(inside_points)
我有一个格式如下的多边形矢量 (x1,y1,x2,y2, ...., xn,yn)
。例如,考虑这个多边形数组:
polyPoints = [3,5,7,8,9,5]
如何获取从这些点生成的多边形内的所有索引(或坐标)?
到目前为止,我看到的答案要求您先创建 2D 蒙版,然后才能获得多边形内的索引。
你可以使用 scikit-image:
import numpy as np
from skimage.draw import polygon
points = [3,5,7,8,9,5]
r, c = polygon(points[1::2], points[::2])
print(r, c)
输出是:
[5 5 5 5 5 5 6 6 6 6 7 7] [3 4 5 6 7 8 5 6 7 8 6 7]
使用口罩可能是您所能达到的最有效的方法。这是一些效率相当低的算法,但可能可以优化以接近掩码方法。这本质上做了一个掩码,但在线上。
方法是:
求所有边的直线方程
查找边界框
对于边界框内的每个 y(或 x,以较小者为准),计算在该 y 处与水平线 (y=yi) 相交的边,并找出它们与 x 相交的边.
对于边界框内的每个 x,找出 x 右侧与直线 y=yi 相交的边数。如果边数为奇数,则点 (x,y) 位于多边形内部。
它确实适用于简单的正方形几何体。
import numpy as np
# taken from:
def line(p1, p2):
A = (p1[1] - p2[1])
B = (p2[0] - p1[0])
C = (p1[0]*p2[1] - p2[0]*p1[1])
return A, B, -C
def intersection(L1, L2):
D = L1[0] * L2[1] - L1[1] * L2[0]
Dx = L1[2] * L2[1] - L1[1] * L2[2]
Dy = L1[0] * L2[2] - L1[2] * L2[0]
if D != 0:
x = Dx / D
y = Dy / D
return x,y
else:
return False
# polyPoints = np.array([0, 0, 4, 0,4, 4, 0, 4])
polyPoints = np.array([[3,5,7,8,9,5]])
polyPoints = polyPoints.reshape(-1, 2)
npoints = polyPoints.shape[0]
polyEgdes = []
for i in range(npoints):
point1, point2 = polyPoints[i, :], polyPoints[(i+1) % npoints, :]
polyEgdes.append(line(point1, point2))
# bounding box
boundingBox = np.vstack((polyPoints.min(axis=0), polyPoints.max(axis=0)))
inside_points = []
for y in range(boundingBox[0, 1], boundingBox[1, 1]):
x_intersect = []
for l in polyEgdes:
# y_ins should be same as y
insect_point = intersection(l, [0, y, 0])
if insect_point:
x_intersect.append(insect_point[0])
x_intersect = np.array(x_intersect)
for x in range(boundingBox[0, 0]+1, boundingBox[1, 0]-1):
x_int_points = x_intersect[(x_intersect - x) >= 0]
if len(x_int_points) % 2 == 1:
inside_points.append((x, y))
print(inside_points)