%% calculate the second starting field at the bottom

function z1=cal_z1(R,W,P,kn,Ch,Ch1,H,h)
kh=W/Ch;
A=4*H*kh*sqrt(1-(kn/kh)^2);
D=H^2*((W/Ch)^2-kn*kn);
W=H^2*((W/Ch1)^2-kn*kn);
z1=R*(1+h/R*sqrt(A)-h^2/3*D-h^2/6*W-h^3/(6*R)*D*sqrt(A));
