插值后向 matplotlib 图添加特征叠加的问题

Problem adding features overlay to matplotlib plot after interpolation

我不得不承认,我仍然无法理解绘图的正确设置和关系以及它与 matplotlib 的部分,仍然混淆 fig 与 plt 与 ax 的相互关系,所以我只是反复试验,文档有时让我更困惑。 :-(

我正在根据 json 绘制天气值并获得积分。我可以用下面的代码绘制如下图

fig=plt.figure(figsize=(10,8))
ax=fig.add_subplot(1,1,1,projection=mapcrs)
ax.set_extent([-93,-86,13,19],datacrs)
ax.add_feature(cfeature.COASTLINE)
ax.add_feature(cfeature.BORDERS, linestyle=':')
ax.scatter(lon,lat,c=dat,transform=datacrs)

而且我可以绘制地图

然后我用这段代码使用 metpy 生成插值

gridx, gridy, gridz = interpolate_to_grid(lon, lat, dat, interp_type='rbf', hres=.1, rbf_func='linear', rbf_smooth=0)

fig=plt.figure(figsize=(15,15))
ax=fig.add_subplot(1,1,1,projection=mapcrs)
#ax = fig.add_axes([0, 0, 1, 1], projection=mapcrs)
#ax.set_extent([-93,-86,13,19])
#ax.add_feature(cfeature.COASTLINE)
#ax.add_feature(cfeature.BORDERS, linestyle=':')
ax.contourf(gridx,gridy,gridz,levels=np.arange(10,60,2),cmap='viridis')
plt.plot(lon,lat,'k.',color='white')

按要求插值得到点,但无法显示特征,请问如何实现?如果我取消对 ax.extent 的注释,我看到的只是一个空白的白色图形。如果我取消注释 ax.features,插值显示为下图而不是地图。

感谢您的帮助和指导。

您在 contourf 函数中缺少 transform 关键字参数以给出插值数据的坐标系。这是一个使用随机数据的最小工作示例,获得的输出如下:

import numpy as np

from cartopy import crs, feature
from matplotlib import pyplot as plt
from scipy.interpolate import griddata

# figure
fig = plt.figure(figsize=(5, 5))

# coordinate systems
crs_map = crs.Mercator()
crs_data = crs.PlateCarree()

# random data
np.random.seed(42)  # for repro.
n = 100
lon = -89 + 2 * np.random.randn(n)
lat = 16 + 2 * np.random.randn(n)
dat = np.random.rand(n)

# interpolated data
ilon = np.linspace(-93, -86, 200)
ilat = np.linspace(13, 19, 200)
ilon, ilat = np.meshgrid(ilon, ilat)
idat = griddata((lon, lat), dat, (ilon, ilat), method="linear")

# show up
ax = fig.add_subplot(1, 1, 1, projection=crs_map)
ax.set_extent([-93, -86, 13, 19], crs_data)
ax.add_feature(feature.COASTLINE)
ax.add_feature(feature.BORDERS, ls=":", lw=0.5)
ax.scatter(lon, lat, c=dat, transform=crs_data)  # this is invisible with contour
ax.plot(lon, lat, "k.", transform=crs_data)  # in order to see the points
ax.contourf(ilon, ilat, idat, levels=np.linspace(0, 1, 10), transform=crs_data)