c$Header: /weru/cvs/weps/weps.src/manage/setbds.for,v 1.1.1.1 1999-03-12 17:05:28 wagner Exp $

c
c
      real function setbds 
     i                    (clay, sand, om)
c
c     + + + PURPOSE + + +
c     The following function estimates settled soil bulk density from
c     intrinsic properties. see Rawls (1983) Soil Science 135, 123-125.
c
c     Should eventually be called by MAIN to initialize the values
c     for each subregion (unless soil composition changes).
c
c
c     + + + KEYWORDS + + +
c     bulk density, initialization
c
c     + + + PARAMETERS AND COMMON BLOCKS + + +
c
c     none at this time
c
c     + + + ARGUMENT DECLARATIONS + + +
c
      real clay, sand, om
c
c     + + + ARGUMENT DEFINITIONS + + +
c
c     setbd - settled bulk density
c     clay  - fraction of soil clay content
c     sand  - fraction of soil sand content
c     om    - organic matter
c
c     + + + LOCAL VARIABLES + + +
c
      integer i, j
      real mbdtv (0:10,0:10), mbd
c     real mbd1, mbd2
c
c     mbd   - mineral bulk density without organic matter
c
c     + + + SUBROUTINES CALLED + + +
c
c     + + + FUNCTION DECLARATONS + + +
c
c     + + + DATA INITIALIZATIONS + + +
c
      data mbdtv /1.48,1.25,1.00,1.06,1.16,1.22,1.30,1.39,1.45,1.51,1.52
     &           ,1.52,1.40,1.19,1.25,1.32,1.40,1.52,1.58,1.63,1.65,0.
     &           ,1.52,1.40,1.25,1.35,1.45,1.53,1.60,1.64,1.72,0.,0.
     &           ,1.52,1.40,1.29,1.41,1.50,1.57,1.63,1.68,0.,0.,0.
     &           ,1.50,1.40,1.35,1.43,1.53,1.61,1.64,0.,0.,0.,0.
     &           ,1.46,1.40,1.40,1.43,1.53,1.62,0.,0.,0.,0.,0.
     &           ,1.45,1.40,1.38,1.42,1.50,0.,0.,0.,0.,0.,0.
     &           ,1.42,1.37,1.33,1.33,0.,0.,0.,0.,0.,0.,0.
     &           ,1.33,1.32,1.20,0.,0.,0.,0.,0.,0.,0.,0.
     &           ,1.23,1.18,0.,0.,0.,0.,0.,0.,0.,0.,0.
     &           ,1.15,0.,0.,0.,0.,0.,0.,0.,0.,0.,0./
c
c     + + + END SPECIFICATIONS + + +
c

      j = nint(clay*100.0/10.0)
      i = nint(sand*100.0/10.0)

c     i = 1 + ( ( clay*100.0 ) - 5 ) / 10
c     j = 1 + ( ( sand*100.0 ) - 5 ) / 10

      mbd = mbdtv(i,j)

c     mbd1  = mbdtv (i,j) + ( mbdtv (i, j+1)
c    &        - mbdtv(i,j) ) * ( (sand*100.0) - 10.0 * (j-i) ) 
c    &        / 10.0
c     mbd2  = mbdtv (i+1,j) + ( mbdtv (i+1, j+1)
c    &         - mbdtv(i+1,j) ) * ( (sand*100.0) -   10.0 * (j-i) )
c    &        / 10.0

c     mbd   = mbd1 + (mbd2-mbd1) * ( (clay * 100.0)- 10.0 * (i-1) )
c    &        /10.0

      setbds = 100.0/
     &        (((om*100.0) / 0.224)  + (100.0 -(om * 100.0) )/ mbd)

      return
      end