为什么我要插值？因为我的日值场和气候态场的网格间隔不一样（前者是0.75°×0.75°，后者是0.25°×0.25°），没办法相减，所以想把日值场插值到气候态场。
如何插值?参考了钦天监摸鱼记事的文章：https://www.jianshu.com/p/ff0888f7add8
遇到的问题？插值前后的图对不上。
解决过程： 画了一下全球的图（我原来只画了我关注的区域），发现是插值后变量的南北纬反了。

反思： 处理数据、计算、画图的代码都找不出问题，可以试试画出全球的图找问题。某女士画500hPa高度图画不出台风也是因为这个（弱智）问题。

代码

插值代码：

# -*- coding: utf-8 -*-
"""
Created on Sun Apr  3 09:20:23 2022

@author: Onefire Jung
"""

import xarray as xr 
import numpy as np 
from scipy.interpolate import griddata #插值函数

import matplotlib.pyplot as plt 
import cartopy.crs as ccrs                                               
import cartopy.feature as cf                                             
from cartopy.mpl.ticker import LongitudeFormatter,LatitudeFormatter
import cartopy.mpl.ticker as cticker
from matplotlib.ticker import MaxNLocator

#c场为气候场；d场为日值场
filec = xr.open_dataset(r'.\part2\hgt1981_2010_whole.nc')
filed = xr.open_dataset(r'.\part2\hgt_2012.nc')

#挑出需要的变量数据：气候场中跳出xxxhPa做时间平均，日值场中挑
hgt_July = np.average(filec.z[:,1,:,:],axis = 0)
hgt_d = filed.z.loc["2012-07-21"]

#c场的经纬度（0.25×0.25）
lon = filec.longitude
lat = filec.latitude

#d场的原经、纬度（0.75×0.75）
lons = filed.longitude.values 
lats = filed.latitude.values


############           插值            ##########

#构造point（被插场的原网格）
x = np.tile(lons,lats.shape[0]) #重复完整的lons n次
y = np.repeat(lats,lons.shape[0]) #重复lats中每个元素m次
point = np.stack([x,y],axis = -1) #拼接

#构造points（目标网格）：官方建议插值使用numpy的mgrid导入网格
lon1,lon2,lat1,lat2 = np.min(lon),np.max(lon),np.min(lat),np.max(lat)
step = 0.25
lat_new1, lon_new1 = np.mgrid[-90.0:90.0+0.25:0.25,0.0:359.75+0.25:0.25]
lon_new= lon_new1.ravel() #展成一维
lat_new= lat_new1.ravel() #展成一维
points= np.stack([lon_new,lat_new],axis = -1)

#把values转成一维（插值函数要求values是一维的）
hgt_d_1dimon = hgt_d.values.ravel()

#插值函数，插完得到的是一维数组，还需要转成二维
result =griddata(point,hgt_d_1dimon,points,method='cubic',fill_value=0,rescale=True) 
hgt_d_resize = result.reshape((721, 1440)) #转成二维

#单位华为位势什米
hgt_July = hgt_July/100.
hgt_d_resize = hgt_d_resize /100.
#异常场
hgt_a = hgt_d_resize - hgt_July

'''
7月的气候场为：hgt_July(全域）
721当天：hgt_d_resize（全域）

'''


############           画图           ##########
plt.rcParams['font.sans-serif']=['SimHei'] #用来正常显示中文标签

fig = plt.figure(dpi=200)
ax1 = fig.add_subplot(1,1,1,projection=ccrs.PlateCarree())

# x,y是横纵坐标的索引，左闭右开
x1,x2 = 240,601
y1,y2 = 160,321 
leftlon, rightlon, lowerlat, upperlat = (60,150,10,50)
#mycrs = ccrs.PlateCarree(central_longitude=180)
mycrs = ccrs.PlateCarree()
ax1.set_global() #使得轴域（Axes即两条坐标轴围城的区域）适应地图的大小
ax1.coastlines() #画出海岸线
ax1.set_extent([leftlon, rightlon, lowerlat, upperlat], crs=mycrs)
ax1.set_xticks(np.arange(leftlon,rightlon+10,10), crs=mycrs)
ax1.set_yticks(np.arange(lowerlat,upperlat+10,10), crs= mycrs)
# # #使横纵向坐标含WESN
lon_formatter = cticker.LongitudeFormatter()
lat_formatter = cticker.LatitudeFormatter()
ax1.xaxis.set_major_formatter(lon_formatter)
ax1.yaxis.set_major_formatter(lat_formatter)

