% solve interpolation with polynomials in Lagrange  form
clear; clf;
t = [0 0.5 1 6 7 9];
y = [0 1.6 2 2 1.5 0];
h = plot(t,y,'b');
set(h,'Linewidth',3);
hold on;
pts = 100;
%x = 0:1/pts:9;
x = linspace(t(1),t(length(t)),pts);
dg = 5;  % degree= number of data points -1
% vector of evaluated intrerpolating polynomial
p = zeros(1,pts);
for i = 0:dg,
   lg = ones(1,pts);
   for j = 0:dg,
      if i==j, continue; end;
      h = (x-t(j+1))/(t(i+1)-t(j+1));
      %if i==0, size(lg), size(h) end;
      lg = lg.*h;
   end;
   p = p + y(i+1)*lg; 
end;
h = plot(x,p,'g');
set(h,'Linewidth',3);
