!*==sbgrid.spg  processed by SPAG 6.70Rc at 15:33 on 10 Dec 2012
!*------------------   SPAG Configuration Options   --------------------
!*--0323,76 000101,-1 000000102011332010100002000000210211210,136 10 --
!*--1100000011112111000000000000,10,10,10,10,10,10,900,100 200000000 --
!*--000000010000000000000,72,72 73,42,38,33 00011112110000100000000  --
!*----------------------------------------------------------------------
 
!++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 
 
 
 
 
 
!$Author: joelevin $
!$Date: 2011-03-24 11:33:26 -0500 (Thu, 24 Mar 2011) $
!$Revision: 11724 $
!$HeadURL: https://eweru-dev1.eweru.ksu.edu/svn/code/weps1/branches/WEPS_F90_update/weps.src/src/lib_erosion/sbgrid.for $
!**********************************************************************
!     subroutine sbgrid
!**********************************************************************
      subroutine sbgrid
      use i_p1werm
      use i_m1geo
      use i_m1subr
      use i_m2geo
      use i_e2grid
      implicit none
!*--SBGRID29
!
!*** Start of declarations rewritten by SPAG
!
! Local variables
!
      real :: dxmin,lx,ly
      integer :: i,icsr,j,ngdpt
!
!*** End of declarations rewritten by SPAG
!
!
!     +++ PURPOSE +++
!     to calculate grid size and spacing for EROSION.
!     grid size assumes outer points are outside field boundary
!      a distance ix/2
!     to calculate number of grid points for EROSION.
!     A max 'interior' square grid of 29X29 is assigned-no barriers
!     A max 'interior' rectangular grid of 59X59 is assigned barriers
!     to assign subregion index no. to each grid point.
!
!     +++ ARGUMENT DECLARATION +++
 
!
 
!     +++ LOCAL DEFINITIONS +++
!     imax  - no. grid intervals in x-direction
!     jmax  - no. grid intervals in y-direction.
!     ix    - grid interval in x-direction (m)
!     jy    - grid interval in y-direction (m)
!     dxmin - minimum grid interval (m)
!     csr   - current subr. index at grid point i,j
!     icsr  - same as csr but not an array
!     i,j   - do loop indexes
!
!     + + + GLOBAL COMMON BLOCKS + + +
 
!
!     + + + LOCAL COMMON BLOCKS + + +
!
!     +++ LOCAL VARIABLES +++
!
!     +++ END SPECIFICATIONS +++
!
!     set min grid spacing
      dxmin = min_grid_sp
!     set max no. of grid points with no barrier
 
      ngdpt = n_g_dpt
!     barriers?
      if (nbr>0) then
!        find shortest barrier to determine dxmin
         do i = 1,nbr
            if (amzbr(i)>0.0) dxmin = min(dxmin,5.0*amzbr(i))
                                        !Check for zero height barriers
         end do
         ngdpt = b_g_dpt  !default to this value if a barrier exists
      end if
 
!     calculate max grid intervals
!       calc. lx and ly sides of field
      lx = amxsim(1,2) - amxsim(1,1)
      ly = amxsim(2,2) - amxsim(2,1)
!
!^^^tmp out
!      write(*,*) 'tmp out from sbgrid, line 69'
!     write (*,*) 'lx=', lx, 'ly=',ly
! ^^^end tmp
!
!     increase grid points on large field
!      if((lx .gt. 200) .or. (ly  .gt. 200)) then
!         ngdpt = B_G_DPT
!      endif
!
!        case where lx > ly
      if (lx>ly) then
         imax = int(lx/dxmin)
         imax = min(imax,ngdpt)
         imax = max(imax,2)
!     calculate spacing for square or with barriers a rectangular grid
         ix = lx/(imax-1)
 
         if (nbr>0) then
            jmax = int(ly/dxmin)
            jmax = min(jmax,ngdpt)
         else
            jmax = anint(ly/ix) + 1
         end if
 
         jmax = max(jmax,2)
         jy = ly/(jmax-1)
 
!        case where lx = ly or lx < ly
      else
         jmax = int(ly/dxmin)
         jmax = min(jmax,ngdpt)
         jmax = max(jmax,2)
         jy = ly/(jmax-1)
 
         if (nbr>0) then
            imax = int(lx/dxmin)
            imax = min(imax,ngdpt)
         else
            imax = anint(lx/jy) + 1
         end if
         imax = max(imax,2)
         ix = lx/(imax-1)
 
      end if
!
!     determine subregion of each grid point
!     for a single subregion now
      icsr = 1
      do j = 0,jmax
         do i = 0,imax
!
!     for multiple subregions
!     do 5  icsr = 1, nsubr
!     for multiple subregions
!     if (i*ix .lt. amxsr(1,1,icsr) .or. i*ix .gt. amxsr(1,2,icsr))
!    & go to 10
!     if (j*jy .lt. amxsr(2,1,icsr) .or. j*jy .gt. amxsr(2,2,icsr))
!    & go to 10
!   5 continue
!
            csr(i,j) = icsr
         end do
      end do
      end subroutine sbgrid