ax1.set_title('(a)插值过的day',loc='left',fontsize = 14)
# C = plt.contour(lon,lat,hgt_d_resize/100.,colors='black',transform=ccrs.PlateCarree()) #生成等值线图
# plt.clabel(C,inline=2,fmt='%1.0f',fontsize=10)
# plt.contourf(lon[240:601],lat[0:321],hgt_d_resize[0:321,240:601]., cmap='RdBu',levels = np.arange(520,650,10),transform=mycrs)
plt.contourf(lon[240:601],lat[160:321],hgt_d_resize[160:321,240:601], cmap='RdBu',levels = 20,transform=mycrs)
plt.colorbar(shrink=1,orientation = 'horizontal',pad = 0.1)
plt.show()

画原图的代码

···# -*- coding: utf-8 -*-
"""
Created on Sun Apr  3 12:06:20 2022

@author: Onefire Jung
"""
import xarray as xr 
import numpy as np 
from scipy.interpolate import griddata #插值函数

import matplotlib.pyplot as plt 
import cartopy.crs as ccrs                                               
import cartopy.feature as cf                                             
from cartopy.mpl.ticker import LongitudeFormatter,LatitudeFormatter
import cartopy.mpl.ticker as cticker
from matplotlib.ticker import MaxNLocator

#c场为气候场；d场为日值场
filed = xr.open_dataset(r'.\part2\hgt_2012.nc')
lon = filed.longitude
lat = filed.latitude

#挑出需要的变量数据：气候场中跳出xxxhPa做时间平均，日值场中挑

hgt_d = filed.z.loc["2012-07-21"].values/100.0
plt.rcParams['font.sans-serif']=['SimHei'] #用来正常显示中文标签

fig = plt.figure(dpi=200)
ax1 = fig.add_subplot(1,1,1,projection=ccrs.PlateCarree())

# x,y是横纵坐标的索引，左闭右开
leftlon, rightlon, lowerlat, upperlat = (60,150,10,50)
#mycrs = ccrs.PlateCarree(central_longitude=180)
mycrs = ccrs.PlateCarree()
ax1.set_global() #使得轴域（Axes即两条坐标轴围城的区域）适应地图的大小
ax1.coastlines() #画出海岸线
ax1.set_extent([leftlon, rightlon, lowerlat, upperlat], crs=mycrs)
ax1.set_xticks(np.arange(leftlon,rightlon+10,10), crs=mycrs)
ax1.set_yticks(np.arange(lowerlat,upperlat+10,10), crs= mycrs)
# # #使横纵向坐标含WESN
lon_formatter = cticker.LongitudeFormatter()
lat_formatter = cticker.LatitudeFormatter()
ax1.xaxis.set_major_formatter(lon_formatter)
ax1.yaxis.set_major_formatter(lat_formatter)

ax1.set_title('(a)没插值的',loc='left',fontsize = 14)
# C = plt.contour(lon,lat,hgt_d_resize/100.,colors='black',transform=ccrs.PlateCarree()) #生成等值线图
# plt.clabel(C,inline=2,fmt='%1.0f',fontsize=10)
# plt.contourf(lon[240:601],lat[0:321],hgt_d_resize[0:321,240:601]., cmap='RdBu',levels = np.arange(520,650,10),transform=mycrs)
plt.contourf(lon[81:201],lat[41:108],hgt_d[41:108,81:201], cmap='RdBu',levels = 20,transform=mycrs)
plt.colorbar(shrink=1,orientation = 'horizontal',pad = 0.1)

plt.show()

图片对比

插值后的图.png

没插值的图.png

hgt_d_resize = hgt_d_resize[::-1,:]

对了.png