c***********************************************************************
c*     subroutine sbwind
c***********************************************************************
      subroutine sbwind
     i (wustfl,awu, wind_dir, ntstep, intstep,
     o aicsr)
c
c     +++ PURPOSE +++
c     to update wzzo at each grid point;
c     To update  soil friction velocity on each grid point
c     and modify it for barriers and hills;
c     To initialize en thresh. and cp thresh. fr. velocites on grid;
c     To calculate max ratios of friction velocity to threshold
c     friction velocity
c
c     +++ ARGUMENT DECLARATIONS +++
      integer wustfl,intstep, ntstep
      real awu, aicsr, wind_dir
c
c     +++ ARGUMENT DEFINITIONS +++
c     intstep  - current index of ntstep thru time
c     ntstep   - max. no. of time steps in day
c     icsr     - index of current subregion.
c     wusp     - subregion soil threshold friction vel. trans. cap. (m/s)
c     aicsr    - max ratio of friction velocity to thresh. friction vel.
c     imax     - no. grid intervals in x-direction.
c     jmax     - no. grid intervals in y-direction.
c     wus      - soil friction velocity at grid points corrected for
c                hills and barriers (m/s).
c     wust     - threshold fr. vel. for en. at grid points
c     wusp     - threshold fr. vel. for trans. cap. at grid points
c     wind_dir - direction of wind (degrees from north)
c
c     + + + GLOBAL COMMON BLOCKS + + +
*$noreference
      include  'p1werm.inc'
      include  'm1geo.inc'
      include  'b1glob.inc'
      include  'h1db1.inc'
      include  'p1const.inc'
      include  's1agg.inc'
      include  's1dbh.inc'
      include  's1sgeo.inc'
c
c     + + + LOCAL COMMON BLOCKS + + +
      include  'erosion/m2geo.inc'
      include  'erosion/w2wind.inc'
      include  'erosion/e2grid.inc'
      include  'erosion/e3grid.inc'
      include  'erosion/s2agg.inc'
      include  'erosion/s2sgeo.inc'
      include  'erosion/s2surf.inc'
c
*$reference
c     +++ LOCAL VARIABLES +++
      integer i,j, icsr,k
c
       real at, wzzo, wzzov, rintstep
c
c     + + + END SPECIFICATIONS + + +

c     assign subregion index, currently only one
      icsr = 1
      aicsr = 0.


c     loop through grid interior to update
      do 40 i = 1, imax-1
      do 30 j = 1, jmax-1
c
c     update aerodynamic roughness
c ^^^ tmp out
c     write (*,*) 'in sbwind, call to sbzo'

      call sbzo
     i (wind_dir, asargo(icsr), szrgh(i,j), asxrgs(icsr), slrr(i,j),
     i  wzoflg, abrlai(icsr), abrsai(icsr), abzht(icsr), asxdks(icsr),
     o  wzzo, wzzov, awzzo)
c
c ^^^ tmp out
c      write (*,*) 'in sbwind, call to sbwus'
c     update surface (below canopy) friction velocity
      call sbwus
     i (anemht, awzzo, awu, wzzov, abrlai(icsr), abrsai(icsr),
     o  wus(i,j) )

c     correct friction velocity for hills
c      if (nhill .ne. 0 ) then
c           wus(i,j) = wus(i,j) * w0hill(i,j,kbr)
c        endif
c
c     correct friction velocity for barriers
      if (nbr .ne. 0 ) then
         wus(i,j) = wus (i,j) * w0br (i,j,kbr)
      endif
c

      if (wustfl .eq. 1) then
c     update threshold friction velocities
c       calculate hour k for surface water content
        rintstep = intstep
        k = aint(rintstep*23.75/ntstep) + 1

      call sbwust
     i  (sf84(i,j), asdagd(1,icsr), sfcr(i,j), asvroc(1,icsr),
     i  sflos(i,j), abffcv(icsr),wzzo, ahrwc0(k,icsr), ahrwcw(1,icsr),
     o   wust(i,j), wusp(i,j) )
c
c ^^^ tmp out
c      write (*,*) 'in sbwind, call to sbwust'
c      write (*,*) 'i  j  wust sf84', i, j, wust(i,j),sf84(i,j)
c      write (*,*) 'sfcr sflos wzzo', sfcr(i,J),sflos(i,J),wzzo
c      write (*,*) 'abffcv, ahrwc0, wusp',abffcv(icsr),ahrwc0(12,icsr),
c     i             wusp(i,j)

      endif

       at = wus(i,j)/wust(i,j)
        aicsr = amax1(aicsr, at)
c

   30 continue
   40 continue
c
c     write (*,*) 'at exit sbwind aicsr =', aicsr
c     write (*,*) ' wus(3,3), wust(3,3)', wus(3,3), wust(3,3)
c
      return
      end
c++++ +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++