I'm working with a GeoDataFrame of points and I need to create a grid of rectangular polygons.

I know how to have bounds of the GeoDataFrame with gdf.geometry.bounds. Now using this information and with two variables length and width (L, W) of rectangular polygon, how can I create a grid?

3 Answers 3


Previous answer gives correct result, but I allowed myself to improve the code, to avoid many unnecessary lines, as well as iterating on indexes instead of values on the list.

import geopandas as gpd
from shapely.geometry import Polygon
import numpy as np
points = gpd.read_file('points.shp')

xmin, ymin, xmax, ymax = points.total_bounds

length = 1000
wide = 1200

cols = list(np.arange(xmin, xmax + wide, wide))
rows = list(np.arange(ymin, ymax + length, length))

polygons = []
for x in cols[:-1]:
    for y in rows[:-1]:
        polygons.append(Polygon([(x,y), (x+wide, y), (x+wide, y+length), (x, y+length)]))

grid = gpd.GeoDataFrame({'geometry':polygons})

The main idea might be the same, but we are creating now much less useless variables and whole code is clearer to understand

  • 3
    When you want to work with non-integer coordinates or grid sizes use np.arange instead of range: cols = list(np.arange(xmin, xmax + wide, wide)) ...
    – quassy
    Commented Feb 27, 2021 at 0:46

There are many solutions.

One of them

enter image description here

 import geopandas as gpd
 from shapely.geometry import Polygon
 import numpy as np
 points = gpd.read_file('points.shp')
 xmin,ymin,xmax,ymax =  points.total_bounds
 width = 2000
 height = 1000
 rows = int(np.ceil((ymax-ymin) /  height))
 cols = int(np.ceil((xmax-xmin) / width))
 XleftOrigin = xmin
 XrightOrigin = xmin + width
 YtopOrigin = ymax
 YbottomOrigin = ymax- height
 polygons = []
 for i in range(cols):
    Ytop = YtopOrigin
    Ybottom =YbottomOrigin
    for j in range(rows):
        polygons.append(Polygon([(XleftOrigin, Ytop), (XrightOrigin, Ytop), (XrightOrigin, Ybottom), (XleftOrigin, Ybottom)])) 
        Ytop = Ytop - height
        Ybottom = Ybottom - height
    XleftOrigin = XleftOrigin + width
    XrightOrigin = XrightOrigin + width

grid = gpd.GeoDataFrame({'geometry':polygons})

enter image description here

You can also truncate the grid (convex hull):

enter image description here

But one of the most interesting is to use the module gpd_lite_toolboox)

  • For calculating number of columns and rows I add abs() because of possibility to have negative coordinates with wgs84 : rows = abs(int(np.ceil((ymax-ymin) / height))) and cols = abs(int(np.ceil((xmax-xmin) / width)))
    – Tim C.
    Commented Jan 26, 2018 at 14:01
  • 1
    @gene I disagree with you since number will always be positive since it is max minus min, not needing abs.
    – gcamargo
    Commented May 18, 2018 at 18:45
  • Where is abs in the answer ?
    – gene
    Commented May 19, 2018 at 19:15

One adjustment to the answer from @gene. otherwise, it worked great for my purposes.

I ran this against some points I had to do some binning with, and noticed that the two points that defined the upper and left limits of the extent (used as the origin point for drawing the first and subsequent polygons) do not actually intersect the resulting grid.

To adjust this I just "buffered" out the x and y min and max values by 1 (1 meter), and that did the trick. I understand this is kind of a "Kentucky Windage" solution though, so if there is a better way, let me know!

xmin, ymin, xmax, ymax = gdf.geometry.total_bounds
xmin = xmin - 1
ymin = ymin - 1
xmax = xmax + 1
ymax = ymax + 1

enter image description here

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.