Lagrange Polynomial Interpolation

x = [20,41,58,75];
y = [0.122, 0.361, 0.520, 0.634];

function [ysz] = lagrange(xsz, x, y)
  n = 4;
  ysz = 0;
  A = [];
  for i=1:n
    p = 1;
    for j=1:n
      if j != i
        p = conv(p,(xsz-x(j)) / (x(i)-x(j)));
      endif
    endfor;
    A(i) = p;
    ysz += y(i) * p;
  endfor;
endfunction;
argument = [20, 25, 36, 41, 53, 58, 75, 80]
result_lagrange = [0.122, lagrange(25,x,y),lagrange(36,x,y), 0.361, lagrange(53,x,y), 0.520, 0.634, lagrange(80,x,y)]
plot(argument, result_lagrange);