function nevilluv_algoritmus(zobraz_vse,x,f,alfa)
%
% Vypocte hodnotu Newtonova interpolacniho polynomu v bode alfa 
% pro funkci zadanou tabulkou [x;f(x)]
%
% nejsou-li zadane vektory x,f => zadaji se mysi 
% neni-li zadano alfa => zada se rucne
% pro zobraz_vse=1 zobrazi vsechny dilci ukony

clc;
clear texty;
barva=['r','m','g','c','k'];
styl=['- ';'--';'-.';': '];
indikator=0;

if ~exist('zobraz_vse')
  zobraz_vse=0;
end

if ~exist('x') || ~exist('f')
  [x,f]=zadej_body;
  indikator=1;
  clc;
end
if ~exist('alfa')
  alfa=input(' Zadej alfa = ');
  clc;
end

presah=10;
pocet_deleni=1000;
[sortx,~]=sort(x);
n=max(size(x));
deleni=(sortx(n)-sortx(1))/pocet_deleni;
[~,ind]=sort(abs(x-alfa));
%[v,ind]=sort(abs(x-alfa),'descend');
x=x(ind);
f=f(ind);

disp('-----------------------------------------------------------------------');
fprintf(' Nevilluv algoritmus pro vypocet hodnoty funkce f v bode alfa=%f\n',alfa);
disp(' ');
fprintf(' pro funkci zadanou tabulkou \n');
disp(' ');
fprintf(' |    x(k)    |    f(k)    |\n');
fprintf(' ===========================\n');
for i=1:n
  fprintf(' | %9.4f  | %9.4f  |\n',x(i),f(i));
end  

A=zeros(n);
A(:,1)=f';
for k=2:n
  for i=k:n
    A(i,k)=A(i,k-1)+(alfa-x(i))*(A(i,k-1)-A(i-1,k-1))/(x(i)-x(i-k+1));
  end
end

if zobraz_vse==1 | indikator==1
 clear texty;
 if alfa<min(x) || alfa>max(x)
   pres(1)=220*(min(abs(alfa-x)))/(max(x)-min(x));
 else
   pres(1)=10;
 end  
  
 [xk,fk]=polynom(x,f,deleni,pres(1)/2);

 xmin=min(min(xk),alfa);
 xmax=max(max(xk),alfa);
 fmin=min(fk);
 fmax=max(fk);
 for i=2:n
  pomf=f(1:i);
  pomx=x(1:i);
  if alfa<min(pomx) || alfa>max(pomx)
    pres(i)=220*(min(abs(alfa-pomx)))/(max(pomx)-min(pomx));
  else
    pres(i)=10;
  end 
  [~,fk]=polynom(pomx,pomf,deleni,pres(i)/2);
  if fmin>min(fk)
    fmin=min(fk);
  end
  if fmax<max(fk)
    fmax=max(fk);
  end
 end
 disp(' ');
 disp('Stiskni klavesu');
 disp(' ');
 pause;
 fprintf('    x(k)   |   f(k)   |   Aproximace f(alfa)  \n');
 fprintf(' ================================================\n');
 figure(1);
 clf;
  set(figure(1),'Units','centimeters');
 set(figure(1),'Position',[7 3 28 15]);
 for i=1:n
  clf;
  hold on;
  axis([xmin-presah/100*(xmax-xmin) xmax+presah/100*(xmax-xmin) ...
        fmin-presah/100*(fmax-fmin) fmax+presah/100*(fmax-fmin)]); 
  zoom on;
  for k=1:min(i,n-1)
   pomf=f(1:k);
   pomx=x(1:k);
   [xk,fk]=polynom(pomx,pomf,deleni,pres(k)/2);
   color=strcat(barva(mod(k-1,5)+1),styl(mod(floor((k-1)/5),5)+1,:));
   if k==1
    plot([min(xk,alfa)-abs(alfa-xk)/10, max(xk,alfa)+abs(alfa-xk)/10],[fk,fk],color);
   else
    plot(xk,fk,color);
   end  
   popis=sprintf('polynom stupne %d            ', k-1);
   texty(k,:)=popis(1:25); 
  end
  if i==n
   popis=sprintf('Interpolacni polynom            ');
   texty(i,:)=popis(1:25); 
   [xk,fk,~]=polynom(x,f,deleni,pres(i)/2);  
   plot(xk,fk,'b-');
  end 
  legend(texty,'Location','bestoutside');
  bbb=plot(x,f,'ro');
  set(get(get(bbb,'Annotation'),'LegendInformation'),'IconDisplayStyle','off');
  for j=1:i
    bbb=plot(alfa,A(j,j),'ko');
    set(get(get(bbb,'Annotation'),'LegendInformation'),'IconDisplayStyle','off');
    text(alfa+deleni*10,A(j,j),num2str(j));
  end

  tiskni=sprintf(' %9.4f |%9.4f ',x(i),f(i));
  for j=2:i
    tiskni=[tiskni,sprintf('|%9.4f ',A(i,j))];
  end
  disp(tiskni);
  
  if i~=n
    pause;
  end
 end
 else
   disp(' ');
   disp(' ');
   fprintf('    x(k)   |   f(k)   |   Aproximace f(alfa)  \n');
   fprintf(' ================================================\n');
   for i=1:n
     tiskni=sprintf(' %9.4f |%9.4f ',x(i),f(i));
     for j=2:i
       tiskni=[tiskni,sprintf('|%9.4f ',A(i,j))];
     end
     disp(tiskni);
   end
end

disp(' ');
disp(' ');
fprintf(' Priblizna hodnota  f(alfa) = %9.4f  \n',A(n,n));

