CODEWITHSUNDEEP
pythonmaskrasterioextent
justin fiedler
justin fiedler

Reputation: 23

Mask raster by extent in Python using rasterio

I want to clip one raster based on the extent of another (smaller) raster. First I determine the coordinates of the corners of the smaller raster using

import rasterio as rio
import gdal
from shapely.geometry import Polygon

src = gdal.Open(smaller_file.tif)
ulx, xres, xskew, uly, yskew, yres  = src.GetGeoTransform()
lrx = ulx + (src.RasterXSize * xres)
lry = uly + (src.RasterYSize * yres)
geometry = [[ulx,lry], [ulx,uly], [lrx,uly], [lrx,lry]]

This gives me the following output geometry = [[-174740.0, 592900.0], [-174740.0, 2112760.0], [900180.0, 2112760.0], [900180.0, 592900.0]]. (Note that the crs is EPSG: 32651). Now I would like to clip the larger file using rio.mask.mask(). According to the documentation, the shape variable should be GeoJSON-like dict or an object that implements the Python geo interface protocol (such as a Shapely Polygon). Therefore I create a Shapely Polygon out of the variable geometry, using

roi = Polygon(geometry)

Now everything is ready to use the rio.mask() function.

output = rio.mask.mask(larger_file.tif, roi, crop = True)

But this gives me the following error

TypeError: 'Polygon' object is not iterable

What do I do wrong? Or if someone knows a more elegant way to do it, please let me know.

(Unfortunately I cannot upload the two files since they're too large)

Upvotes: 1

Views: 2619

Answers (2)

ciranzo
ciranzo

Reputation: 82

You can perform the clip entirely with rasterio package. Just open the larger raster using the bounding box of the smaller one.

import rasterio
from rasterio.windows import from_bounds

small_img = rasterio.open("small_img.tif")

# Open the larger raster
with rio.open("big_img.tif") as big_src:
    # Get a window that corresponds to the smaller raster's bounds
    small_window = from_bounds(*small_img.bounds, transform=big_src.transform)
    
    # Read the data from the large raster using the bbox of small raster
    output = big_src.read(window=small_window)

Upvotes: 0

Julian_P
Julian_P

Reputation: 101

I found your question when I needed to figure out this kind of clipping myself. I got the same error and fixed it the following way:

rasterio.mask expects a list of features, not a single geometry. So the algorithm wants to run masking over several features bundled in an iterable (e.g. list or tuple) so we need to pass it our polygon within a list (or tuple) object. The code you posted works after following change:

roi = [Polygon(geometry)]

All we have to do is to enclose the geometry in a list/tuple and then rasterio.mask works as expected.

Upvotes: 0

Related Questions