anchor.diff:  c    +                    - 0.3606d0)
               rectchan = 0.5d0 * width(s) * depth(s)
               transslope(s) = transfactor * curvature(s) 
  -            xmagtransslope = sign(transslope(s), 1.d0)
  +            xmagtransslope = abs(transslope(s))
               if (xmagtransslope .le. corner) then
  -               Acs(s) = rectchan - width(s) ** 2.d0 * 
  +               Acs(s) = rectchan - width(s) ** 2 * 
        +                  xmagtransslope / 8.d0
                  rightdepth(s) = depth(s) + width(s) * 
        +                         transslope(s) / 2.d0
                  leftdepth(s) = depth(s) - width(s) * transslope(s) / 2.d0
               else
  -               Acs(s) = 0.5d0 * depth(s) ** 2.d0 / xmagtransslope
  +               Acs(s) = 0.5d0 * depth(s) ** 2 / xmagtransslope
                  if (curvature(s) .lt. 0.d0) then
                     rightdepth(s) = 0.d0
                     leftdepth(s) = 2.d0 * depth(s)
  
positive.code:       real*8 delx(*), dely(*), dels(*), curvature(*), mag, sn, a, b, c,
       +       carg
        intrinsic acos
  c
        do s = 2, stnserod - 1
           a = dels(s) 
  
negative.code:       else
           do s = 1, stations - 1
              if (width(s) .ne. 0.d0) then
                 forcefactor = rho * vel(s) ** 2 / depth(s) 
       +              * (1.d0 - 2.d0 * depth(s) / width(s))
                 latforce(s) = 0.d0
                 lag(s) = 0.d0
                 xmagcurvep1 = abs(curvature(s + 1)) 
                 xmagcurve = abs(curvature(s)) 
                 if ((curvature(s + 1) * curvature(s) .gt. 0.d0 .and.
       +           xmagcurvep1 .gt. xmagcurve .and.
       +           xmagcurvep1 * width(s) .le. 2.d0) .or.
  