Example - Clip

[1]:
import rioxarray

%matplotlib inline

Load in xarray dataset

See docs for rioxarray.open_rasterio

Notes:

  • chunks=True only works if you have a dask installed. Otherwise, you can skip this.
  • masked=True will convert from integer to float64 and fill with NaN. If this behavior is not desired, you can skip this.
[2]:
xds = rioxarray.open_rasterio(
    "../../test/test_data/compare/small_dem_3m_merged.tif",
    masked=True,
    chunks=True,
)
[3]:
xds
[3]:
<xarray.DataArray (band: 1, y: 245, x: 574)>
dask.array<shape=(1, 245, 574), dtype=float64, chunksize=(1, 245, 574)>
Coordinates:
  * band         (band) int64 1
  * y            (y) float64 4.616e+06 4.616e+06 ... 4.615e+06 4.615e+06
  * x            (x) float64 4.25e+05 4.251e+05 ... 4.268e+05 4.268e+05
    spatial_ref  int64 0
Attributes:
    transform:     (3.0, 0.0, 425047.68381405267, 0.0, -3.0, 4615780.040546387)
    scales:        (1.0,)
    offsets:       (0.0,)
    grid_mapping:  spatial_ref
[4]:
xds.plot()
[4]:
<matplotlib.collections.QuadMesh at 0x7f1c4807acc0>
../_images/examples_clip_geom_5_1.png

Clip using a geometry

By default, it assumes that the CRS of the geometry is the same as the CRS of the dataset. If it is different, make sure to pass in the CRS of the geometry.

API Reference for rio.clip:

[5]:
geometries = [
    {
        'type': 'Polygon',
        'coordinates': [[
            [425499.18381405267, 4615331.540546387],
            [425499.18381405267, 4615478.540546387],
            [425526.18381405267, 4615478.540546387],
            [425526.18381405267, 4615331.540546387],
            [425499.18381405267, 4615331.540546387]
        ]]
    }
]
[6]:
clipped = xds.rio.clip(geometries)
[7]:
clipped
[7]:
<xarray.DataArray (band: 1, y: 49, x: 9)>
dask.array<shape=(1, 49, 9), dtype=float64, chunksize=(1, 49, 9)>
Coordinates:
  * band         (band) int64 1
  * y            (y) float64 4.615e+06 4.615e+06 ... 4.615e+06 4.615e+06
  * x            (x) float64 4.255e+05 4.255e+05 ... 4.255e+05 4.255e+05
    spatial_ref  int64 0
Attributes:
    transform:     (3.0, 0.0, 425500.68381405267, 0.0, -3.0, 4615477.04054638...
    scales:        (1.0,)
    offsets:       (0.0,)
    grid_mapping:  spatial_ref
[8]:
clipped.plot()
[8]:
<matplotlib.collections.QuadMesh at 0x7f1c39904400>
../_images/examples_clip_geom_10_1.png
[9]:
clipped.rio.to_raster("clipped.tif", compress='LZMA', tiled=True, dtype="int32")

Clip using a GeoDataFrame

API Reference for rio.clip:

[10]:
import geopandas
from shapely.geometry import box, mapping
[11]:
geodf = geopandas.GeoDataFrame(
    geometry=[
        box(425499.18381405267, 4615331.540546387, 425526.18381405267, 4615478.540546387)
    ],
    crs=xds.rio.crs.to_dict()
)
[12]:
clipped = xds.rio.clip(geodf.geometry.apply(mapping), geodf.crs, drop=False, invert=True)
[13]:
clipped.plot()
[13]:
<matplotlib.collections.QuadMesh at 0x7f1c3938fcf8>
../_images/examples_clip_geom_16_1.png
[14]:
clipped.rio.to_raster("clipped_invert.tif", compress='LZMA', tiled=True, dtype="int32")