subroutine crop_read_input(bn1, bn2, bn3, bnu1, bnu2, bp1, bp2, bp3, bpu1, bpu2, cvm)
!**** PURPOSE: THIS SUBROUTINE READS CROP PARAMETERS from crop.dat
!**** CALLED IN: MAIN
!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
! PARAMETERS & VARIABLES
!
! >>>>> COMMON PARAMETERS & VARIABLES
! ald1(ic) = shape parameter for the LAI developement equation
! ald2(ic) = shape parameter for the LAI developement equation
! almn(ic) = minimum Leaf Area Index (for forest and natural vegetation)
! be(ic) = biomass-energy ratio
! blai(ic) = maximum Leaf Area Index for crop
! bn1(ic) = nitrogen uptake parameter #1: normal fraction of N
! in crop biomass at emergence, kg N/kg biomass
! bn2(ic) = nitrogen uptake parameter #2: normal fraction of N
! in crop biomass at 0.5 maturity, kg N/kg biomass
! bn3(ic) = nitrogen uptake parameter #3: normal fraction of N
! in crop biomass at maturity, kg N/kg biomass
! bnu1(ic) = shape parameter to calc optimal N fraction in crop biomass
! bnu2(ic) = shape parameter to calc optimal N fraction in crop biomass
! bp1(ic) = phosphorus uptake parameter #1: normal fraction of P
! in crop biomass at emergence, kg P/kg biomass
! bp2(ic) = phosphorus uptake parameter #2: normal fraction of P
! in crop biomass at 0.5 maturity, kg P/kg biomass
! bp3(ic) = phosphorus uptake parameter #3: normal fraction of P
! P in crop biomass at maturity, kg P/kg biomass
! bpu1(ic) = shape parameter to calc optimal P fraction in crop biomass
! bpu2(ic) = shape parameter to calc optimal P fraction in crop biomass
! cnyld(ic) = fraction of nitrogen in yield, kg N/kg yield
! cpyld(ic) = fraction of phosphorus in yield, kg P/kg yield
! cvm(ic) = minimum value of C factor for water erosion
! dlai(ic) = fraction of growing season when leaf area declines
! dlp1 = complex number: fraction of grow. season, max corresp. LAI
! dlp2 = complex number: fraction of grow. season, max corresp. LAI
! hi(ic) = harvest index
! hun(ic) = potential heat units required for maturity of crop
! icnum(ic) = crop number
! ilcc(ic) = land cover category:
! (1) annual crop (row crop / small grain)
! (2) annual winter crop (row crop / small grain)
! (3) perennial (grass, brush, urban, water)
! (4) woods
! (5) annual legumes (row crop)
! (6) annual winter legumes (row crop)
! (7) perennial legumes (grass)
! pt2(ic) = 2nd point on radiation use efficiency curve:
! complex number: The value to the left of (.) is a CO2
! atm. concentration higher than the ambient (in units of
! microliters CO2/liter air). The value to the right of (.)
! is the corresponding biomass-energy ratio divided by 100
! rdmx(ic) = maximum root depth (m), then converted to mm
! sla(ic) = specific leaf area, m2/kg
! tb(ic) = base temperature for plan growth, degrees C
! to(ic) = optimal temperature for plant growth, degrees C
! >>>>>
! >>>>> STATIC PARAMETERS
! b1 = intermediate parameter
! b2 = intermediate parameter
! b3 = intermediate parameter
! ic = count parameter
! titl = text
! >>>>>
!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
!**** Include common parameters
use vegetation, only : &
vegetation_s_curve_parameters, &
ald1, &
ald2, &
blai, &
almn, &
ilcc, &
sla, &
to, &
tb, &
dlai
use input, only : &
read_integer_column, &
read_string_column, &
read_real_column, &
input_count_rows
real(dp), dimension(:), intent(inout) :: bn1
real(dp), dimension(:), intent(inout) :: bn2
real(dp), dimension(:), intent(inout) :: bn3
real(dp), dimension(:), intent(inout) :: bnu1
real(dp), dimension(:), intent(inout) :: bnu2
real(dp), dimension(:), intent(inout) :: bp1
real(dp), dimension(:), intent(inout) :: bp2
real(dp), dimension(:), intent(inout) :: bp3
real(dp), dimension(:), intent(inout) :: bpu1
real(dp), dimension(:), intent(inout) :: bpu2
real(dp), dimension(:), intent(inout) :: cvm
real(dp), dimension(:), allocatable :: dlp1_array
real(dp), dimension(:), allocatable :: dlp2_array
integer ic
real(dp) b1, b2, b3
real(dp) x3, x4
x3 = .5
x4 = 1.
allocate(dlp1_array(mcrdb))
allocate(dlp2_array(mcrdb))
call read_integer_column(crop_input_file_id, "icnum", icnum, 0)
call read_string_column(crop_input_file_id, "cpnm", cnam, "")
call read_real_column(crop_input_file_id, "be", be, 0.0_dp)
call read_real_column(crop_input_file_id, "hi", hi, 0.0_dp)
call read_real_column(crop_input_file_id, "blai", blai, 0.0_dp)
call read_real_column(crop_input_file_id, "dlp1", dlp1_array, 0.0_dp)
call read_real_column(crop_input_file_id, "dlp2", dlp2_array, 0.0_dp)
call read_real_column(crop_input_file_id, "bn1", bn1, 0.0_dp)
call read_real_column(crop_input_file_id, "bn2", bn2, 0.0_dp)
call read_real_column(crop_input_file_id, "bn3", bn3, 0.0_dp)
call read_real_column(crop_input_file_id, "bp1", bp1, 0.0_dp)
call read_real_column(crop_input_file_id, "bp2", bp2, 0.0_dp)
call read_real_column(crop_input_file_id, "bp3", bp3, 0.0_dp)
call read_real_column(crop_input_file_id, "cnyld", cnyld, 0.0_dp)
call read_real_column(crop_input_file_id, "cpyld", cpyld, 0.0_dp)
call read_real_column(crop_input_file_id, "rdmx", rdmx, 0.0_dp)
call read_real_column(crop_input_file_id, "cvm", cvm, 0.0_dp)
call read_real_column(crop_input_file_id, "pt2", pt2, 0.0_dp)
call read_integer_column(crop_input_file_id, "hun", hun, 0)
! Vegetation variables
call read_real_column(crop_input_file_id, "almn", almn, 0.0_dp)
call read_integer_column(crop_input_file_id, "ird", ilcc, 0)
call read_real_column(crop_input_file_id, "sla", sla, 0.0_dp)
call read_real_column(crop_input_file_id, "to", to, 0.0_dp)
call read_real_column(crop_input_file_id, "tb", tb, 0.0_dp)
call read_real_column(crop_input_file_id, "dlai", dlai, 0.0_dp)
!**** READ 5 - crop parameter database (unit 5-crop.dat)
do ic = 1, mcrdb
dlp1 = dlp1_array(ic)
dlp2 = dlp2_array(ic)
cvm(ic) =log(cvm(ic))
rdmx(ic) = rdmx(ic) * 1000.
!**** RECALCULATE: dlp1, dlp2 need to be broken into two variables
!**** Isolate number to left of decimal: b1, b2 - fraction of the grow. season
b1 = .01 * Int(dlp1)
b2 = .01 * Int(dlp2)
!**** Isolate number to right of decimal, dlp1, 2 is now max LAI corr. to b1, 2
dlp1 = dlp1 - Int(dlp1)
dlp2 = dlp2 - Int(dlp2)
!#### determine shape parameters for the LAI developement equation
call vegetation_s_curve_parameters(dlp1, dlp2, b1, b2, ald1(ic), ald2(ic))
!**** Calc nitrogen uptake parameters
b1 = bn1(ic) - bn3(ic)
b2 = 1. - (bn2(ic) - bn3(ic)) / b1
b3 = 1. - .00001 / b1
!#### CALL ASCRV
call vegetation_s_curve_parameters(b2, b3, x3, x4, bnu1(ic), bnu2(ic))
!**** CORRECT bp3 & Calc phosphorus uptake parameters
if (bp2(ic) - bp3(ic) < .0001) bp3(ic) = .75 * bp3(ic)
b1 = bp1(ic) - bp3(ic)
b2 = 1. - (bp2(ic) - bp3(ic)) / b1
b3 = 1. - .00001 / b1
!#### CALL ASCRVg
if (b2 > 0. .AND. b2 < 1.) then
call vegetation_s_curve_parameters(b2, b3, x3, x4, bpu1(ic), bpu2(ic))
else
call log_error("crop_read_input", &
"Value of b2 in call vegetation_s_curve_parameters() is causing an error "// &
"(ic, b2, b1, bp2, bp3):", int=ic, reals=(/b2, b1, bp2, bp3/))
end if
end do
deallocate(dlp1_array)
deallocate(dlp2_array)
call log_debug("crop_read_input", 'Crop parameters READ!')
end subroutine crop_read_input