function RES=InterPol1(X,Y)
% function RES=InterPol1(X,Y)
% Знаходження коефіцієнтів інтерполяційного полінома
% "Numerical recipes in C", 3.5, p.121, функція polcoe
% Арументи:
%   X - вектор-рядок аргументів функції
%   Y - вектор-рядок значень функції
% Результати:
%   RES - вектор-рядок коефіцієнтів інтерполюючого полінома, a0, a1, a2,...
% Обмеження:
%   1) X та Y повинні мати ОДНАКОВУ кількість елементів

% визначаємо, чи має X один рядок
[m,n]=size(X);
% якщо ні - транспонуємо
if (m>1) X=X'; [m,n]=size(X); end

% визначаємо, чи має Y один рядок
[m,n]=size(Y);
% якщо ні - транспонуємо
if (m>1) Y=Y'; [m,n]=size(Y); end

S=zeros(1,n);
RES=zeros(1,n);

S(1,n)=-X(1,1);
for i=2:n                               % коефіцієнти Si головного полінома P(x)
    for j=n-i+1:n-1                     % знаходяться рекурентно
        S(1,j)=S(1,j)-X(1,i)*S(1,j+1);
    end
    S(1,n)=S(1,n)-X(1,i);
end

for j=1:n
    phi=n;
    for k=n:-1:2                        % визначення phi=П    (Xj-Xk) знаходиться як похідна
        phi=(k-1)*S(1,k)+X(1,j)*phi;    %                 j<>k
    end                                 % від P(x)
    disp(fprintf('Phi = %f',phi));
    ff=Y(1,j)/phi;
    b=1;                                % коефіцієнти поліномів на кожному кроці формули Лагранжа
    for k=n:-1:1                        % знаходяться штучним діленням P(x) на (x-xj).
        RES(1,k)=RES(1,k)+b*ff;         % Розв"язок Ck накопичується
        b=S(1,k)+X(1,j)*b;
    end
end

clear m n S i j k phi b;