A1=0.6890;
z=[0 10 20 30 40 50 60]; % Data points taken along z - axis (independent)
d=[3.2 3.1 2.1 3.8*(1+A1) 4 4.5 3.9]; %[units:mm] meaured values of d at respective points along z

%=====================================================
%Polynomial Interpolation (Lagrange)
%=====================================================

syms z
L=0;
for ii=1:length(z)
   Lj=1; % resets Lj to 1 for each iteration
   for jj=1:length(z)
       if ii~=jj %ends iteration when ii=jj
           Lj=Lj*(z-z(jj))/(z(ii)-z(jj));%Formulating the Lj term using the lagrange equation
       end
   end
   L=L + d(ii)*Lj;
end

%=====================================================
%Plot of Lagrange equation
%=====================================================

ezplot(L, [0 60])
xlabel('Surface Coordinate z (mm)')
ax = gca;
ax.FontWeight = ('bold');
ylabel('Thickness d (mm)')
title('Polynomial (Lagrange) Interpolation')
hold on
plot(z,d,'o')
hold off