- b1(1,i-2)=bnew1(1,1,iti)*sin(theta(i-1)/2.0)
- & +bnew1(2,1,iti)*sin(theta(i-1))
- & +bnew1(3,1,iti)*cos(theta(i-1)/2.0)
- gtb1(1,i-2)=bnew1(1,1,iti)*cos(theta(i-1)/2.0)/2.0
- & +bnew1(2,1,iti)*cos(theta(i-1))
- & -bnew1(3,1,iti)*sin(theta(i-1)/2.0)/2.0
-c & +bnew1(3,1,iti)*sin(alpha(i))*cos(beta(i))
+ b1(1,i-2)=bnew1(1,1,iti)*dsin(theta(i-1)/2.0d0)
+ & +bnew1(2,1,iti)*dsin(theta(i-1))
+ & +bnew1(3,1,iti)*dcos(theta(i-1)/2.0d0)
+ gtb1(1,i-2)=bnew1(1,1,iti)*dcos(theta(i-1)/2.0d0)/2.0d0
+ & +bnew1(2,1,iti)*dcos(theta(i-1))
+ & -bnew1(3,1,iti)*dsin(theta(i-1)/2.0d0)/2.0d0
+c & +bnew1(3,1,iti)*dsin(alpha(i))*cos(beta(i))