
`regulargridinterpolator` 不支持含 nan 的输入数据,需改用 `griddata` 等无结构插值器:先剔除 nan 网格点及其对应坐标,再以 `(x_flat, y_flat, z_flat)` 形式调用 `griddata(..., method='cubic')`,即可安全完成三次插值并天然处理边界与空洞区域。
当面对带大量 NaN 的二维规则地理网格(如 1°×1° 的纬度-经度场)并对任意散点坐标(如浮标位置或模型输出点)执行高阶平滑插值时,scipy.interpolate.RegularGridInterpolator 的 'cubic' 模式会直接报错 ValueError: Array must not contain infs or nans.——这是因为其底层依赖 make_interp_spline,而该函数严格要求输入值数组完全有限(finite),不接受任何缺失值。
根本原因在于:RegularGridInterpolator 是为“完整规则网格”设计的;它假设每个 (lat_i, lon_j) 都有定义良好的值,从而构建张量积样条基函数。一旦出现 NaN,不仅破坏了网格完整性,更导致沿任一维度做一维样条拟合时失败(如错误栈中 _do_spline_fit 所示)。
✅ 正确解法是切换到无结构(unstructured)插值范式:将原始网格展平为有效点云,再使用 scipy.interpolate.griddata。该函数原生支持 method='cubic'(基于分片双三次 Hermite 插值),且自动忽略无效点、鲁棒处理空洞邻域,并对查询点外推时返回 NaN(符合预期行为)。
以下是推荐实现流程(适配您的地理数据场景):
import numpy as np
from scipy.interpolate import griddata
# 假设 lat0 (shape: M), lon0 (shape: N), data_nan (shape: M×N)
# 其中 data_nan 含约 40% NaN,但所有查询点 (lat, lon) 均在网格范围内
# 1. 构建完整坐标网格并展平
latG, lonG = np.meshgrid(lat0, lon0, indexing='ij') # shape: (M, N)
lat_flat = latG.ravel()
lon_flat = lonG.ravel()
data_flat = data_nan.ravel()
# 2. 筛选非 NaN 有效点(保留坐标与值的一致性)
mask_valid = ~np.isnan(data_flat)
points_valid = np.column_stack((lat_flat[mask_valid], lon_flat[mask_valid]))
values_valid = data_flat[mask_valid]
# 3. 对目标散点 (lat, lon) 执行三次插值
# 注意:griddata 要求查询点为 (N, 2) 形状的数组
xi = np.column_stack((lat, lon)) # shape: (K, 2)
interped_cub = griddata(
points=points_valid,
values=values_valid,
xi=xi,
method='cubic',
fill_value=np.nan # 显式指定越界/空洞处返回 NaN(默认即此行为)
)
# interped_cub.shape == (len(lat),) —— 与 linear 插值结果维度一致!? 关键优势说明:
- ✅ 自动容错:griddata 在 'cubic' 模式下会基于 Delaunay 三角剖分,仅利用局部有效邻点构造插值多项式,天然跳过 NaN 区域;
- ✅ 语义一致:对靠近 NaN 区域的查询点,若缺乏足够邻近有效支撑点,结果自动为 NaN,与 RegularGridInterpolator(method='linear') 的失败行为逻辑统一;
- ✅ 无需后处理:避免了“先线性 + 再最近邻填补”的两阶段繁琐流程;
- ⚠️ 注意性能:griddata 的 'cubic' 比 'linear' 计算开销大,但对数千个查询点(如您的 ~2000)仍非常高效;若需更高性能,可考虑 scipy.interpolate.CloughTocher2DInterpolator(同样支持 fill_value,且缓存三角剖分)。
? 补充建议:
- 若地理跨度大(如跨极区),建议先将经纬度转换为等距投影坐标(如 pyproj 中的 EPSG:3857 或 EPSG:4326 下的球面距离加权),以提升插值几何合理性;
- 对于超大规模网格(>10⁵ 有效点),可启用 rescale=True 参数缓解坐标尺度差异导致的数值不稳定;
- 验证时,可用小范围人工 NaN 区域对比 griddata 与手动掩膜后 RegularGridInterpolator('linear') 结果,确认空洞边缘过渡平滑性。
至此,您即可在保持代码简洁性的同时,获得比线性插值更光滑、比最近邻更物理合理的三次插值结果,且完全兼容真实地球科学数据中常见的稀疏、不规则缺失模式。









