Spline Toolbox

csapi

Cubic spline interpolation

Syntax

values = csapi(x,y,xx)
pp = csapi(x,y)

Description

A cubic spline s with knot sequence x is constructed that satisfies
s(x(j))=y(j) for all j, as well as the not-a-knot end conditions,
jump_x(2) D³s = 0 = jump_x(end-1) D³s (with D³s the third derivative of s).

The call csapi(x,y,xx) returns the values s(xx) of this interpolating cubic spline at the given argument sequence xx.

The alternative call csapi(x,y) returns instead the ppform of the cubic spline, for later use with fnval, fnder, etc.

If x is a cell array, containing sequences x1, ..., xm, of lengths n1, ..., nm respectively, then y is expected to be an array, of size [n1,...,nm] (or of size [d,n1,...,nm] if the interpolant is to be d-vector-valued), and the output will be an m-cubic spline interpolant to such data. Precisely, if there are only two input arguments, then the output will be the ppform of this interpolant, while, if there is a third input argument, xx, then the output will be the values of the interpolant at the sites specified by xx. If xx is a cell array with m sequence entries, then the corresponding m- (or (m+1)-)dimensional array of grid values is returned. Otherwise, xx must be a list of m-vectors and, the corresponding list of values of the interpolant at these sites is returned.

This command is essentially the MATLAB function spline, which, in turn, is a stripped-down version of the Fortran routine CUBSPL in PGS, except that csapi (and now also spline) accepts vector-valued values and can handle gridded data.

Examples

See the demo csapidem for various examples.

Up to rounding errors, and assuming that x has at least four entries, the statement pp = csapi(x,y) should put the same spline into pp as the statements  

n = length(x);
pp = fn2fm(spapi(augknt(x([1 3:(n-2) n]),4),x,y),'pp');

except that the description of the spline obtained the second way will use no break at x(2) and x(n - 1).

Here is a simple bivariate example, a bicubic spline interpolant to the Mexican Hat function being plotted:

x =.0001+[-4:.2:4]; y = -3:.2:3;
[yy,xx] = meshgrid(y,x); r = pi*sqrt(xx.^2+yy.^2); z = sin(r)./r;
bcs = csapi( {x,y}, z ); fnplt( bcs ), axis([-5 5 -5 5 -.5 1])

Note the reversal of x and y in the call to meshgrid, needed since MATLAB likes to think of the entry z(i,j) as the value at (x(j),y(i)) while this toolbox follows the Approximation Theory standard of thinking of z(i,j) as the value at (x(i),y(j)). Similar caution has to be exerted when values of such a bivariate spline are to be plotted with the aid of MATLAB's mesh, as is shown here (note the use of the transpose of the matrix of values obtained from fnval):

xf = linspace(x(1),x(end),41); yf = linspace(y(1),y(end),41);
mesh(xf, yf, fnval( bcs, {xf, yf}).')

Algorithm

The relevant tridiagonal linear system is constructed and solved, using MATLAB's sparse matrix capability.

The not-a-knot end condition is used, thus forcing the first and second polynomial piece of the interpolant to coincide, as well as the second-to-last and the last polynomial piece.

See Also

csape, spapi, spline, tspdem

Cautionary Note

If the sequence x is not nondecreasing, both x and y will be reordered in concert to make it so.

csape

csaps