Interpolate large irregular grid onto another irregular grid in 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  = g开发者_如何学Goriddata(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.