subroutine crop_operation(j, je, k, alai, dm, frar, g, ida, iy, olai, rd, rsd, ws)
!**** PURPOSE: TO DEFINE PLANT OPERATIONS
!**** CALLED IN: CRPMD
!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
! PARAMETERS & VARIABLES
!
! >>>>> COMMON PARAMETERS & VARIABLES
! alai(j, je) = leaf area index
! aryld(j, k, icr) = fraction of area by crop per sub, soil
! arylda(iy, icr) = fraction of area by crop per year
! aryldc(icr) = fraction of area by crop
! arylds(k, icr) = fraction of area by crop per soil
! avyld(j, k, icr) = av. yld per sub, soil, crop, kg/ha
! avylda(iy, icr) = av yld per year, crop, kg/ha
! avyldc(icr) = av yld per crop, kg/ha
! avylds(k, icr) = av. yld per soil, crop, kg/ha
! cva(j, je) = vegetation cover, kg/ha
! dm(j, je) = total biomass, kg/ha
! fon(j, je, l) = fresh org N, kg/ha
! fop(j, je, l) = fresh org P, kg/ha
! frar(j, je) = fractional area of hydrotope in subbasin
! g(j, je) = fraction of heat units to maturity accumulated
! hi(icr) = harvest index for crop (database), for maize & potat.
! hia(j, je) = harvest index
! hiad(j, je) = harvest index, adjusted
! huharv(j, je) = harvest index heat unit
! icc = index for cover crop corr. number in crop database
! ida = current day
! idayx = par = ida, to calc ndgro - number of growth days
! idop(5, iop) = day of operation
! igro(j, je) = vegetation index, =1 if yes
! iopc(5, iop) = opeartion code: 1 - planting, ...
! ipo = index for potatoes corr. number in crop database
! isba = index for s. barley corr. number in crop databas
! istyr = starting year
! iwb = index for w. barley corr. number in crop databas
! iwr = index for w. rye corr. number in crop databas
! iww = index for w. wheat corr. number in crop databas
! iy = current year as counter (1, ..., nbyr)
! ncrp(iop) = crop number
! ndgro = number of growth days
! ndpri = day to write crop yield for GIS output
! nucr(j, je) = crop number (database)
! olai(j, je) = alai(j, je) - leaf area index
! rsd(j, je, 2) = residue, kg/ha
! rwt(j, je) = fraction of root weight
! sbar(j) = subbasin area, m2
! snup(j, je) = N uptake, kg/ha
! spup(j, je) = P uptake, kg/ha
! swh(j, je) = actual transp. by plants, mm
! swp(j, je) = potent. transp. by plants, mm
! ws(j, je) = water stress
! yld(j, je) = crop yield, kg/ha
! ylda(j, k) = crop yield for subbasin and soil, kg/ha
! >>>>>
! >>>>> STATIC PARAMETERS
! icr = local par
! ii = local par
! ioper = local par
! xx = local par
! yield = current yield
! >>>>>
!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
use output, only: output_store_hydrotope_value
!**** Include common parameters
real(dp), dimension(:, :), intent(inout) :: alai
real(dp), dimension(:, :), intent(inout) :: dm
real(dp), dimension(:, :), intent(in) :: frar
real(dp), dimension(:, :), intent(inout) :: g
integer, intent(in) :: ida
integer, intent(in) :: iy
real(dp), dimension(:, :), intent(out) :: olai
real(dp), dimension(:, :), intent(inout) :: rd
real(dp), dimension(:, :, :), intent(inout) :: rsd
real(dp), dimension(:, :), intent(inout) :: ws
integer j, je, k
integer icr, ii, ioper
real(dp) xx, yield
!*********************************************************** START IF (IGRO=0)
!**** CHECK if day of planting, then goto 10 - planting
if (igro(j, je) .eq. 0) then
do ii = 1, mop
if (ida .eq. idop(ii)) then
nucr(j, je) = ncrp(ii)
ioper = iopc(ii)
if (ioper .eq. 1) then
igro(j, je) = 1
g(j, je) = 0.
dm(j, je) = 0.01
snup(j, je) = 0.
spup(j, je) = 0.
swh(j, je) = 0.
swp(j, je) = 0.
huharv(j, je) = 0.
hia(j, je) = 0.
olai(j, je) = 0.
rwt(j, je) = 0.
idayx = - 99
ndgro = 0
EXIT
end if
endif
end do
endif
!*********************************************************** END IF (IGRO=0)
!*********************************************************** START IF(IGRO=1)
!**** CHECK if day of harvest and kill
if (igro(j, je) .eq. 1) then
do ii = 1, mop
if (ida .eq. idop(ii)) then
nucr(j, je) = ncrp(ii)
ioper = iopc(ii)
if (ioper .eq. 2) then
!**** CALC HARVEST AND KILL
igro(j, je) = 0
if (hiad(j, je) .gt. hi(nucr(j, je))) hiad(j, je) = hi(nucr(j, je))
if (hiad(j, je) .gt. 1.) hiad(j, je) = 1.
!**** CALC residue & fresh org N and P (no residue, fon, fop for cover crop)
if (nucr(j, je) .ne. icc) then
rsd(j, je, 1) = 0.25 * (1. - rwt(j, je)) * (1. - hiad(j, je)) * dm(j, je) &
+ rsd(j, je, 1)
rsd(j, je, 2) = 0.75 * (1. - rwt(j, je)) * (1. - hiad(j, je)) * dm(j, je) &
+ rsd(j, je, 2)
if (rsd(j, je, 1) .le. 0.) rsd(j, je, 1) = 1.e-6
if (rsd(j, je, 2) .le. 0.) rsd(j, je, 2) = 1.e-6
fon(j, je, 1) = rsd(j, je, 1) * .008
fon(j, je, 2) = rsd(j, je, 2) * .008
fop(j, je, 1) = rsd(j, je, 1) * .0011
fop(j, je, 2) = rsd(j, je, 2) * .0011
endif
!**** CALC yield
! ATTN! Harvest index used for grains: hia(), for maize, potatoes: hi()
! No yield for cover crop (icc)
if (nucr(j, je) .eq. iww .or. nucr(j, je) .eq. iwb .or. &
nucr(j, je) .eq. iwr .or. nucr(j, je) .eq. isba) then
yield = 0.85 * dm(j, je) * hia(j, je)
else if (nucr(j, je) .eq. ipo) then
yield = 1.00 * dm(j, je) * hi(nucr(j, je))
else
yield = 0.85 * dm(j, je) * hi(nucr(j, je))
endif
if (nucr(j, je) .eq. icc) yield = 0.
ylda(j, k) = yield
yld(j, je) = yld(j, je) + yield
call output_store_hydrotope_value(crop_yield_output_id, j, je, yield)
dm(j, je) = 0.
rd(j, je) = 0.
ws(j, je) = 1.
alai(j, je) = 0.
hia(j, je) = 0.
cva(j, je) = rsd(j, je, 1) + rsd(j, je, 2)
g(j, je) = 0.
idayx = 0
icr = nucr(j, je)
!**** CALC Day to write crop yield for GIS output (except cover crop)
if (icr .ne. icc) ndpri = ida + 3
!**** CALC average yield
! avyld(j, k, icr) av yld per sub, soil, crop & aryld(j, k, icr): frac. area
! avylds(k, icr) av yld per soil, crop & arylds(k, icr): frac. area
! avyldc(icr) av yld per crop & aryldc(icr): frac. area
! avylda(iy, icr) av yld per year, crop & arylda(iy, icr): frac. area
if (icrop == 1) then
avyld(j, k, icr) = avyld(j, k, icr) + ylda(j, k) * frar(j, je) / 100.
aryld(j, k, icr) = aryld(j, k, icr) + frar(j, je)
endif
avylds(k, icr) = avylds(k, icr) + ylda(j, k) * frar(j, je) / 100.
arylds(k, icr) = arylds(k, icr) + frar(j, je)
avyldc(icr) = avyldc(icr) + ylda(j, k) * frar(j, je) / 100.
aryldc(icr) = aryldc(icr) + frar(j, je)
avylda(iy, icr) = avylda(iy, icr) + ylda(j, k) * frar(j, je) / 100.
arylda(iy, icr) = arylda(iy, icr) + frar(j, je)
return
!**** END CALC HARVEST AND KILL
end if
if (ioper .eq. 3) then
!**** CALC HARVEST ONLY - CUTTING, NO KILL
if (hiad(j, je) .gt. hi(nucr(j, je))) &
hiad(j, je) = hi(nucr(j, je))
yield = (1. - rwt(j, je)) * dm(j, je) * hiad(j, je)
yld(j, je) = yld(j, je) + yield
xx = dm(j, je)
dm(j, je) = dm(j, je) - yield
alai(j, je) = alai(j, je) * dm(j, je) / xx
g(j, je) = g(j, je) * dm(j, je) / xx
return
!**** END CALC HARVEST ONLY
end if
if (ioper .eq. 4) then
!**** CALC KILL
igro(j, je) = 0
dm(j, je) = 0.
ws(j, je) = 1.
alai(j, je) = 0.
cva(j, je) = 0.
g(j, je) = 0.
rd = 0.
!**** END CALC KILL
end if
endif
end do
endif
!*********************************************************** END IF (IGRO=1)
return
end subroutine crop_operation