fitting a spline to titanium test data -凯发k8网页登录
this example shows how to use commands from curve fitting toolbox™ to fit a spline to titanium test data with manual and automatic selection of knots.
manual knot choice for spline interpolation
here are some data that record a certain property of titanium, measured as a function of temperature. we'll use it to illustrate some issues with spline interpolation.
[xx,yy] = titanium;
a plot of the data shows a rather sharp peak.
plot(xx,yy,'bx');
frame = [-10 10 -.1 .3] [min(xx),max(xx),min(yy),max(yy)];
axis(frame);
we pick a few data points from these somewhat rough data, since we want to interpolate. here is a picture of the data, with the selected data points marked.
pick = [1 5 11 21 27 29 31 33 35 40 45 49]; tau = xx(pick); y = yy(pick); hold on plot(tau,y,'ro'); hold off
since a spline of order k with n k knots has n degrees of freedom, and we have 12 data points, a fit with a fourth order spline requires 12 4 = 16 knots. moreover, this knot sequence t must be such that the i-th data site lies in the support of the i-th b-spline. we achieve this by using the data sites as knots, but add two simple knots at either end.
dl = tau(2) - tau(1); dr = tau(end) - tau(end-1); t = [tau(1)-dl*[2 1] tau tau(end) dr*[1 2]]; % construct the knot sequence plot(tau,y,'ro'); hold on axis(frame [-2*dl 2*dr 0 0]) plot(t,repmat(frame(3) .03,size(t)),'kx') hold off legend({'data values' 'knots'},'location','nw')
we use this knot sequence to construct an interpolating cubic spline.
sp = spapi(t,tau,y);
now, for the plot. since we do not care about the part of the spline outside the data interval, we restrict the plot to that interval.
plot(tau,y,'ro') axis(frame) hold on fnplt(sp,[tau(1) tau(end)], 'k') hold off
a closer look at the left part of the spline fit shows some undulations.
xxx = linspace(tau(1),tau(5),41); plot(xxx, fnval(sp, xxx), 'k', tau, y, 'ro'); axis([tau(1) tau(5) 0.6 1.2]);
the unreasonable bump in the first interval stems from the fact that our spline goes smoothly to zero at its first knot. to see that, here is a picture of the entire spline, along with its knot sequence and the data points.
fnplt(sp,'k'); hold on plot(tau,y,'ro', t,repmat(.1,size(t)),'kx'); hold off legend({'spline interpolant' 'data values' 'knots'},'location','nw')
here is a simple way to enforce a more reasonable boundary behavior. we add two more data points outside the given data interval and choose as our data there the values of the straight line through the first two data points.
tt = [tau(1)-[4 3 2 1]*dl tau tau(end) [1 2 3 4]*dr]; xx = [tau(1)-[2 1]*dl tau tau(end) [1 2]*dr]; yy = [y(1)-[2 1]*(y(2)-y(1)) y y(end) [1 2]*(y(end)-y(end-1))]; sp2 = spapi(tt,xx,yy); plot(tau,y,'ro', xx([1 2 end-1 end]),yy([1 2 end-1 end]),'bo'); axis(frame [-2*dl 2*dr 0 0]); hold on fnplt(sp2,'b',tau([1 end])) hold off legend({'original data' 'data added for end conditions' ... 'fit with added data'},'location','nw')
here is a comparison of the two spline fits, to show the reduction of the undulation in the first and last interval.
hold on fnplt(sp,'k',tau([1 end])) hold off legend({'original data' 'data added for end conditions' ... 'fit with added data' 'original fit'},'location','nw')
finally, here is a closer look at the first four data intervals that shows more clearly the reduction of the undulation near the left end.
plot(tau,y,'ro',xxx,fnval(sp2,xxx),'b',xxx,fnval(sp,xxx),'k'); axis([tau(1) tau(5) .6 1.2]); legend({'original data' 'fit with added data' ... 'original fit'},'location','nw')
automatic knot choice for interpolation
if all this detail turns you off, let curve fitting toolbox choose the knots for you. specify the desired order of the interpolant as the first input argument to the spline interpolation command spapi, rather than a knot sequence.
autosp = spapi(4, tau, y); knots = fnbrk(autosp,'knots'); plot(tau, y, 'ro') hold on fnplt(autosp,'g') plot(knots, repmat(.5,size(knots)),'gx') hold off legend({'data values' 'fit with knots chosen by spapi' ... 'knots chosen by spapi'}, 'location','nw')
below is the result of a much better knot choice, obtained by shifting the knot at 842 slightly to the right and the knot at 985 slightly to the left.
knots([7 12]) = [851, 971]; adjsp = spapi(knots, tau, y); hold on fnplt(adjsp,'r',2) plot(knots, repmat(.54,size(knots)),'rx') hold off legend({'data values' 'fit with knots chosen by spapi' ... 'knots chosen by spapi' 'fit with knots adjusted' ... 'adjusted knots'}, 'location','nw')
or else, simply try the standard cubic spline interpolant, supplied by csapi. this amounts to a slightly different choice of knots.
autocs = csapi(tau, y); plot(tau, y, 'ro') hold on fnplt(autocs,'c') hold off
with such rapidly-varying data, it is hard to get agreement among all reasonable interpolants, even if each of them is a cubic spline. the plot below shows all five interpolants, for comparison.
plot(tau, y, 'ro') hold on fnplt(sp,'k',tau([1 end])) % black: original fnplt(sp2,'b',tau([1 end])) % blue: with special end conditions fnplt(autosp,'g') % green: automatic knot choice by spapi fnplt(autocs,'c') % cyan: automatic knot choice by csapi fnplt(adjsp,'r',2) % red: knot choice by spapi slightly changed hold off legend({'data values' 'original fit' 'special end conditions' ... 'with knots chosen by spapi' 'with knots chosen by csapi' ... 'with adjusted knots'},'location','nw')
select a web site
choose a web site to get translated content where available and see local events and offers. based on your location, we recommend that you select: .
you can also select a web site from the following list:
americas
- (español)
- (english)
- (english)
europe
- (english)
- (english)
- (deutsch)
- (español)
- (english)
- (français)
- (english)
- (italiano)
- (english)
- (english)
- (english)
- (deutsch)
- (english)
- (english)
- switzerland
- (english)
asia pacific
- (english)
- (english)
- (english)
- 中国
- (日本語)
- (한국어)