matplotlib.mlab.griddata非常慢，当有效数据为inpu时返回nan数组

Question

我试图将一个不规则网格数据集（原始卫星数据）与相关的纬度和经度映射到由basemap.makegrid()给出的一组规则网格化的纬度和经度集合。我正在使用安装了mpl_toolkits.natgrid的matplotlib.mlab.griddata。下面是ipython中whos用作输出的变量列表以及变量的一些统计信息：

Variable   Type       Data/Info
-------------------------------
datalat    ndarray    666x1081: 719946 elems, type `float32`, 2879784 bytes (2 Mb)
datalon    ndarray    666x1081: 719946 elems, type `float32`, 2879784 bytes (2 Mb)
gridlat    ndarray    1200x1000: 1200000 elems, type `float64`, 9600000 bytes (9 Mb)
gridlon    ndarray    1200x1000: 1200000 elems, type `float64`, 9600000 bytes (9 Mb)
var        ndarray    666x1081: 719946 elems, type `float32`, 2879784 bytes (2 Mb)

In [11]: var.min()
Out[11]: -30.0

In [12]: var.max()
Out[12]: 30.0

In [13]: datalat.min()
Out[13]: 27.339874

In [14]: datalat.max()
Out[14]: 47.05302

In [15]: datalon.min()
Out[15]: -137.55658

In [16]: datalon.max()
Out[16]: -108.41629

In [17]: gridlat.min()
Out[17]: 30.394031556984299

In [18]: gridlat.max()
Out[18]: 44.237140350357713

In [19]: gridlon.min()
Out[19]: -136.17646180595321

In [20]: gridlon.max()
Out[20]: -113.82353819404671

datalat和{}是原始数据坐标

gridlat和gridlon是要插值到的坐标

{cd9>包含实际数据^

使用这些变量，当我调用griddata(datalon, datalat, var, gridlon, gridlat)时，它需要20分钟才能完成并返回一个nan数组。从数据来看，纬度和经度似乎是正确的，原始坐标重叠了新区的一部分，而一些数据点位于新区之外。有人有什么建议吗？nan值表明我在做蠢事。。。在

Answer 1

更可能的是，griddata太难了。它是用来处理随机采样数据的。您的数据几乎肯定是定期采样的——只是与目标输出网格不在同一网格上。在

如果地球的拓扑或曲率会影响你的结果，那么看看一个更简单的方法，比如仿射变换或者一系列小芯片上的仿射变换。在

有一些现成的解决方案可能会有所帮助。GDAL就是一个很好的例子。在

此外，这类问题经常在GIS中讨论。参见：

https://gis.stackexchange.com/questions/10430/changing-image-projection-using-python

Answer 2

如果您的数据在网格上，使得点(datalon[i], datalat[j])的数据点在data[i,j]中，那么可以使用scipy.interpolate.RectBivariateSpline，而不是{}。不过，一些特定于地理位置的库可能提供更多的功能。在

Answer 3

看起来mlab.griddata例程可能会在输出数据上引入不必要的附加约束。虽然输入位置可以是任何位置，但是输出位置must be a regular grid-因为您的示例是在lat/lon空间中，所以您选择的地图投影可能会违反这一点（即x/y中的规则网格不是lat/lon中的规则网格）。在

您可以尝试来自SciPy的interpolate.griddata例程作为替代-但是您需要将您的位置变量组合到一个数组中，因为调用签名是不同的：类似于

import scipy.interpolate
data_locations = np.vstack(datalon.ravel(), datalat.ravel()).T
grid_locations = np.vstack(gridlon.ravel(), gridlat.ravel()).T
grid_data      = scipy.interpolate.griddata(data_locations, val.ravel(),
                                            grid_locations, method='nearest')

对于最近邻插值。这会将位置放入一个数组中，数组中有2列对应于您的2个维度。还可能需要在贴图投影的变换空间中执行插值。在