如何获取Python多边形内的网格点坐标?
How to get the coordinates of grid points inside a polygon in Python?
我要输入一个多边形的4个顶点的坐标,指定多边形的边(等分)的点数,目标是生成一个矩阵,里面有网格点的坐标多边形。
下图可以更好地解释我的目标:
所以在那种情况下,我将输入点 (P) 的坐标并指定我希望网格为 3 x 2。输出将是一个 3 x 2 矩阵,其坐标 (x,y) 为网格点数 (N).
我已经搜索了很多,但仍然找不到方法,老实说,我对 Python 一点经验都没有。我发现结合使用 numpy 的 meshgrid 和 matplotlib.path 的 contains_points 在多边形内创建网格,但我不知道如何获取网格点的坐标。我看到 shapely 在类似这样的问题中被大量使用,但同样,我没有这方面的经验,所以一些帮助将不胜感激!
提前谢谢大家!
为了解释这个方法,我们分为三个步骤:
- 找到 4 边形每边的刻度线
- 找到网格线
- 找到网格线的交点
def det(a, b):
return a[0] * b[1] - a[1] * b[0]
ticks = []
A = (1, 2) # wlog., suppose the polygon is ABCD;
B = (3, 2)
C = (3, 4)
D = (1, 4)
polygon = [A, B, C, D]
n = 3 # number of parts on each side of the grid
# we first find ticks on each side of the polygon
for j in range(4): # because we are talking about 4-gons
temp_ticks = []
for i in range(n-1):
t = (i+1)*(1/n)
Ex = polygon[j][0] * (1-t) + polygon[(j+1) % 4][0] * t
Ey = polygon[j][1] * (1-t) + polygon[(j+1) % 4][1] * t
temp_ticks.append((Ex, Ey))
if j < 2:
ticks.append(temp_ticks)
else: # because you are moving backward in this part
temp_ticks.reverse()
ticks.append(temp_ticks)
# then we find lines of the grid
h_lines = []
v_lines = []
for i in range(n-1):
h_lines.append((ticks[0][i], ticks[2][i]))
v_lines.append((ticks[1][i], ticks[3][i]))
# then we find the intersection of grid lines
for i in range(len(h_lines)):
for j in range(len(v_lines)):
line1 = h_lines[i]
line2 = v_lines[j]
xdiff = (line1[0][0] - line1[1][0], line2[0][0] - line2[1][0])
ydiff = (line1[0][1] - line1[1][1], line2[0][1] - line2[1][1])
div = det(xdiff, ydiff)
if div == 0:
raise Exception('lines do not intersect')
d = (det(*line1), det(*line2))
x = det(d, xdiff) / div
y = det(d, ydiff) / div
print(x,y) # this is an intersection point that you want
注:求交点部分代码来自here
我要输入一个多边形的4个顶点的坐标,指定多边形的边(等分)的点数,目标是生成一个矩阵,里面有网格点的坐标多边形。
下图可以更好地解释我的目标:
所以在那种情况下,我将输入点 (P) 的坐标并指定我希望网格为 3 x 2。输出将是一个 3 x 2 矩阵,其坐标 (x,y) 为网格点数 (N).
我已经搜索了很多,但仍然找不到方法,老实说,我对 Python 一点经验都没有。我发现结合使用 numpy 的 meshgrid 和 matplotlib.path 的 contains_points 在多边形内创建网格,但我不知道如何获取网格点的坐标。我看到 shapely 在类似这样的问题中被大量使用,但同样,我没有这方面的经验,所以一些帮助将不胜感激!
提前谢谢大家!
为了解释这个方法,我们分为三个步骤:
- 找到 4 边形每边的刻度线
- 找到网格线
- 找到网格线的交点
def det(a, b):
return a[0] * b[1] - a[1] * b[0]
ticks = []
A = (1, 2) # wlog., suppose the polygon is ABCD;
B = (3, 2)
C = (3, 4)
D = (1, 4)
polygon = [A, B, C, D]
n = 3 # number of parts on each side of the grid
# we first find ticks on each side of the polygon
for j in range(4): # because we are talking about 4-gons
temp_ticks = []
for i in range(n-1):
t = (i+1)*(1/n)
Ex = polygon[j][0] * (1-t) + polygon[(j+1) % 4][0] * t
Ey = polygon[j][1] * (1-t) + polygon[(j+1) % 4][1] * t
temp_ticks.append((Ex, Ey))
if j < 2:
ticks.append(temp_ticks)
else: # because you are moving backward in this part
temp_ticks.reverse()
ticks.append(temp_ticks)
# then we find lines of the grid
h_lines = []
v_lines = []
for i in range(n-1):
h_lines.append((ticks[0][i], ticks[2][i]))
v_lines.append((ticks[1][i], ticks[3][i]))
# then we find the intersection of grid lines
for i in range(len(h_lines)):
for j in range(len(v_lines)):
line1 = h_lines[i]
line2 = v_lines[j]
xdiff = (line1[0][0] - line1[1][0], line2[0][0] - line2[1][0])
ydiff = (line1[0][1] - line1[1][1], line2[0][1] - line2[1][1])
div = det(xdiff, ydiff)
if div == 0:
raise Exception('lines do not intersect')
d = (det(*line1), det(*line2))
x = det(d, xdiff) / div
y = det(d, ydiff) / div
print(x,y) # this is an intersection point that you want
注:求交点部分代码来自here