x=linspace(0,1,11); y=sin(x); xx=linspace(0,1,101); yy = spline(x,y,xx); plot(x,y,'*',xx,yy) pp = spline(x,y); pp pp = form: 'pp' breaks: [1x11 double] coefs: [10x4 double] pieces: 10 order: 4 dim: 1 pp.coefs ans = -0.1652 -0.0003 1.0000 0 -0.1652 -0.0499 0.9950 0.0998 -0.1615 -0.0994 0.9801 0.1987 -0.1567 -0.1479 0.9553 0.2955 -0.1501 -0.1949 0.9211 0.3894 -0.1422 -0.2399 0.8776 0.4794 -0.1325 -0.2826 0.8253 0.5646 -0.1227 -0.3223 0.7648 0.6442 -0.1073 -0.3592 0.6967 0.7174 -0.1073 -0.3913 0.6216 0.7833 pp.breaks ans = Columns 1 through 5 0 0.1000 0.2000 0.3000 0.4000 Columns 6 through 10 0.5000 0.6000 0.7000 0.8000 0.9000 Column 11 1.0000 ppval(pp,xx) ans = Columns 1 through 5 0 0.0100 0.0200 0.0300 0.0400 Columns 6 through 10 0.0500 0.0600 0.0699 0.0799 0.0899 Columns 11 through 15 0.0998 0.1098 0.1197 0.1296 0.1395 Columns 16 through 20 0.1494 0.1593 0.1692 0.1790 0.1889 Columns 21 through 25 0.1987 0.2085 0.2182 0.2280 0.2377 Columns 26 through 30 0.2474 0.2571 0.2667 0.2764 0.2860 Columns 31 through 35 0.2955 0.3051 0.3146 0.3240 0.3335 Columns 36 through 40 0.3429 0.3523 0.3616 0.3709 0.3802 Columns 41 through 45 0.3894 0.3986 0.4078 0.4169 0.4259 Columns 46 through 50 0.4350 0.4439 0.4529 0.4618 0.4706 Columns 51 through 55 0.4794 0.4882 0.4969 0.5055 0.5141 Columns 56 through 60 0.5227 0.5312 0.5396 0.5480 0.5564 Columns 61 through 65 0.5646 0.5729 0.5810 0.5891 0.5972 Columns 66 through 70 0.6052 0.6131 0.6210 0.6288 0.6365 Columns 71 through 75 0.6442 0.6518 0.6594 0.6669 0.6743 Columns 76 through 80 0.6816 0.6889 0.6961 0.7033 0.7104 Columns 81 through 85 0.7174 0.7243 0.7311 0.7379 0.7446 Columns 86 through 90 0.7513 0.7578 0.7643 0.7707 0.7771 Columns 91 through 95 0.7833 0.7895 0.7956 0.8016 0.8076 Columns 96 through 100 0.8134 0.8192 0.8249 0.8305 0.8360 Column 101 0.8415 unmkpp(pp) ans = Columns 1 through 5 0 0.1000 0.2000 0.3000 0.4000 Columns 6 through 10 0.5000 0.6000 0.7000 0.8000 0.9000 Column 11 1.0000 [x,P]=unmkpp(pp) x = Columns 1 through 5 0 0.1000 0.2000 0.3000 0.4000 Columns 6 through 10 0.5000 0.6000 0.7000 0.8000 0.9000 Column 11 1.0000 P = -0.1652 -0.0003 1.0000 0 -0.1652 -0.0499 0.9950 0.0998 -0.1615 -0.0994 0.9801 0.1987 -0.1567 -0.1479 0.9553 0.2955 -0.1501 -0.1949 0.9211 0.3894 -0.1422 -0.2399 0.8776 0.4794 -0.1325 -0.2826 0.8253 0.5646 -0.1227 -0.3223 0.7648 0.6442 -0.1073 -0.3592 0.6967 0.7174 -0.1073 -0.3913 0.6216 0.7833 help csape CSAPE Cubic spline interpolation with various end conditions. PP = CSAPE(X,Y) returns the cubic spline interpolant (in ppform) to the given data (X,Y) using Lagrange end conditions (the default in table below). The interpolant matches, at the data site X(j), the given data value Y(:,j), j=1:length(X). The data values may be scalars, vectors, matrices, or even ND-arrays. Data points with the same site are averaged. For interpolation to gridded data, see below. PP = CSAPE(X,Y,CONDS) uses the end conditions specified by CONDS, with corresponding end condition values endcondvals . If there are two more data values than data sites, then the first (last) data value is taken as the value for the left (right) end condition, i.e., endcondvals = Y(:,[1 end]). Otherwise, default values are used. CONDS may be a *string* whose first character matches one of the following: 'complete' or 'clamped', 'not-a-knot', 'periodic', 'second', 'variational', with the following meanings: 'complete' : match endslopes (as given, with default as under *default*) 'not-a-knot' : make spline C^3 across first and last interior break (ignoring given end condition values if any) 'periodic' : match first and second derivatives at first data point with those at last data point (ignoring given end condition values if any) 'second' : match end second derivatives (as given, with default [0 0], i.e., as in variational) 'variational' : set end second derivatives equal to zero (ignoring given end condition values if any) The *default* : match endslopes to the slope of the cubic that matches the first four data at the respective end. By giving CONDS as a 1-by-2 matrix instead, it is possible to specify *different* conditions at the two endpoints, namely CONDS(i) with value endcondvals(:,i), with i=1 (i=2) referring to the left (right) endpoint. CONDS(i)=j means that the j-th derivative is being specified to be endcondvals(:,i) , j=1,2. CONDS(1)=0=CONDS(2) means periodic end conditions. If CONDS(i) is not specified or is different from 0, 1 or 2, then the default value for CONDS(i) is 1 and the default value of endcondvals(:,i) is taken. If no end condition values are specified, then the default value for endcondvals(:,i) is taken to be deriv. of cubic interpolant to nearest four points, if CONDS(i)=1; 0 if CONDS(i)=2. For backward compatibility, it is also possible to specify endcondvals as an optional fourth input argument, but only for univariate data. It is also possible to handle gridded data, by having X be a cell array containing m univariate meshes and, correspondingly, having Y be an m-dimensional array (or an (m+length(d))-dimensional array if the function is to have values of size d ). Correspondingly, CONDS is a cell array with m entries, and, consistent with the tensor product procedure used, also the tensor product of the univariate end conditions is enforced and must be supplied with values; see the example below. For example, fnplt(csape( [0:4], [1 0 -1 0 1;0 1 0 -1 0], 'periodic')), axis equal plots a circle, while x = linspace(0,2*pi,21); pp = csape( x, [1 sin(x) 0], [1 2] ); gives a good approximation to the sine function on the interval [0 .. 2*pi] (matching its slope 1 at the left endpoint, x(1) = 0, and its second derivative 0 at the right endpoint, x(21) = 2*pi, in addition to its value at every x(i), i=1:21). As an illustration of the specification of end conditions in a multivariate setting, here is complete bicubic interpolation, with the data explicitly derived from the bicubic polynomial g(x,y) = x^3y^3, to make it easy for you to see exactly where the slopes and slopes of slopes (i.e., cross derivatives) must be placed in the data values supplied. Since our g is a bicubic polynomial, its interpolant, f , must be g itself. We test this. sites = {[0 1],[0 2]}; coefs = zeros(4,4); coefs(1,1) = 1; g = ppmak(sites,coefs); Dxg = fnval(fnder(g,[1 0]),sites); Dyg = fnval(fnder(g,[0 1]),sites); Dxyg = fnval(fnder(g,[1 1]),sites); f = csape(sites,[Dxyg(1,1), Dxg(1,:), Dxyg(1,2); ... Dyg(:,1), fnval(g,sites), Dyg(:,2) ; ... Dxyg(2,1), Dxg(2,:), Dxyg(2,2)], ... {'complete','complete'}); if any(squeeze(fnbrk(f,'c'))-coefs), 'something went wrong', end As a multivariate vector-valued example, here is a sphere, done as a parametric bicubic spline, using prescribed slopes in one direction and periodic end conditions in the other: x = 0:4; y=-2:2; s2 = 1/sqrt(2); clear v v(3,:,:) = [0 1 s2 0 -s2 -1 0].'*[1 1 1 1 1]; v(2,:,:) = [1 0 s2 1 s2 0 -1].'*[0 1 0 -1 0]; v(1,:,:) = [1 0 s2 1 s2 0 -1].'*[1 0 -1 0 1]; sph = csape({x,y},v,{'clamped','periodic'}); values = fnval(sph,{0:.1:4,-2:.1:2}); surf(squeeze(values(1,:,:)),squeeze(values(2,:,:)),squeeze(values(3,:,:))) % the previous two lines could have been replaced by: fnplt(sph) axis equal, axis off See also csapi, spapi, spline. Reference page in Help browser doc csape x=[1,2,3,4]; y=[2,4,5,7]; xx=linspace(1,4,100); yy=spline(x,y,xx); plot(xx,yy,x,y,'*') t=[1,2,3,4]; x=[2,3,4,6,4,3]; y=[1,0,-1,0,2,4]; plot(x,y,'*') tt=linspace(1,4,100); xx=spline(t,x,tt); yy=spline(t,y,tt); plot(xx,yy) plot(xx,yy,x,y,'*') x=[2,3,4,6,4,3,2]; y=[1,0,-1,0,2,4,1]; xx=spline(t,x,tt); ??? Error using ==> chckxy at 89 The number of sites, 4, is incompatible with the number of values, 7. Error in ==> spline at 55 [x,y,sizey,endslopes] = chckxy(x,y); t=[1,2,3,4,5,6]; x=[2,3,4,6,4,3,2]; t=[1,2,3,4,5,6,7]; x=[2,3,4,6,4,3,2]; y=[1,0,-1,0,2,4,1]; xx=spline(t,x,tt); yy=spline(t,y,tt); plot(xx,yy,x,y,'*') tt=linspace(1,7,100); xx=spline(t,x,tt); yy=spline(t,y,tt); plot(xx,yy,x,y,'*') x = [1,2,3,4]; y = sin(x); pp = spline(x,y); f = @(x) ppval(pp,x); quad(f,1,4) ans = 1.2135 diary off