xc = (1/2) .* (1-cos(theta))
   zc = zeros(1,21);
   i=1;
   while xc(i)<=0.2025
       zc(i) = 2.6595*xc(i)*(xc(i)^2 - 0.6075*xc(i) + 0.1147);
       i = i+1;
   end
   while i<=21
       zc(i) = 0.02208*(1-xc(i));
       i = i+1;
   end