function pxx=lagrange_interp2(x,f,xx)
% PXX=LAGRANGE_INTERP2(X,F,XX)
% Deadas unas abscisas X=[X(1),X(2),...,X(n)], los correspondientes
% valores de la funcion F=[F(1),F(2),...,F(n)], y un conjunto de
% puntos XX=[XX(1),XX(2),...,XX(M)]; calcula los valores del polinomio
% interpolador P(1), P(2),...,P(M) en cada uno de los puntos
% XX(1),XX(2),...,XX(M).
%
% INPUT:
%  x=[x(1),x(2),..,x(n)] : conjunto de abcisas para calcular el
%                          polinomio interpolador.
%  f=[f(1),f(2),...,f(n)]: valores de la funcion f en las abcisas
%                          x(1),x(2),...,x(n) respectivamente.
% xx=[xx(1),xx(2),...,xx(m)] : conjunto de puntos sobre el eje X en los
%                          que se evaluara el polinomio interpolador.
% OUTPUT:
%  pxx=[pxx(1),pxx(2),...,pxx(m)]:
%                          valor del polinomio interpolador en los
%                          puntos xx(1),xx(2),...,xx(m).
n=length(x);
m=length(xx);
pxx=0;
for i=1:n
  L=ones(1,m);
  for j=1:n
    if j ~= i
      L = L .* (xx-x(j))/(x(i)-x(j));
    end
  end
  pxx = pxx + f(i)*L;
end

return
end
