Python三维数组的二维插值

Question

我有一个维度为（时间=410，纬度=72，经度=144）的全球降水数据。我希望将lat插值到180，将lon插值到360（即，将数据从2.5度缩小到1度）。

在Matlab中，我曾经做过：

%LAT,LON,XI,YI are old and new meshgridded lat-lons
for t = 1:size(precip,1)
    newPrecip(t,:,:)=interp2(LON,LAT,squeeze(precip(t,:,:)),XI,YI);
end

在python中，我尝试过interp2d、map_coordinates和interp（basemap）函数，但都没有成功。

这是我的map_coordinates代码，它使我最接近于获得结果

new_indicies = np.mgrid[0:410, -89.5:89.5:180J, -179.5:179.5:360J]
newPrecip = ndimage.map_coordinates(precip, new_indicies, order=1)

新数据确实具有我想要的维度（410180360），但它实际上并没有插值，而是用0值填充新添加的网格（围绕原始数据）。

我是python的新手，如果能为我找到这个问题的解决方案，我将不胜感激。

Answer 1

我通常使用Basemap附带的插值例程。

import mpl_toolkits.basemap as mp
Zg = mp.interp(dataIN,lonIN,latIN,lonOUT,latOUT,
               checkbounds=False, masked=False, order=1)

这里，lonIN和latIN是原始网格的1D经度-纬度（假设为常规网格），dataIn是原始数据的2D数组，lonOUT和latOUT是要插值到的2D网格，而Zg是新网格上的输出。可以使用创建输出网格

x=np.arange(-180,180,1)
y=np.arange(-90,90,1)
lonOUT,latOUT=np.meshgrid(x,y)

希望这能有所帮助。干杯，Trond

Answer 2

我会研究scipy.interpolate。看起来你可能需要网格数据。

Answer 3

多亏了Trond，我意识到自己做错了什么。这是适用于我的最后一个代码。希望它能为其他人提供参考。这是我第一次在Stackoverflow上发帖，很高兴我的查询得到了快速正确的回答！

#this code assumes a input data 'precip' of dimensions (410,72,144)
import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.basemap import interp
LAT=np.arange(-89.5,90.5,1)
LON=np.arange(-179.5,180.5,1)
LON,LAT=np.meshgrid(LON,LAT)
lat=np.arange(-88.75,91.25,2.5)
lon=np.arange(-178.75,181.25,2.5)
newPrecip=np.zeros((410,180,360), dtype='float')
for i in range(410):
    newPrecip[i,:,:]=interp(np.squeeze(precip[i,:,:]),lon,lat,LON,LAT,order=1)
plt.figure(1)
plt.pcolor(lon,lat,precip.mean(axis=0))
plt.figure(2)
plt.pcolor(LON,LAT,newPrecip.mean(axis=0))

Answer 4

如果您想要普通克里格插值，另一个选项是PyKrige库！