Cartopy绘制二维经纬度坐标下的RFnet地理数据

针对非规则二维经纬度网格的RFnet数据，可直接通过Cartopy结合matplotlib原生支持二维坐标的绘图接口实现地理空间可视化，不需要提前对数据做插值重格。

完整实现代码 #

import xarray as xr
import numpy as np
import matplotlib.pyplot as plt
import cartopy.crs as ccrs
import cartopy.feature as cfeature

# 读取NetCDF数据集
ds = xr.open_dataset('R90C_20P_MSG_201403250100.nc', engine='netcdf4')

# 提取目标变量为numpy数组
lon = ds['LONGITUDE'].data
lat = ds['LATITUDE'].data
rfnet = ds['RFnet'].data

# 掩膜无效像元：过滤经纬度超出合理范围的卫星边缘填充值
rfnet = np.where(
    (lon >= -180) & (lon <= 180) & (lat >= -90) & (lat <= 90),
    rfnet,
    np.nan
)

# 创建画布，指定绘图投影（可替换为其他Cartopy支持的投影）
fig = plt.figure(figsize=(12, 6))
ax = fig.add_subplot(111, projection=ccrs.PlateCarree())

# 绘制填色图：传入二维经纬度，指定坐标参考系
im = ax.pcolormesh(
    lon, lat, rfnet,
    cmap='seismic',
    vmin=-10,
    vmax=10,
    transform=ccrs.PlateCarree(),
    shading='auto'
)

# 添加地理辅助要素
ax.add_feature(cfeature.COASTLINE.with_scale('10m'), linewidth=0.8)
ax.add_feature(cfeature.BORDERS.with_scale('10m'), linewidth=0.5)

# 自动适配数据范围设置绘图边界，也可手动传入[minlon, maxlon, minlat, maxlat]自定义范围
ax.set_extent(
    [np.nanmin(lon), np.nanmax(lon), np.nanmin(lat), np.nanmax(lat)],
    crs=ccrs.PlateCarree()
)

# 添加经纬度网格
grid = ax.gridlines(draw_labels=True, linestyle='--', alpha=0.5)
grid.top_labels = False
grid.right_labels = False

# 添加色标
cbar = fig.colorbar(im, location='bottom', shrink=0.6, pad=0.1, ax=ax)
cbar.set_label('RFnet')

plt.tight_layout()
plt.show()

关键注意事项 #

• 不要使用imshow绘图：imshow仅支持等间距规则矩形网格，静止卫星原始产品的经纬度为随像元变化的二维弯曲网格，必须使用支持二维坐标输入的pcolormesh或contourf接口。
• transform参数不可省略：该参数用于告知Cartopy传入的经纬度为WGS84地理坐标系下的坐标，Cartopy会自动完成坐标到目标绘图投影的转换，无需手动处理投影变换。
• 维度匹配要求：传入pcolormesh的经度、纬度、数据数组三者维度必须完全一致，本案例中三个数组均为(image_y, image_x)形状，可直接传入，若出现维度报错可检查数组是否被意外转置。
• 无效值处理：卫星扫描边缘的填充像元通常存在异常经纬度值，提前掩膜为NaN可避免绘图出现跨区域的异常条纹。
• 投影自定义：如果需要更换绘图投影（如兰伯特共形投影、墨卡托投影等），仅需修改add_subplot时的projection参数即可，其余代码无需调整。

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