The task at hand is seemingly simple: I have a 2D grid of data. The data is available in 2D arrays for X and Y coordinates, as well as the input variable which I want to interpolate. This means I can plot the data using rectangular cells, which means it is possible to use bilinear interpolation. Unfortunately, the data is not precisely aligned with the coordinates, and also not precisely spaced. There were some numerics involved in creating the data, which means that all sampling locations are a little off the mark, and the cell spacing is smooth but not uniform.
I would like to interpolate from this input grid to a set of predefined sample coordinates (as opposed to simply refining the mesh).
In short, an example for my type of input is:
# a nice, regular grid
Xs, Ys = np.meshgrid(np.linspace(0, 1, num=3), np.linspace(0, 1, num=5))
# ...perturbed by some systematic and some random noise ...
X_in = Xs + np.random.normal(scale=0.03, size=(5, 3))
# ...and some systematic deviation
Y_in = (Ys + np.random.normal(scale=0.03, size=(5, 3)))* (1 + Xs**1.5)
# and some variable at each node to interpolate
Z_in = np.random.normal(scale=1, size=(5, 3))
So (X_in, Y_in) are arrays of shape (n, m) which define a mesh with quadrilateral cells, and Z_in another array of ther same shape which provides a value at each node in that mesh. I am looking for some Python library that performs bilinear interpolation of Z_in across those cells.
However, all methods I have found so far either ignore the rectangular structure (and triangulate the data, or fit some 2D spline through arbitrary point clouds), or require a perfectly rectangular and equally-spaced grid as input (which mine is not).
Examples of answers/methods that seem not to be appliccable:
This answer recommends using scipy.ndimage.map_coordinates -- but that effectively uses the indices of the 2D input data array as coordinates, which won't work for me.
scipy.interpolate.interp2d requires either a regular grid (node locations provided by 1D X and Y arrays), or an irregular one, which is flattened, which means that the algorithm cannot know which nodes form a cell. This means it either fits some spline through unstructured data, or triangulates it. And it only interpolates onto regular grids or individual points.
scipy.interpolate.RectBivariateSpline is recommended for interpolation from gridded data but only accepts input points which are perfectly aligned with the coordinate system.
There's also a Matplotlib toolkit for interpolation, which I had thought should be able to do this sort of thing, as it also does interpolated contour plots of rectangular meshes, but as it turns out, even though mpl_toolkits.basemap.interp accepts arbitrary quadrilateral meshes as target for interpolation, it cannot use them as inputs ...
Upon closer inspection, it turns out that even matplotlib.plt.contour() does not seem to perform bilinear interpolation when plotting the input data:
plt.contour(X_in, Y_in, Z_in, levels=np.linspace(Z_in.min(), Z_in.max(), 50))
plt.plot(X_in, Y_in, 'k-')
plt.plot(X_in.T, Y_in.T, 'k-')
As you can see, the contour lines within the cells are straight, but with bilinear interpolation, they should not be, and there should not be those empty quadrilateral areas in the mittle of some cells. I suspect that Matplotlib only finds the contour values on the cell edges and simply draws straight lines between them.
I have found two explanations of the maths of bilinear interpolation from grids which are not perfectly aligned, but I was hoping to come across a ready-made implementation somewhere because I'm sure that this kind of task is not so rare, and a numpy or scipy implementation (if it exists) is probably way faster than whatever I'd implement myself.
