xarray插值CMIP6气候数据后绘图出现本初子午线及极点空白问题

问题原因 #

这个空白缺口是两个原因共同导致的：

1. 经度边界不闭合：CMIP6不同模式的经度范围存在0_360°和-180180°两种格式，插值后360°（对应本初子午线）位置的数值没有和0°位置的数值衔接，绘图时就会出现纵向空白线。
2. 极点无数据：大部分模式的格点纬度最高只到±89°左右，没有覆盖南北纬90°的极点位置，绘图时极点区域就会出现空白。

解决方法 #

你可以通过补全循环经度点、补充极点数值两步解决，修改后的代码如下：

import xarray as xr
import matplotlib.pyplot as plt
import cartopy.crs as ccrs
from cartopy.util import add_cyclic_point

cesm2_waccm = xr.open_dataset('pr_day_CESM2-WACCM_ssp245_r2i1p1f1_gn_20750101-20841231.nc')
mpi = xr.open_dataset('pr_day_MPI-ESM1-2-LR_ssp245_r1i1p1f1_gn_20750101-20941231.nc')

cesm2_waccm_subset = cesm2_waccm.sel(time=slice('2075-01-01', '2075-12-31')).mean(dim='time')
mpi_subset = mpi.sel(time=slice('2075-01-01', '2075-12-31')).mean(dim='time')

map_proj = ccrs.PlateCarree()

# 插值操作
mpi_interp = mpi_subset.interp(lat=cesm2_waccm_subset['lat'], lon=cesm2_waccm_subset['lon'])

# 第一步：补全循环经度点，消除本初子午线空白
pr_cyclic, lon_cyclic = add_cyclic_point(mpi_interp.pr, coord=mpi_interp.lon)
# 生成补全后的数据集
interp_cyclic = xr.DataArray(pr_cyclic, dims=['lat', 'lon'], coords={'lat': mpi_interp.lat, 'lon': lon_cyclic}).to_dataset(name='pr')

# 第二步：可选，补充极点数据消除极点空白
# 北极点（北纬90°）赋值为北纬最高纬度圈的平均值
north_pole_val = interp_cyclic.pr.isel(lat=-1).mean()
# 南极点（南纬90°）赋值为南纬最低纬度圈的平均值
south_pole_val = interp_cyclic.pr.isel(lat=0).mean()
# 补充极点坐标
interp_with_pole = interp_cyclic.reindex(lat=[-90] + interp_cyclic.lat.values.tolist() + [90], fill_value=0)
interp_with_pole.pr.loc[dict(lat=90)] = north_pole_val
interp_with_pole.pr.loc[dict(lat=-90)] = south_pole_val

# 绘图就不会有空白了
plot = interp_with_pole.pr.plot(subplot_kws={'projection': map_proj})
plot.axes.coastlines()
plt.show()

补充说明 #

如果你不需要保留插值后的数据集仅用于绘图，也可以直接在xarray的plot方法中添加参数add_cyclic=True直接补全循环点，代码更简洁：

plot = mpi_interp.pr.plot(subplot_kws={'projection': map_proj}, add_cyclic=True)

这个方法可以直接解决本初子午线的空白问题，如果你不需要展示极地范围，不需要额外处理极点数据。

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