Plotting N number of grid points inside a polygon/shapefile Python

40 views Asked by At

So my goal is to place points or coordinates with a calculated spacing inside a rectangle polygon (however it could also be an irregular polygon). This is the code im working with but it seems to not be able to check if the points are inside the polygon. Is there something im doing wrong? Thank you!

import numpy as np
import geopandas as gpd
from shapely.geometry import Point


df = gpd.read_file(r'C:\Users\User\Desktop\PyTrials\Rectangle.shp')

def grid(row, pointcount):


    finished = 0
    maxiter = 20 
    iterations = 0
    #prep_polygon =prep(row)
    #pointcount=2
    while finished != 1 or iterations<maxiter: #Until there are n points inside the polygon or maxiter is reached
        iterations+=1
        xmin, ymin, xmax, ymax = row.geometry.bounds #Find polygon bounds


        #Find suitable spacing. The 50 buffer is to prevent points at polygon border.
        spacing=(((row.geometry.buffer(50).area/pointcount)**0.5)*0.99)*(1.3-np.random.random()*0.6)
       

        nspacingx=np.ceil((xmax-xmin)/spacing) #Number of points in x
        nspacingy=np.ceil((ymax-ymin)/spacing) #and y

        randomstart=[xmin+spacing*np.random.random(),ymin+spacing*np.random.random()] #A random start coordinate

        xlist=[randomstart[0]+(x*spacing) for x in range(int(nspacingx)+1)] #Create a list of x coordinates
        ylist=[randomstart[1]+(y*spacing) for y in range(int(nspacingy)+1)] #And y

        #Using meshgrid, list all points combinations [Point(x1,y1), Point(x1,y2), ... and check if they are inside the polygon NOT WORKING
        points = [Point(coords) for coords in np.array(np.meshgrid(xlist, ylist)).T.reshape(-1, 2).tolist() if Point(coords).intersects(row.geometry)] #This step is when it is unable to evaluate if #points are inside the polygon

        numberpoints=(len(points))
        print (numberpoints) #when I print this, it gives me 0 even though #the generated points are inside the polygon

        if len(points)==pointcount:
            #print(sum(points))
            return points
            finished =1

df['gridlist'] = df.apply(grid, args=(10,), axis=1) #Apply the function with 10 m spacing
df = df.explode('gridlist') #Explode so each point in gridlist become a row
df['geometry'] = df['gridlist']
df = df.drop('gridlist', axis=1)

df.to_file(r'C:\Users\User\Desktop\PyTrials\PointsRectangle.shp')
0

There are 0 answers