我已经设法写了一个符合我的目的的函数。它通过沿网格线插值，然后在x和y方向插值平面，并取两者的平均值，从坐标网格中插值(填充)平面。在

通过将坐标重塑为一维矢量，一次性插值平面，然后再重新塑造为二维，应该可以稍微加快这一速度。但是，对于合理的平面尺寸来说，这个代码已经足够快了。在

如果坐标也在平面外，似乎也可以工作。

如果网格近似规则，则外推法也有效。不管怎样，它都会外推，但是随着栅格不规则度的增加，你会看到一些尖锐的折痕远离边缘。在

这是密码。docstring中提供了一个示例。在def interlin2d(x,y,z,fsize):

"""

Linear 2D interpolation of a plane from arbitrary gridded points.

:param x: 2D array of x coordinates

:param y: 2D array of y coordinates

:param z: 2D array of z coordinates

:param fsize: Tuple of x and y dimensions of plane to be interpolated.

:return: 2D array with interpolated plane.

This function works by interpolating lines along the grid point in both dimensions,

then interpolating the plane area in both the x and y directions, and taking the

average of the two. Result looks like a series of approximately curvilinear quadrilaterals.

Note, the structure of the x,y,z coordinate arrays are such that the index of the coordinates

indicates the relative physical position of the point with respect to the plane to be interpoalted.

Plane is allowed to be a subset of the range of grid coordinates provided.

Extrapolation is accounted for, however sharp creases will start to appear

in the extrapolated region as the grid of coordinates becomes increasingly irregular.

Scipy's interpolation function is used for the grid lines as it allows for proper linear extrapolation.

However Numpy's interpolation function is used for the plane itself as it is robust against gridlines

that overlap (divide by zero distance).

Example:

#set up number of grid lines and size of field to interpolate

nlines=[3,3]

fsize=(100,100,100)

#initialize the coordinate arrays

x=np.empty((nlines[0],nlines[1]))

y=np.empty((nlines[0],nlines[1]))

z=np.random.uniform(0.25*fsize[2],0.75*fsize[2],(nlines[0],nlines[1]))

#set random ordered locations for the interior points

spacings=(fsize[0]/(nlines[0]-2),fsize[1]/(nlines[1]-2))

for k in range(0, nlines[0]):

for l in range(0, nlines[1]):

x[k, l] = round(random.uniform(0, 1) * (spacings[0] - 1) + spacings[0] * (k - 1) + 1)

y[k, l] = round(random.uniform(0, 1) * (spacings[1] - 1) + spacings[1] * (l - 1) + 1)

#fix the edge points to the edge

x[0, :] = 0

x[-1, :] = fsize[1]-1

y[:, 0] = 0

y[:, -1] = fsize[0]-1

field = interlin2d(x,y,z,fsize)

"""

from scipy.interpolate import interp1d

import numpy as np

#number of lines in grid in x and y directions

nsegx=x.shape[0]

nsegy=x.shape[1]

#lines along the grid points to be interpolated, x and y directions

#0 indicates own axis, 1 is height (z axis)

intlinesx=np.empty((2,nsegy,fsize[0]))

intlinesy=np.empty((2,nsegx,fsize[1]))

#account for the first and last points being fixed to the edges

intlinesx[0,0,:]=0

intlinesx[0,-1,:]=fsize[1]-1

intlinesy[0,0,:]=0

intlinesy[0,-1,:]=fsize[0]-1

#temp fields for interpolation in x and y directions

tempx=np.empty((fsize[0],fsize[1]))

tempy=np.empty((fsize[0],fsize[1]))

#interpolate grid lines in the x direction

for k in range(nsegy):

interp = interp1d(x[:,k], y[:,k], kind='linear', copy=False, fill_value='extrapolate')

intlinesx[0,k,:] = np.round(interp(range(fsize[0])))

interp = interp1d(x[:, k], z[:, k], kind='linear', copy=False, fill_value='extrapolate')

intlinesx[1, k, :] = interp(range(fsize[0]))

intlinesx[0,:,:].sort(0)

# interpolate grid lines in the y direction

for k in range(nsegx):

interp = interp1d(y[k, :], x[k, :], kind='linear', copy=False, fill_value='extrapolate')

intlinesy[0, k, :] = np.round(interp(range(fsize[1])))

interp = interp1d(y[k, :], z[k, :], kind='linear', copy=False, fill_value='extrapolate')

intlinesy[1, k, :] = interp(range(fsize[1]))

intlinesy[0,:,:].sort(0)

#interpolate plane in x direction

for k in range(fsize[1]):

tempx[k, :] = np.interp(range(fsize[1]),intlinesx[0,:,k], intlinesx[1,:,k])

#interpolate plane in y direction

for k in range(fsize[1]):

tempy[:, k] = np.interp(range(fsize[0]), intlinesy[0, :, k], intlinesy[1, :, k])

return (tempx+tempy)/2