Last active
July 22, 2022 12:39
-
-
Save robbibt/216d01b92aab9b5f92bd4473e374b270 to your computer and use it in GitHub Desktop.
Rasterio geometry mask example for Stack Exchange
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
import xarray as xr | |
import geopandas as gpd | |
import rasterio | |
# Open your shapefile and xarray object | |
ds = raster_mask | |
gdf = vector_mask | |
# Select shapefile feature you want to analyse | |
# and reproject to same CRS as xarray | |
gdf = gdf.iloc[[0]].to_crs(ds.crs) | |
# create polygon mask | |
mask = rasterio.features.geometry_mask( | |
gdf.geometry, | |
out_shape=ds.geobox.shape, | |
transform=ds.geobox.affine, | |
all_touched=False, | |
invert=False) | |
mask = xr.DataArray(mask, dims=("y", "x")) | |
# mask ds with rasterized gdf | |
ds = ds.where(~mask) |
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
Thanks for the answer @robbibt. Unfortunately, using
odc.geo.xr
does not add thegeobox
function onto my xarray dataset. Is it ok to simply open it like this or do I have to link it to old somehow?