Wow. I guess I've been working with those tools for so many years, as well as building variations of them, that I think of them as second nature.
What is often called bilinear interpolation, when done on a lattice is what interp2 or griddedInterpolant calls linear. And, I guess it may not be obvious that what one tool calls linear, is what someone else might call bilinear (for example, photoshop, as I recall has a bilinear interpolation tool for image interpolation) or someone else might call tensor product linear interpolation. But the formula you have written is essentially that tensor product linear interpolant. So in any dimension, the result will be a piecewise linear function.
I'll show how to use griddedInterpolant, which I prefer myself. This allows you to learn to use one tool, in various numbers of dimensions.
In the help for griddedInterpolant, we see:
F = griddedInterpolant(...
'linear' - (default) linear, bilinear, trilinear,...
So when we ask griddedInterpolant for 'linear' interpolation in multiple dimenssions, it does a bilinear or trilinear interpolation, as you wish.
x = 0:0.5:1;
y = linspace(0,1,4);
[xx,yy] = ndgrid(x,y)
xx =
0 0 0 0
0.5 0.5 0.5 0.5
1 1 1 1
yy =
0 0.33333 0.66667 1
0 0.33333 0.66667 1
0 0.33333 0.66667 1
Now create an array z, of the form z(x,y). So two independent variables, and thus a bilinear interpolant.
zz = xx.^2 + yy.^2
zz =
0 0.11111 0.44444 1
0.25 0.36111 0.69444 1.25
1 1.1111 1.4444 2
z(x,y) is pretty boring, I'll admit. Sorry about that. It is a paraboloid of revolution, in two dimensions.
GI = griddedInterpolant(xx,yy,zz,'linear')
GI =
griddedInterpolant with properties:
GridVectors: {[0 0.5 1] [0 0.33333 0.66667 1]}
Values: [3×4 double]
Method: 'linear'
ExtrapolationMethod: 'linear'
Can I convince you that GI is a bilienar interpolant? For example, fix x at some value.
yint = linspace(0,1,100);
xint = reshape(0.2,size(yint));
plot(yint,GI(xint,yint),'-r')
xlabel 'yint'
ylabel 'z(xint,yint)'
grid on
So 3 linear segments. This is a tensor product linear (bilinear) interpolant, as you asked about.
With 3 independent variables, so trilinear interpolation, you would create xx,yy,zz using ndgrid. Then create a 3-d array of values that defines the relation w(x,y,z). Now griddedInterpolant will do the trilinear interpolation you want.
If you prefer to use interp2 interp3, or interpn, you can also do so, but these tools require arrays that are ordered as meshgrid would create, instead of ndgrid.
Be careful though. One of the things people need to understand about these interpolants is they are not truly linear when viewed along any path through space. For example,
[xx,yy] = ndgrid([0 1],[0 1])
xx =
0 0
1 1
yy =
0 1
0 1
zz = [0 .5 ; .5 2]
zz =
0 0.5
0.5 2
Now, consider the path that moves acrosss the diagonal of this cell, running from the points (0,0) to (1,1), with x==y along that path.
GI = griddedInterpolant(xx,yy,zz,'linear');
xint = linspace(0,1,100);
yint = xint;
zint = GI(xint,yint);
plot(xint,GI(xint,yint),'-r')
xlabel 'xint == yint'
ylabel 'z(xint,yint)'
grid on
So even though this is a bilinear interpolant, when you vary two variables at once, the shape of the function can exhibit a quadratic polynomial shape.
In the case of a triinear interpolant, you will see cubic polynomials when traversed along a diagonal. This is completely as expected. Except of course, it seems to come as a surprise to many, probably because this is a "linear" interpolation. I will admit the first time I saw this behavior, it seemed surprising to me. Once I learned how to manipulate these tools, my next trick was to build a 3-d array, such that when interpolated along the diagonal using trilinear interpolation had the complete behavior of a user supplied cubic spline. Totally useless of course, unless you wanted exactly that behavior, which sometimes, we did. :)
Best Answer