How do I find the derivative of a spline curve in MATLAB 7.9 (R2009b)?

108 views (last 30 days)
I am trying to find a derivative of a spline using a POLYDER function. The result looks wrong.
x=[0 600 1000 1381 1887 3087 4265 5367 6267 6987 ; 40 100 100 200 200 800 700 349 400 40];
pp=spline(x(1,:),x(2,:));
xx=linspace(x(1,1),x(1,10),2000);
erg=ppval(pp,xx);
plot(xx,erg,'k');
hold on
[breaks,coefs,l,k,d] = unmkpp(pp);
pp2 = mkpp(breaks,repmat(k-1:-1:1,d*l,1).*coefs(:,1:k-1),d);
y_der=polyder(ppval(pp2,linspace(x(1,1),x(1,10),2000)));
plot(linspace(x(1,1),x(1,10),2000-1),y_der,'r');

Accepted Answer

MathWorks Support Team
MathWorks Support Team on 28 Oct 2010
The POLYDER function is meant to be used to differentiate polynomials, and not piece-wise polynomials, such as splines.
If you have Spline Toolbox, you can use FNDER function to find the derivative of a spline:
pp=spline(x,y);
p_der=fnder(pp,1);
y_prime=ppval(p_der,x);
For more information on the FNDER function, please refer to the documentation by executing the following in the MATLAB command prompt:
doc fnder
If you do not have access to the Spline Toolbox, you can use UNMKPP function to break down your polynomial and then use MKPP function to assemble a new polynomial that will be a derivative of the first polynomial as in the following example:
% create points we will use to construct the spline
x = 0:10;
y = (x-5).^3+3*3+x+5;
% points the spline will be plotted at
xx = linspace(0,10,20);
% create the spline
pp = spline(x,y);
% plot the spline and the initial points used to create it
figure
plot(x,y,'o',xx,ppval(pp,xx))
hold on
% extract details from piece-wise polynomial by breaking it apart
[breaks,coefs,l,k,d] = unmkpp(pp);
% make the polynomial that describes the derivative
pp2 = mkpp(breaks,repmat(k-1:-1:1,d*l,1).*coefs(:,1:k-1),d);
% plot the derivative of the polynomial
plot(xx,ppval(pp2,xx),'-r')
% to calculate 2nd derivative differentiate the 1st derivative
[breaks,coefs,l,k,d] = unmkpp(pp2);
pp3 = mkpp(breaks,repmat(k-1:-1:1,d*l,1).*coefs(:,1:k-1),d);
plot(xx,ppval(pp3,xx),'-g')
For more information on the UNMKPP and MKPP functions please refer to the documentation by executing the following in the MATLAB command prompt:
doc unmkpp
doc mkpp

More Answers (0)

Products


Release

R2009b

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!