使用rioxarray裁剪NetCDF：如何缩小输出文件至流域区域尺寸？

我尝试用一份流域Shapefile裁剪全球降水NetCDF文件，所用代码如下：

import geopandas as gpd
import rioxarray
import xarray as xr
import matplotlib.pyplot as plt
from shapely.geometry import mapping

# Load shapefile
shapefile_path = "my_shapefile.shp"
shapefile = gpd.read_file(shapefile_path, crs="epsg:4326")

# Load NetCDF file
netcdf_path = "my_netcdffile.nc"
ds = xr.open_dataset(netcdf_path)

pre = ds["pre"] # selecting precipitation

# clipping
pre.rio.set_spatial_dims(x_dim="lon", y_dim="lat", inplace=True)
pre.rio.write_crs("epsg:4326", inplace=True)
clipped_nc = pre.rio.clip(
    shapefile.geometry.apply(mapping), pre.rio.crs, drop=False
)

clipped_nc.to_netcdf("clipped_netcdf.nc") # saving nc file

裁剪后输出的clipped_netcdf.nc仍保留原始全球文件的尺寸，仅流域区域有有效数据。如何让输出文件仅包含流域所在区域，从而缩小文件体积？

解决方案 #

1. 修改rio.clip的drop参数 #

问题出在你调用rio.clip时设置了drop=False，这个参数会强制保留原始数据的所有空间维度（即全球范围的经纬度），即使大部分区域是NaN。只需将其改为drop=True，rioxarray就会自动删除全为NaN的行和列，只保留包含有效数据的区域：

修改后的裁剪代码片段：

clipped_nc = pre.rio.clip(
    shapefile.geometry.apply(mapping), pre.rio.crs, drop=True
)

2. 先预切片NetCDF再裁剪（推荐大文件使用） #

如果你的全球NetCDF文件非常大，先通过Shapefile的边界范围对数据进行预切片，再执行裁剪操作，能大幅提升效率并进一步缩小输出文件体积：

完整代码示例：

import geopandas as gpd
import rioxarray
import xarray as xr
from shapely.geometry import mapping

# Load shapefile
shapefile_path = "my_shapefile.shp"
shapefile = gpd.read_file(shapefile_path, crs="epsg:4326")

# 获取流域的经纬度边界
min_lon, min_lat, max_lon, max_lat = shapefile.total_bounds

# Load NetCDF file并预切片到流域范围
netcdf_path = "my_netcdffile.nc"
ds = xr.open_dataset(netcdf_path)
pre = ds["pre"].sel(
    lon=slice(min_lon - 0.1, max_lon + 0.1),  # 留少量缓冲避免边界数据丢失
    lat=slice(min_lat - 0.1, max_lat + 0.1)
)

# 设置空间维度和CRS
pre.rio.set_spatial_dims(x_dim="lon", y_dim="lat", inplace=True)
pre.rio.write_crs("epsg:4326", inplace=True)

# 执行裁剪并丢弃无效区域
clipped_nc = pre.rio.clip(
    shapefile.geometry.apply(mapping), pre.rio.crs, drop=True
)

# 保存结果
clipped_nc.to_netcdf("clipped_netcdf.nc")

这里给边界加了0.1的缓冲，是为了避免因为Shapefile和NetCDF的网格对齐问题，导致流域边缘的有效数据被误删，你可以根据数据的分辨率调整缓冲值。

内容的提问来源于stack exchange，提问作者jei L