function  retval = bernpoly(f,x)
% usage: retval = bernpoly(f,x)
% description: Returns nth degree Bernstein polynomial for
% f(x) evaluated at vector of abscissas x, given input of
% x and vector f of values f(k/n), k=0,...,n.

n = length(f)-1;
m = length(x);
retval = zeros(size(x));
for j = 1:m
  xj = x(j);
    if (xj == 0)
      ytmp = f(1);
    elseif (xj == 1)
      ytmp = f(n+1);
    else
      onemx = 1-xj;
      ytmp =  f(1);
      coef = onemx^n;
      for k = 1:n
        coef = coef*(n-k+1)/k*xj/onemx;
        ytmp = ytmp + coef*f(k+1);
      end
    end
    retval(j) = ytmp;
end