如何对xarray数据集旋转n度？大尺寸数据集适配方案问询

针对你用scipy.ndimage.rotate处理4200×4200数据集时返回全NaN的问题，核心原因是超大数组一次性加载计算导致内存过载，以下是几种可行的替代方案：

方法1：结合Dask分块优化scipy旋转 #

通过Dask将大数组拆分为小分块并行处理，避免内存溢出，同时保留scipy的旋转逻辑：

import xarray as xr
from scipy import ndimage
import dask.array as da
import numpy as np

# 加载你的大尺寸数据集（示例模拟4200x4200规模）
airtemps = xr.tutorial.open_dataset("air_temperature")
da = airtemps['air'][0].repeat(20, dim='lat').repeat(20, dim='lon')
da = da.chunk({"lat": 500, "lon": 500})  # 分块大小可按需调整

# 执行旋转
rotated_da = xr.apply_ufunc(
    ndimage.rotate,
    da,
    10,  # 旋转角度
    axes=(1, 0),
    reshape=False,
    dask="parallelized",
    output_dtypes=[da.dtype],
    kwargs={"mode": "constant", "cval": np.nan}  # 边界填充方式可自定义
)

# 计算并可视化
rotated_da[::5, ::5].compute().plot(vmin=220, vmax=300)

方法2：用Rasterio做地理空间旋转 #

Rasterio专为地理空间数据设计，天然支持分块处理大文件，通过自定义旋转投影实现数据旋转：

import xarray as xr
import rasterio
from rasterio.warp import reproject, Resampling
from affine import Affine
import numpy as np

# 加载数据
airtemps = xr.tutorial.open_dataset("air_temperature")
da = airtemps['air'][0].repeat(20, dim='lat').repeat(20, dim='lon')

# 原数据的地理变换与CRS（假设为WGS84）
src_res_lon = da.lon[1] - da.lon[0]
src_res_lat = da.lat[1] - da.lat[0]
src_transform = Affine.translation(da.lon[0], da.lat[0]) * Affine.scale(src_res_lon, src_res_lat)
src_crs = 'EPSG:4326'

# 定义以数据中心为原点旋转10度的目标CRS
center_lon = da.lon.mean().item()
center_lat = da.lat.mean().item()
dst_crs = f'+proj=omerc +lat_0={center_lat} +lonc={center_lon} +alpha=10 +k=1 +x_0=0 +y_0=0 +datum=WGS84 +units=m +no_defs'

# 执行重投影（实现旋转）
dst_data = np.empty_like(da)
reproject(
    source=da.data,
    destination=dst_data,
    src_transform=src_transform,
    src_crs=src_crs,
    dst_transform=src_transform,
    dst_crs=dst_crs,
    resampling=Resampling.bilinear
)

# 转回Xarray DataArray
rotated_da = xr.DataArray(dst_data, dims=da.dims, coords=da.coords)
rotated_da[::5, ::5].plot(vmin=220, vmax=300)

方法3：用OpenCV高效旋转 #

OpenCV的warpAffine处理大数组效率更高，适合对性能要求高的场景：

import xarray as xr
import cv2
import numpy as np

# 加载数据
airtemps = xr.tutorial.open_dataset("air_temperature")
da = airtemps['air'][0].repeat(20, dim='lat').repeat(20, dim='lon')
data = da.data

# 生成旋转矩阵，计算旋转后尺寸避免裁剪
angle = 10
h, w = data.shape
center = (w // 2, h // 2)
rot_mat = cv2.getRotationMatrix2D(center, angle, 1.0)

abs_cos = abs(rot_mat[0, 0])
abs_sin = abs(rot_mat[0, 1])
new_w = int(h * abs_sin + w * abs_cos)
new_h = int(h * abs_cos + w * abs_sin)

rot_mat[0, 2] += new_w / 2 - center[0]
rot_mat[1, 2] += new_h / 2 - center[1]

# 执行旋转
rotated_data = cv2.warpAffine(
    data, rot_mat, (new_w, new_h),
    flags=cv2.INTER_LINEAR,
    borderMode=cv2.BORDER_CONSTANT,
    borderValue=np.nan
)

# 转回Xarray DataArray（需重新生成坐标）
rotated_da = xr.DataArray(
    rotated_data,
    dims=['new_lat', 'new_lon'],
    coords={
        'new_lat': np.linspace(da.lat.min(), da.lat.max(), new_h),
        'new_lon': np.linspace(da.lon.min(), da.lon.max(), new_w)
    }
)

rotated_da[::5, ::5].plot(vmin=220, vmax=300)

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