python三维插值从规则网格到散点,将大型不规则网格插值到Python中的另一个不规则网格...

I am trying to interpolate complex values from one irregular grid to another irregular grid using Python. The grids are in 2D and there are 103,113 data points. I am using Python 2.6.6, Scipy 0.7.2, Numpy 1.3.0, Matplotlib 0.99.3

In Matlab using griddata this is achieved in roughly 5 seconds.

BnGRID2 = griddata(R_GRID1,Z_GRID1,BnGRID1,R_GRID2,Z_GRID2) (MATLAB)

(Note all arrays are 201 x 513)

However, if I try using matplotlib.mlab.griddata I get a memoryError even if I try to work with the real part only:

mlab.griddata(R_GRID1.flatten(),Z_GRID1.flatten(),num.real(BnGRID1.flatten()),R_GRID2.flatten(),Z_GRID2.flatten())

If I try using interp2d I get a segmentation fault and Python exits:

a = interp.interp2d(R_GRID1,Z_GRID1,num.real(BnGRID1))

I have tried using KDTree and this seems to work ok, however, it takes a few minutes compared with the few seconds for Matlab, but I haven't explored this option too much yet.

Was wondering if anyone has any ideas how I can get this done as quickly as Matlab seems to? I noticed that the newer version of Scipy also has griddata, does anyone know if this can handle large irregular grids?

解决方案

Scipy's griddata seems to be able to deal with data sets of this size without problems:

import numpy as np

import scipy.interpolate

# old grid

x, y = np.mgrid[0:1:201j, 0:1:513j]

z = np.sin(x*20) * (1j + np.cos(y*3))**2 # some data

# new grid

x2, y2 = np.mgrid[0.1:0.9:201j, 0.1:0.9:513j]

# interpolate onto the new grid

z2 = scipy.interpolate.griddata((x.ravel(), y.ravel()), z.ravel(), (x2, y2), method='cubic')

The griddata step takes about 5s on an old AMD Athlon.

If your data is on a grid (i.e., the coordinates corresponding to value z[i,j] are (x[i], y[j])), you can get more speed by using scipy.interpolate.RectBivariateSpline

z3 = (scipy.interpolate.RectBivariateSpline(x[:,0], y[0,:], z.real)(x2[:,0], y2[0,:])

+ 1j*scipy.interpolate.RectBivariateSpline(x[:,0], y[0,:], z.imag)(x2[:,0], y2[0,:]))

which takes 0.05s. It's much faster, because even if your grid spacings are irregular, a more efficient algorithm can be used as long as the grid is rectangular.