Bilinear Interpolation for Data on a Rectangular grid

This is an implementation of a bilinear interpolating function.

For a point (x0,y0) contained in a rectangle (x1,y1),(x2,y1), (x2,y2),(x1,y2) and x1<x2, y1<y2, the first step is to get z() at locations (x0,y1) and (x0,y2) as convex linear combinations z(x0,y*)=a*z(x1,y*)+(1-a)*z(x2,y*) where a=(x2-x1)/(x0-x1) for y*=y1,y2. In a second step z(x0,y0) is calculated as convex linear combination between z(x0,y1) and z(x0,y2) as z(x0,y1)=b*z(x0,y1)+(1-b)*z(x0,y2) where b=(y2-y1)/(y0-y1).

Finally, z(x0,y0) is a convex linear combination of the z values at the corners of the containing rectangle with weights according to the distance from (x0,y0) to these corners.

The grid lines can be unevenly spaced.

Usage

bilinear(x, y, z, x0, y0)
BiLinear(x, y, z, x0, y0)

Arguments

x

a vector containing the x coordinates of the rectangular data grid.

y

a vector containing the y coordinates of the rectangular data grid.

z

a matrix containing the z[i,j] data values for the grid points (x[i],y[j]).

x0

vector of x coordinates used to interpolate at.

y0

vector of y coordinates used to interpolate at.

Note

This Fortran function was part of the akima package but not related to any of Akimas algorithms and under GPL. So it could be transfered into the interp package without changes.

BiLinear is a C++ reimplementation, maybe it will replace the Fortran implementation later, so its name may change in futrure versions.

Value

This function produces a list of interpolated points:

x

vector of x coordinates.

y

vector of y coordinates.

z

vector of interpolated data z.

If you need an output grid, see bilinear.grid.

References

Pascal Getreuer, Linear Methods for Image Interpolation, Image Processing On Line, 2011, http://www.ipol.im/pub/art/2011/g_lmii/article.pdf

Note

Use interpp for the general case of irregular gridded data!

See also

interp, bilinear.grid

Examples

data(akima474)
# interpolate at the diagonal of the grid [0,8]x[0,10]
akima.bil <- bilinear(akima474$x,akima474$y,akima474$z,
                     seq(0,8,length=50), seq(0,10,length=50))
plot(sqrt(akima.bil$x^2+akima.bil$y^2), akima.bil$z, type="l")