format long
x = [cos(0.6)-0.6 cos(0.7)-0.7 cos(0.8)-0.8 cos(0.9)-0.9];
y = [0.6 0.7 0.8 0.9];
z = 0;
for k = 1:4
  L(k) = 1
  for i = 1:4
    if i != k 
      L(k) = L(k)*(z-x(i))/(x(k)-x(i));
      end
      end
      end
p = 0;
for j = 1:4
  p = p+y(j)*L(j)
  end