function [der_fce,D_P,D_L,D_C]=derivace_vypocet(f, x0, h, zobraz)
%
% Vypocte numericky prvni derivaci symbolicky zadane funkce f 
% v bode x0 s krokem h pomoci zakladnich pomernych diferenci
% D_P ... prava
% D_L ... leva
% D_C ... centralni

is_symbolic=~isempty(ver('symbolic'));
%is_toolbox('symbolic');

if ~is_symbolic
  disp(' ');
  disp('Lituji, ale neni k dispozici Symbolic Math Toolbox !');
  disp(' '); 
else

if zobraz==1
 clc;
 disp('----------------------------------------------------------------');
 fprintf('  Numericky vypocet derivace funkce f(x)=%s\n',f);
 fprintf('  v bode x0=%f s krokem h=%f\n',x0,h);
 disp('----------------------------------------------------------------');

 figure(1);
 clf;
 hold on;
 zoom on;
end

x=sym('x');

if verzematlabu==0
 % verze do 2017a
 f=sym(f);
 sym_derivace=diff(sym(f),x);
else
 % verze od 2017b
 sym_derivace=diff(str2sym(f),x);
end

x=x0;
%derivace_f_x0=eval(sym_derivace)
derivace_f_x0=double(subs(sym_derivace,x));
der_fce=derivace_f_x0;
f_x0=eval(f);
x=x0+h;
f_x0_plus_h=eval(f);
x=x0-h;
f_x0_minus_h=eval(f);

D_P=(f_x0_plus_h-f_x0)/h;
D_L=(f_x0-f_x0_minus_h)/h;
D_C=(f_x0_plus_h-f_x0_minus_h)/h/2;

if zobraz==1
 okoli=1.5*h;
 deleni=h/100;
 i=1;
 for x=x0-3*h:deleni:x0+3*h
   ff(i)=eval(f);
   i=i+1;
 end

 xx=x0-3*h:deleni:x0+3*h;
 pomf1=f_x0-derivace_f_x0*okoli;
 pomf2=f_x0+derivace_f_x0*okoli;
 pomf3=f_x0-D_P*okoli;
 pomf4=f_x0+D_P*(h+okoli);
 pomf5=f_x0-D_L*(okoli+h);
 pomf6=f_x0+D_L*okoli;
 pomf7=f_x0_minus_h-D_C*okoli;
 pomf8=f_x0_plus_h+D_C*okoli;

 max_y=max([ff,pomf1,pomf2,pomf3,pomf4,pomf5,pomf6,pomf7,pomf8]);
 min_y=min([ff,pomf1,pomf2,pomf3,pomf4,pomf5,pomf6,pomf7,pomf8]);
 sy=max_y-min_y;

 axis([x0-3.5*h x0+3.5*h min_y-.15*sy max_y+.15*sy]);
 plot(xx,ff,'b-');
 plot(x0,f_x0,'ro');
 texty(1,:)=sprintf('zadana funkce                            ');
 texty(2,:)=sprintf('bod [x0,f(x0)]                           ');
 legend(texty);
 pause;
 
 plot([x0-okoli x0+okoli],[pomf1, pomf2],'k-'); 
 texty(3,:)=sprintf('tecna ........ presna hodnota derivace   ');
 legend(texty);
 disp(' ');
 fprintf(' Presna hodnota derivace je                 f`(x_0) = %f\n',derivace_f_x0);
 disp(' ');
 pause;
  
 plot([x0-okoli x0+h+okoli],[pomf3, pomf4],'g-'); 
 plot(x0+h,f_x0_plus_h,'r>');
 texty(4,:)=sprintf('secna(1) ... prava pomerna diference     ');
 texty(5,:)=sprintf('bod [x0+h,f(x0+h)]                       ');
 legend(texty);
disp(sprintf([' Priblizna hodnota derivace',...
                '                                                  chyba ']));
 fprintf(...
 '     - pomoci prave pomerne diference       D_P(f,x0) = %f          %f\n',D_P,D_P-derivace_f_x0);
 pause;              
 
 plot([x0-okoli-h x0+okoli],[pomf5, pomf6],'c-'); 
 plot(x0-h,f_x0_minus_h,'r<');
 texty(6,:)=sprintf('secna(2) ... leva pomerna diference      ');
 texty(7,:)=sprintf('bod [x0-h,f(x0-h)]                       ');
 legend(texty);
 fprintf(...
 '     - pomoci leve pomerne diference        D_L(f,x0) = %f           %f\n',D_L,D_L-derivace_f_x0);
 pause;               

 plot([x0-okoli-h x0+okoli+h],[pomf7, pomf8],'m-'); 
 texty(8,:)=sprintf('secna(3) ... centralni pomerna diference ');
 legend(texty);
 fprintf(...
 '     - pomoci centralni pomerne diference   D_C(f,x0) = %f          %f\n',D_C,D_C-derivace_f_x0);
 disp(' ');
end

end







