% plot Lagrange polynomials on [a,b]

a = -1;
b =  1;
N =  7;



% points for interpolation
Xb = linspace(a,b,N);
% plot the interpolation points
figure(1); clf;
fs = 24; %font size
plot(Xb,zeros(size(Xb)),'x', ...
     'LineWidth',3,'MarkerSize',10);
grid on;
xlabel('x', 'FontSize', fs); ylabel('y', 'FontSize', fs);
hold on;


% points for visualization
Xf = linspace(a,b,100);
ColorSet = lines(N);
for i = 1:N,
    % form the Lagrange polynomial i
    Li = @(x) prod(x-Xb([1:i-1,i+1:N]))/prod(Xb(i)-Xb([1:i-1,i+1:N]));
    % plot it
    plot(Xf,arrayfun(Li,Xf),'-',...
         Xb,arrayfun(Li,Xb),'o',...
         'LineWidth',3,'Color',ColorSet(i,:));
    if(i<N), pause; end
end
hold off;