function yy=BarycentricInterpolation(x,y,lambda,xx)
% BARYCENTRICINTERPOLATION interpolate using the barycentric formula
%   yy=BarycentricInterpolation(x,y,lambda,xx) uses the points (x,y)
%   and the precalculated barycentric coefficients lambda for the
%   barycentric form of the interpolation polynomial P and iterpolates
%   the values yy=P(xx)

x=x(:); y=y(:); xx=xx(:);
nn=length(xx);
for i=1:nn,
  z=(xx(i)-x)+1e-30;  % prevents a division by zero
  mue=lambda'./z;
  yy(i)=mue'*y/sum(mue);
end;
