c$Author: fredfox $
c$Date: 2001-12-06 00:22:15 $
c$Revision: 1.9 $
c$Source: /weru/cvs/weps/weps.src/mproc/remove.for,v $

      subroutine remove
     i           (pyieldf,pstalkf,prootf,rstandf,rflatf,rbgf,rrootf,
     m                  cstemht, cstems, cyld, cstand,
     m                  croot, crootz,
     m                  tstemht, tstems, tyld, tstand, tflat,
     m                  tbg, tbgz, troot, trootz,
     m                  dstemht, dstems, dstand, dflat,
     m                  dbgz, drootz,
     i                  nlay, tot_mass_rem, sel_mass_left)

c     + + + PURPOSE + + +
c     This subroutine performs the biomass manipulation of removing
c     biomass.  Removal can be from one or more of the following
c     processes: harvest, graze, burns and in the harvest process
c     the following can be removed: grain, crop, or all.  The amount of
c     each component removed is determined by the fraction passed into
c     this subroutine for each crop component.

c     Possible future enhancements
c     a)  bioflg - selects which age pools will be processed.  Probably the
c     same definition as other biomass manipulation process effects use.
c     b)  xxlocflg - selects the individual mass component pools that are
c     being effected (material being removed in this case).  There would
c     likely need to be more than one of these flags, possibly one for each
c     "age" pool.  Example settings could be:
c   
c     crlocflg (st,yld,flt,bg,rt)	decomp1locflg (st,flt,bg,rt)
c   
c     bit             val                    bit              val
c     x    st+yld+flt  0                      x    st+flt      0
c     0    yld*fract   1                      0    -           1
c     1    st*fract    2                      1    st*fract    2
c     2    fl*fract    4                      2    fl*fract    4
c     3    bg*fract    8                      3    bg*fract    8
c     4    rt*fract    16                     4    rt*fract    16
c     5    st*cutht    32                     5    st*cutht    32
c   
c   2. It appears likely to me that many harvest operations may require
c   multiple "remove" effect calls to adequately simulate their effects.
c   This would help simplify the logic within the "remove" effect, but
c   will require more code in the user interface to interact with the
c   user and to display all the "remove" effects with their appropriate
c   parameters and values. 

c     + + + KEYWORDS + + +
c     remove, biomass manipulation

c     + + + COMMON BLOCKS + + +
      include 'p1werm.inc'
      include 's1layr.inc'

c     + + + ARGUMENT DECLARATIONS + + +
      real   pyieldf,pstalkf,prootf,rstandf,rflatf,rbgf,rrootf
      real   cstemht, cstems, cyld, cstand
      real   croot, crootz(mnsz)
      real   tstemht, tstems, tyld, tstand, tflat
      real   tbg, tbgz(mnsz), troot, trootz(mnsz)
      real   dstemht(mnbpls), dstems(mnbpls)
      real   dstand(mnbpls), dflat(mnbpls)
      real   dbgz(mnsz,mnbpls), drootz(mnsz,mnbpls)
      integer nlay
      real   tot_mass_rem, sel_mass_left

c     + + + ARGUMENT DEFINITIONS + + +
c     pyieldf   - fraction of yield (reproductive components) removed (kg/kg)
c     pstalkf   - fraction of standing crop removed (kg/kg)
c     prootf    - fraction of plant roots removed (kg/kg)
c     rstandf   - fraction of standing residue removed (kg/kg)
c     rflatf    - fraction of flat residue removed (kg/kg)
c     rbgf      - fraction of below ground biomass material removed (kg/kg)
c     rrootf    - fraction of residue roots removed (kg/kg)

