function der=nd4_13(Fun,x0,h);
%%  function der=nd4_13(Fun,x0,h);
%%
%%  Numerical differentiation using the formula:
%%       f'(x_0)=(f(x0+h)-f(x0-h))/2h
%%
%%  INPUT:
%%    Fun - the function data which may be entered EITHER as
%%           a string, enclosed in ' ', with the name the m-file of the 
%%           function ,  
%%           OR as a row vector of data [f(x0-h),f(x0+h)]
%%     x0 -  the point at which the derivative is to be estimated
%%     h  -  the step size h (several h may be given in a row vector
%%            when Fun is given by a m-file).
%%  OUTPUT:
%%    der - the approximate value of the derivative (or column of
%%          approximate values if more than one h is given).


if isstr(Fun)==1;  % test whether from a function file or from given data
   xh=[x0-h;x0+h];
   values(:,1:length(h))=feval(Fun,xh); %function values 
else
    values(:,1)=Fun'; % given data as 2nd column
end;

der=(values(2,:)-values(1,:))./(2*h);
