function [maxx,maxerr]=GenericInterp(x)

%% Function values (nodes given as input, a row vector)
x=sort(x); n=length(x)-1; a=x(1); b=x(n+1); y=f(x);

%% Barycentric weights
dd=(x-x')+eye(n+1);
wt=1./prod(dd);

%% Evaluation points
m=4; hh=(b-a)/(m*n); xx=a:hh:b;

figure(1);
yf=f(xx);
yp=p(xx,x,y,wt);
plot(xx,yp,'o-',xx,yf,'v-');
err=abs(yf-yp);
figure(2);
plot(xx,err,'s-');
[maxerr,maxi]=max(err); maxx=xx(maxi);

end

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function val=f(x)

val=1./(1+25*x.*x);

end

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function val=p(xx,x,y,wt)

dd=1./(xx-x');
[I,J]=find(isinf(dd));
denom=wt*dd;
numer=(wt.*y)*dd;
val=numer./denom;
val(J)=y(I);

end