% Example:
% n = 10; x = [1:(3 - 1)/(n - 1):3]; y = 10*rand(1, n);
% [A, c] = my_polyfit(x, y);
%
% http://www-hm.ma.tum.de/archiv/in1/ws0102/links/Interpol/Lagrange.html
function [A, c] = my_polyfit(x, y)
[m n] = size(x);
A = zeros(n, n);
for k=1:n
    A(:,k)=x'.^(k - 1); 
end
c = A\y';
xb = x(1, 1);
xe = x(1, n);
xp = xb:(xe - xb)/100:xe;
yp = polyval((c(n:-1:1))', xp);
plot(x, y, 'o', xp, yp, '-')