c     cstemht   - current crop height (m)
c     cstems    - crop population (# stems/m^2)
c     cyld      - (crop pool) yield biomass (kg/m^2)
c     cstand    - (crop pool) standing biomass (kg/m^2)
c     croot     - (crop pool) total root biomass (kg/m^2)
c     crootz    - (crop pool) root biomass by layer (kg/m^2)

c     tstemht   - (temporary crop pool) current (dead) crop height (m)
c     tstems    - (temporary crop pool) stem population (# stems/m^2)
c     tyld      - (temporary crop pool) yield biomass (kg/m^2)
c     tstand    - (temporary crop pool) standing biomass (kg/m^2)
c     tflat     - (temporary crop pool) flat biomass (kg/m^2)
c     tbg       - (temporary crop pool) total below ground biomass (kg/m^2)
c     tbgz      - (temporary crop pool) below ground biomass by layer (kg/m^2)
c     troot     - (temporary crop pool) total root biomass (kg/m^2)
c     trootz    - (temporary crop pool) root biomass by layer (kg/m^2)

c     dstemht   - (decomp pool) standing residue height (m)
c     dstems    - (decomp pool) number of standing residue stems (#/m^2)
c     dstand    - (decomp pool) standing biomass by age pool (kg/m^2)
c     dflat     - (decomp pool) surface biomass by age pool  (kg/m^2)
c     dbgz      - (decomp pool) buried biomass by age pool & layer (kg/m^2)
c     drootz    - (decomp pool) root biomass by age pool & layer (kg/m^2)

c     nlay      - number of layer from which below ground biomass is removed
c     tot_mass_rem - mass of material removed by this harvest operation (kg/m^2)
c     sel_mass_left - mass of material left in pools from which mass is removed
c                     by this harvest operation (kg/m^2)

c     + + + ACCESSED COMMON BLOCK VARIABLE DEFINITIONS + + +
c     mnbpls        - max number of decomposition pools (currently=3)
c     mnsz          - max number of soil layers
c     nslay         - number of soil layers used in a subregion

c     + + + PARAMETERS + + +

c     + + + LOCAL VARIABLES + + +
      character cropname*80
      integer  i,lay, pool_flag
      real mass_rem

c     + + + LOCAL VARIABLE DEFINITIONS + + +
c     i         - loop variable for decomp pools (3 pools total)
c     lay       - number of layers in a specified subregion
c     nlay      - number of layers in the tillage zone

c     + + + END SPECIFICATIONS + + +

      tot_mass_rem = 0.0
      sel_mass_left = 0.0
      pool_flag = 0

c    (note that all pools are updated)
c*****************************************************************
c     crop pool
c     Remove fraction (yldfr) of crop yield mass
      mass_rem = cyld * pyieldf
      if( mass_rem.gt.0.0 ) then 
          pool_flag = 1
          cyld = cyld - mass_rem
          tot_mass_rem = tot_mass_rem + mass_rem
      end if

c     Determine the amount of standing crop mass to be removed.
      mass_rem = cstand * pstalkf
      if( mass_rem.gt.0.0 ) then 
          pool_flag = 1
          cstand = cstand - mass_rem
          tot_mass_rem = tot_mass_rem + mass_rem
      end if

c     Remove fraction of crop root mass from each soil layer
      do lay=1,nlay
          crootz(lay) = crootz(lay) * (1.0 - prootf)
      end do
c     this is the total of root mass removed and we total mass from here
c     since there is no discrimination by layer at this point
      mass_rem = croot * prootf
      if( mass_rem.gt.0.0 ) then 
          pool_flag = 1
          croot = croot - mass_rem
          tot_mass_rem = tot_mass_rem + mass_rem
      end if

c     add biomass to selected mass if biomass was removed from pool
      if( pool_flag.eq.1 ) then
          sel_mass_left = sel_mass_left + cyld + cstand
          pool_flag = 0
      end if

c******************************************************************
c     temporary pool
c     Remove fraction (yldfr) of temporary yield mass
      mass_rem = tyld * pyieldf
      if( mass_rem.gt.0.0 ) then 
          pool_flag = 1
          tyld = tyld - mass_rem
          tot_mass_rem = tot_mass_rem + mass_rem
      end if

c     Determine the amount of standing crop mass to be removed.
      mass_rem = tstand * pstalkf
      if( mass_rem.gt.0.0 ) then 
          pool_flag = 1
          tstand = tstand - mass_rem
          tot_mass_rem = tot_mass_rem + mass_rem
      end if

c     update crop root mass from each soil layer
      do lay=1,nlay
          trootz(lay) = trootz(lay) * (1.0 - prootf)
      end do
c     this is the total of root mass removed and we total mass from here
c     since there is no discrimination by layer at this point
      mass_rem = troot * prootf
      if( mass_rem.gt.0.0 ) then 
          pool_flag = 1
          troot = troot - mass_rem
          tot_mass_rem = tot_mass_rem + mass_rem
      end if

c     Remove fraction of all below ground biomass in temporary pool
      mass_rem = tbg * rbgf
      if( mass_rem.gt.0.0 ) then 
          pool_flag = 1
          tbg = tbg - mass_rem
          tot_mass_rem = tot_mass_rem + mass_rem
      end if

c     update temporary below ground by layer
      do lay=1,nlay
          tbgz(lay) = tbgz(lay) * (1.0 - rbgf)
      end do

c     add biomass to selected mass if biomass was removed from pool
      if( pool_flag.eq.1 ) then
          sel_mass_left = sel_mass_left + tyld + tstand + tflat
          pool_flag = 0
      end if

c****************************************************************
c     residue pools
c     Determine the amount of standing residue mass to be removed.
c    (note that all decomp pools are checked and updated)
      do i=1,mnbpls
          mass_rem = dstems(i) * rstandf
          dstems(i) = dstems(i) - mass_rem
          tot_mass_rem = tot_mass_rem + mass_rem
          mass_rem = dstand(i) * rstandf
          dstand(i) = dstand(i) - mass_rem
          tot_mass_rem = tot_mass_rem + mass_rem
      end do

c     Remove fraction of all flat residue (temporary and decomp pools)
      mass_rem = tflat * rflatf
      tflat = tflat - mass_rem
      tot_mass_rem = tot_mass_rem + mass_rem
      do i=1,(mnbpls)
          mass_rem = dflat(i) * rflatf
          dflat(i) = dflat(i) - mass_rem
          tot_mass_rem = tot_mass_rem + mass_rem
      end do

c     Remove fraction of all below ground residue
c     and residue root fraction by layer
c     take this not the stuff below by layer
      do lay=1,nlay
          do i=1,mnbpls
              mass_rem = dbgz(lay,i) * rbgf
              dbgz(lay,i) = dbgz(lay,i) - mass_rem
              tot_mass_rem = tot_mass_rem + mass_rem
              mass_rem = drootz(lay,i) - mass_rem
              drootz(lay,i) = drootz(lay,i) * (1.0 - rrootf)
              tot_mass_rem = tot_mass_rem + mass_rem
          end do
      end do

      end