如何从含非维度坐标的xarray数据集选取指定区域并保存？

问题根源 #

ORAS5的nav_lat和nav_lon是非维度坐标（不属于数据集的维度），且属于不规则网格——每个(x,y)对应唯一经纬度值，没有规则序列。直接用sel()这类基于维度坐标的切片方法时，xarray无法将经纬度值映射到x/y索引上，所以会出现从起始位置错误切片的情况。

可行解决方法 #

方法1：用xarray.where直接筛选区域 #

这是最直接的方式，通过经纬度条件保留目标网格点，还能自动删除全NaN的无效维度：

# 替换成你需要的经纬度范围
lon_min, lon_max = -10, 30
lat_min, lat_max = 30, 60

# 筛选数据集
ds_filtered = ds.where(
    (ds.nav_lon >= lon_min) & (ds.nav_lon <= lon_max) &
    (ds.nav_lat >= lat_min) & (ds.nav_lat <= lat_max),
    drop=True
)

# 保存结果
ds_filtered.to_netcdf("filtered_oras5_data.nc")

如果数据集特别大，提前计算掩码能优化性能（基于dask延迟计算，不会占用过多内存）：

# 先生成掩码
mask = (ds.nav_lon >= lon_min) & (ds.nav_lon <= lon_max) & (ds.nav_lat >= lat_min) & (ds.nav_lat <= lat_max)
# 先删除全NaN的x/y列，再应用掩码
ds_filtered = ds.isel(x=mask.any(dim='y'), y=mask.any(dim='x'))
ds_filtered = ds_filtered.where(mask, drop=True)

方法2：先找x/y索引再切片（适合近似规则的区域） #

如果目标区域在网格中是连续的，可以先定位对应的x/y索引边界，再用isel()切片：

# 找出符合经纬度范围的x、y索引
x_idx = ((ds.nav_lon >= lon_min) & (ds.nav_lon <= lon_max)).any(dim='y').values.nonzero()[0]
y_idx = ((ds.nav_lat >= lat_min) & (ds.nav_lat <= lat_max)).any(dim='x').values.nonzero()[0]

# 先按索引切片，缩小范围
ds_filtered = ds.isel(x=slice(x_idx.min(), x_idx.max()+1), y=slice(y_idx.min(), y_idx.max()+1))

# 再筛选掉区域内的无效点（可选）
ds_filtered = ds_filtered.where(
    (ds_filtered.nav_lon >= lon_min) & (ds_filtered.nav_lon <= lon_max) &
    (ds_filtered.nav_lat >= lat_min) & (ds_filtered.nav_lat <= lat_max),
    drop=True
)

方法3：用regionmask按地理区域筛选（适合特定海区/国家） #

如果需要按预设地理区域（如海盆、国家）筛选，用regionmask更便捷：

import regionmask

# 示例：用全球海盆掩码，筛选北大西洋（编号6），替换成你需要的区域
mask = regionmask.defined_regions.natural_earth.ocean_basins_50m.mask(ds.nav_lon, ds.nav_lat)
ds_filtered = ds.where(mask == 6, drop=True)

注意事项 #

• ORAS5采用dask分块数组，所有操作均为延迟计算，保存时xarray会自动处理分块计算，无需将全量数据加载到内存。
• 如果保存时出现内存不足，可先调整分块大小：ds_filtered = ds_filtered.chunk({'x': 256, 'y': 256})，再执行保存。

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