clear;
clc;


n=5;
X=[-1.32;-0.84;-0.41;0.16;0.33;1.09];
x_bar=0.25;

function y=f(x);
    y=exp(x);
endfunction;

for i = 1 : n
  Y(i+1)=f(X(i+1));
end;

function FH=test_X(X,n); // FH=[eq;h]
    FH(1)=1; // eq=1
    FH(2)=X(2)-X(1); // h=x1-x0
    for i=2:n
        if X(i+1)-X(i)<>FH(2) then
            FH=[0;0];
            break;
        end;
    end;
endfunction;

function y=omega(x,n,X);
    for i=0:n
        T(i+1)=x-X(i+1);
    end;
    y=prod(T);
endfunction;


D=zeros(n+1,n+1);
D(:,1)=Y;
if test_X(X,n)(1)==1 then
    h=test_X(X,n)(2);
    for k=1:n
        for i=0:(n-k)
            D(i+1,k+1)=D(i+2,k)-D(i+1,k);
        end;
    end;
    Tmp(i+1)=D(1,1);
    for i=1:k
        Tmp(i+1)=(D(1,i+1)*omega(x_bar,i-1,X))/(factorial(i)*h^i);
    end;
    Wn=sum(Tmp);
else
    // punkt 7
    for k=1:n
        for i=0:(n-k)
            D(i+1,k+1)=(D(i+2,k)-D(i+1,k))/(X(i+k+1)-X(i+1));
        end;
    end;
    Tmp(i+1)=D(1,1);
    for i=1:k
        Tmp(i+1)=D(1,i+1)*omega(x_bar,i-1,X);
    end;
    Wn=sum(Tmp);
end;

disp(Wn,"Wn=");