Create point grid inside a shapefile using python

Question

I am working on a shapefile in python using geopandas and gdal. I am looking to create meshgrid (with regular 1000m interval points) inside the polygon shapefile. I have reprojected the file so that units can be meters. However, I could not find any direct way to implement this. Can any one guide in this regard?

I am sharing the code, I have tried so far:

from osgeo import gdal, ogr
import numpy as np
import matplotlib.pyplot as plt
import os
import sys
import pandas as pd
import geopandas as gpd
from shapely.geometry import Polygon
source_ds = ogr.Open(r"E:\review paper\sample tb data for recon\descend\tiffbt\alaska_bound.shp")
boundFile =gpd.read_file(r"E:\review paper\sample tb data for recon\descend\tiffbt\alaska_bound.shp")
bound_project = boundFile.to_crs({'init': 'EPSG:3572'})
print(bound_project.crs)
print(bound_project.total_bounds)

The coordinate system and bounding box coordinates are as below (output of above code):

+init=epsg:3572 +type=crs
[-2477342.73003557 -3852592.48050272  1305143.81797914 -2054961.64359753]

Answer 1

Find xmin, xmax, ymin, ymax. Create ranges from xmin to xmax and ymin to ymax with the interval of your spacing. Combine all possible coordinates using meshgrid:

import numpy as np
import geopandas as gpd

shapefile = r"C:\GIS\data\testdata\ak_riks.shp"
spacing = 100000

df = gpd.read_file(shapefile)
xmin, ymin, xmax, ymax = df.total_bounds

xcoords = [i for i in np.arange(xmin, xmax, spacing)]
ycoords = [i for i in np.arange(ymin, ymax, spacing)]

pointcoords = np.array(np.meshgrid(xcoords, ycoords)).T.reshape(-1, 2) #A 2D array like [[x1,y1], [x1,y2], ...
points = gpd.points_from_xy(x=pointcoords[:,0], y=pointcoords[:,1])
grid = gpd.GeoSeries(points, crs=df.crs)
grid.name = 'geometry'

#If you just want to points inside polygons:
gridinside = gpd.sjoin(gpd.GeoDataFrame(grid), df[['geometry']], how="inner")

#Plot
import matplotlib.pyplot as plt
fig, ax = plt.subplots(figsize=(15, 15))
df.plot(ax=ax, alpha=0.7, color="pink", edgecolor='red', linewidth=1)
grid.plot(ax=ax, markersize=150, color="blue")
gridinside.plot(ax=ax, markersize=100, color="yellow")