clear all; close all

% simple example to illustrate polynomial interpolation
xd=[0 0.5 1 1.5 2]';
yd=[0 .19 .26 .29 .31]';


load peak.dat
xd=peak(:,1);
yd=peak(:,2);
plot(xd,yd,'O:')
% polynomial wiggle, look at fitted curve with axis([-3 3 -10 10])


load census.dat
xd=census(:,1);
yd=census(:,2);
plot(xd,yd,'O:')
% polyfit command fails to find interpolating polynomial
% polyinterp function (Lagrange formula) works


load velocity.dat
xd=velocity(:,1);
yd=velocity(:,2);
plot(xd,yd,'O:')
% polyfit command fails to find interpolating polynomial
% polyinterp function (Lagrange formula) works
% polynomial wiggle, look at fitted curve with axis([0 30 -10 100])


n=length(xd);
xp=linspace(xd(1),xd(n),100);

% solve the linear system for the coefficients
V=vander(xd)
c=V\yd
yp=polyval(c,xp);
plot(xd,yd,'bO',xp,yp,'r')

% using polyfit command
coeff=polyfit(xd,yd,n-1)
yp=polyval(coeff,xp);
plot(xd,yd,'bO',xp,yp,'r')

% using Lagrange interpolating polynomial
yp=polyinterp(xd,yd,xp);
plot(xd,yd,'bO',xp,yp,'r')
% symbolic expression of interpolating polynomial
symx=sym('x');
P=polyinterp(xd,yd,symx)
P=simplify(P)


