0

I want to input the coordinates of 4 vertices of a polygon, specify the number of points to divide the polygon's edges (in equal segments), and the goal is to generate a matrix with the coordinates of the grid points inside the polygon.

The following picture might explain better my goal: enter image description here

So in that case I would input the coordinates of the points (P) and specify that I want grid to be 3 by 2. The output would be a 3 by 2 matrix with the coordinates (x,y) of the grid points (N).

I have searched a lot but still couldn't find a way to do this, and honestly I am not at all experienced with Python. I found something using numpy's meshgrid in combination with matplotlib.path's contains_points to create a grid inside a polygon but I have no idea how to get the grid point's coordinates. I saw shapely being used a lot in problems like this but again, I'm not experience in this so some help would be appreciated!

Thank you all in advance!

Vidya Ganesh
  • 788
  • 11
  • 24
raulfilipe127
  • 13
  • 5
  • It's a bit unclear what you need help with. If you need help with the math this isn't the right community. If you need help implementing an algorithm you should post your code in a way that can reproduce where you are getting stuck. https://stackoverflow.com/help/how-to-ask – gph Jan 07 '21 at 22:33

1 Answers1

1

To explain the approach, we have three steps:

  1. Find the ticks on each side of the 4-gon
  2. Find the grid lines
  3. Find the intersection of the grid lines
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

Note: the part of the code for finding the intersection of lines are from here

Amin Gheibi
  • 639
  • 1
  • 8
  • 15
  • Perfect, that's exactly what I was looking for! I just had to slightly change the code to get a different "number of parts on each side of the grid" for vertical and horizontal lines. Anyway, that's a great solution my problem is solved now, thanks for the help! – raulfilipe127 Jan 08 '21 at 14:17
  • You're welcome. If that's the solution please mark it as the correct answer. That helps others when they come to this question in future. – Amin Gheibi Jan 08 '21 at 17:52