!$Author: wagner $ !$Date: 2007-09-13 02:11:59 $ !$Revision: 1.10 $ !$Source: /weru/cvs/weps/weps_ref/src/ref_friction_vel.f95,v $ !----------------------------------------------------------------------- !! Convert between _reference_ and _non-reference_ wind speeds !! !! Convert from a friction velocity to a _reference_ wind speed (10m ht @ 25mm aerodynamic roughness) !! !! Calculate the friction velocity for a variety of situations. !! - _reference_ friction velocity (10m ht @ 25mm aerodynamic roughness) !! - for a given wind speed under _reference_ conditions !! - for a given wind speed under _non-reference_ conditions !! - friction velocity for a bare surface (RR and/or oriented ridges) !! - for a given _reference_ wind speed !! - for a given _non-reference_ wind speed !! - friction velocity above a standing vegetative canopy !! - for a given _reference_ wind speed !! - for a given _non-reference_ wind speed !! - friction velocity at the surface (below a standing biomass canopy) !! - for a given _reference_ wind speed !! - for a given _non-reference_ wind speed !! - friction velocity at the surface (_with or without_ an existing biomass canopy) !! - for a given _reference_ wind speed !! - for a given _non-reference_ wind speed !! !! Calculate the aerodynamic roughness for a variety of situations. !! - for a randomly rough surface !! - for an oriented (ridge) roughened surface !! - above a standing canopy of vegetative biomass !! - at the surface below a standing canopy of vegetative biomass !! !! Calculate the _representative_ aerodynamic roughness for a given friction velocity MODULE Friction_Velocity IMPLICIT NONE PUBLIC :: WU_ref PRIVATE :: WU_to_WUref ! 3 args ! PRIVATE :: WUF_to_WUref! 1 arg (not sure if I want it here - LEW) PUBLIC :: WU_from_WUref PUBLIC :: WUF_to_WUref ! probably best here - LEW PUBLIC :: WUF_ref PRIVATE :: WUF_ref_0 ! 1 arg PRIVATE :: WUF_ref_1 ! 3 args PUBLIC :: WUS_surf PUBLIC :: WUS_b ! 2 args PUBLIC :: WUS_bv ! 4 args PUBLIC :: WUS_b_bv ! 2 args PUBLIC :: WUS_canopy PUBLIC :: WUS_v ! 4 args PUBLIC :: WZo_surf PUBLIC :: WZo_rr ! 1 arg PUBLIC :: WZo_rg ! 2 args PUBLIC :: WZo_rgrr ! 3 args PUBLIC :: WZo_bv ! 5 args PUBLIC :: WZo_canopy PUBLIC :: WZo_v PUBLIC :: WZo_rep PUBLIC :: BRcd PRIVATE :: BRcd_0 ! 2 args PRIVATE :: BRcd_1 ! 6 args PRIVATE :: BRcd_2 ! 3 args (uses WUS_roughness, WUS_crop and WUS_residue structures) PUBLIC :: BRsai PUBLIC :: BRlai PUBLIC :: Rcrow PUBLIC :: Rrg PUBLIC :: n180_angle PUBLIC :: rel_angle !! This interface allows one to compute a _reference_ wind speed !! when given a wind speed at a specified anemometer height and aerodynamic roughness INTERFACE WU_ref MODULE PROCEDURE WU_to_WUref !! wind speed (m/s) at _reference_ conditions END INTERFACE WU_ref !! This interface allows one to compute a _reference_ friction velocity !! when given a wind speed under _reference_ conditions, i.e: !!> 10m anemometer height !!< 25mm aerodynamic roughness !! or when given a wind speed at a specified anemometer height and aerodynamic roughness INTERFACE WUF_ref MODULE PROCEDURE WUF_ref_0 !! wind speed (m/s) at _reference_ conditions MODULE PROCEDURE WUF_ref_1 !! wind speed (m/s) not at _reference_ conditions END INTERFACE WUF_ref !! This interface allows one to compute a friction velocity (wind speed) !! under _reference_ conditions, i.e: !!> 10m anemometer height !!< 25mm aerodynamic roughness INTERFACE WUF_to_WUref MODULE PROCEDURE WUF_to_WUref !! friction velocity wind speed (m/s) at _reference_ conditions END INTERFACE WUF_to_WUref !! Calculate the friction velocity (m/s) at _non-reference_ surface conditions !!> surface can be bare ground or above standing biomass canopy !! surface can be below a standing biomass canopy, or !!< surface can be with or without standing biomass canopy INTERFACE WUS_surf MODULE PROCEDURE WUS_b !! 2 params (bare ground surface) MODULE PROCEDURE WUS_bv !! 4 params (ground surface below a canopy) ! MODULE PROCEDURE WUS_b_bv END INTERFACE WUS_surf !! Calculate the friction velocity (m/s) at _non-reference_ surface conditions !!> above standing biomass canopy !!< INTERFACE WUS_canopy MODULE PROCEDURE WUS_v !! 4 param (Above standing canopy friction velocity) END INTERFACE WUS_canopy !! Calculate the surface aerodynamic roughness (mm) INTERFACE WZo_surf MODULE PROCEDURE WZo_rr !! 1 param (RR effects only) MODULE PROCEDURE WZo_rg !! 2 param (Oriented roughness effects only) MODULE PROCEDURE WZo_rgrr !! 3 param (Oriented and random roughness effects only) MODULE PROCEDURE WZo_bv !! 5 param (Standing canopy effects on surface WZo, which include RR and ridge effects) END INTERFACE WZo_surf INTERFACE WZo_canopy MODULE PROCEDURE WZo_v !! 5 param (Above standing canopy WZo, which include RR and ridge effects) END INTERFACE WZo_canopy !! Calculate the "effective" biomass drag coefficient due to standing biomass INTERFACE BRcd MODULE PROCEDURE BRcd_0 !! 2 params (_effective_ lai and sai which have already combined the crop and residue pools) MODULE PROCEDURE BRcd_1 !! 6 params (individual standing crop and residue component inputs) MODULE PROCEDURE BRcd_2 !! 3 params (Uses WUS_* structures) MODULE PROCEDURE BRcd_3 !! 8 params (individual standing crop and residue component inputs as well as ridge and crop ht) END INTERFACE BRcd TYPE, PUBLIC :: WUS_info REAL :: WU REAL :: WUF REAL :: WUS REAL :: anemom_ht REAL :: Zo END TYPE WUS_info TYPE, PUBLIC :: WUS_roughness REAL :: RR REAL :: r_ht REAL :: r_sp REAL :: r_dir REAL :: r_width END TYPE WUS_roughness TYPE, PUBLIC :: WUS_crop REAL :: lai REAL :: sai REAL :: ht REAL :: row_sp !! crop row spacing (zero implies broadcast seeding) INTEGER :: ridge_plant_flag !! crop planted in furrow (0) or on ridge (1) END TYPE WUS_crop TYPE, PUBLIC :: WUS_residue REAL :: lai REAL :: sai REAL :: ht END TYPE WUS_residue TYPE, PUBLIC :: WUS_biomass REAL :: biomass_lai REAL :: biomass_sai REAL :: biomass_ht END TYPE WUS_biomass PRIVATE REAL, PARAMETER :: WZo_ref = 25.0 !! "reference" aerodynamic roughness (mm) REAL, PARAMETER :: Anemht_ref = 10.0 !! "reference" anemometer height (m) REAL, PARAMETER :: m_to_mm = 1000.0 !! m to mm conversion factor (mm/m) REAL, PARAMETER :: k = 0.4 !! Von Karman's constant REAL, PARAMETER :: PI = 3.1415927 !! PI CONTAINS !! Calculate the _reference_ friction velocity (m/s) for a given wind speed @ !! standard _reference_ weather station conditions, e.g: !!> 10m anemometer height !!< 25mm aerodynamic roughness REAL FUNCTION WUF_ref_0 (wu_ref) RESULT(frict_vel) REAL, INTENT (IN) :: wu_ref !! wind velocity under _reference_ weather station conditions (m/s) frict_vel = k * wu_ref / log(Anemht_ref*m_to_mm/WZo_ref) END FUNCTION WUF_ref_0 !! Calculate the _reference_ friction velocity (m/s) for a given wind speed @ !! specified weather station conditions, e.g: !!> anemometer height (m) !!< aerodynamic roughness (mm) REAL FUNCTION WUF_ref_1 (wu, anemht, wzo) REAL, INTENT (IN) :: wu !! wind velocity at measurement location (m/s) REAL, INTENT (IN) :: anemht !! anemometer height at measurement location (m) REAL, INTENT (IN) :: wzo !! aerodynamic roughness at measurement location (mm) ! adjust wu to reference conditions, if required !WU_ref(wu, anemht, wzo) WUF_ref_1 = k * WU_ref(wu, anemht, wzo) / log(Anemht_ref*m_to_mm/WZo_ref) END FUNCTION WUF_ref_1 !! Calculate the _reference_ wind velocity (m/s) for !! standard weather station conditions e.g: !!> 10m anemometer height !!< 25mm aerodynamic roughness !! given a wind velocity at a specified anemometer height and aerodynamic roughness REAL FUNCTION WU_to_WUref (wu, anemht, wzo) REAL, INTENT (IN) :: wu !! wind velocity (m/s) at _non-standard_ weather station REAL, INTENT (IN) :: anemht !! anemometer height (m) at non-standard weather station REAL, INTENT (IN) :: wzo !! aerodynamic roughness (mm) at non-standard weather station ! adjust wu to reference conditions IF (wzo <= tiny(1.0)) THEN ! What does a wzo of zero mean? !write(*,*) 'wzo/tiny/huge', wzo, tiny(wzo),huge(wzo), 'log(huge(wzo))', log(huge(wzo)) WU_to_WUref = wu * log(Anemht_ref*m_to_mm/WZo_ref)/log(huge(wzo)) ELSE IF ((anemht*m_to_mm-wzo) > tiny(wzo)) THEN ! anemht must be greater than wzo WU_to_WUref = wu * log(Anemht_ref*m_to_mm/WZo_ref)/log(anemht*m_to_mm/wzo) ELSE ! Error exit if anemht is less than wzo write(0,*) 'Error: WU_to_ref() function' write(0,*) 'WU_to_ref(wu (m/s) =',wu,'anemht(m) =',anemht,'wzo(mm) =',wzo,')' write(0,*) 'anemht cannot be less than or equal to wzo' call exit(1) END IF END FUNCTION WU_to_WUref !! Calculate the _nonreference_ wind velocity (m/s) for !! specified non-standard weather station conditions e.g: !!> anemometer height (m) !!< aerodynamic roughness (mm) !! given a _reference_ wind velocity at standard weather station conditions e.g. !!> 10m anemometer height !!< 25mm aerodynamic roughness !! anemometer height and aerodynamic roughness REAL FUNCTION WU_from_WUref (wu_ref, anemht, wzo) REAL, INTENT (IN) :: wu_ref !! wind velocity (m/s) at _non-standard_ weather station REAL, INTENT (IN) :: anemht !! anemometer height (m) at non-standard weather station REAL, INTENT (IN) :: wzo !! aerodynamic roughness (mm) at non-standard weather station ! adjust wu to reference conditions IF (wu_ref == 0.0) THEN WU_from_WUref = 0.0 ELSE IF (wzo <= tiny(wzo)) THEN !write(*,*) 'wzo is small', wzo WU_from_WUref = wu_ref / (log(Anemht_ref*m_to_mm/WZo_ref)/log(huge(wzo))) ELSE WU_from_WUref = wu_ref / (log(Anemht_ref*m_to_mm/WZo_ref)/log(anemht*m_to_mm/wzo)) END IF END FUNCTION WU_from_WUref !! Calculate the _reference_ wind velocity (m/s) for !! standard weather station conditions e.g: !!> 10m anemometer height !!< 25mm aerodynamic roughness !! given a friction wind velocity REAL FUNCTION WUF_to_WUref (wuf) REAL, INTENT (IN) :: wuf !! friction velocity (m/s) WUF_to_WUref = wuf * log(Anemht_ref*m_to_mm/WZo_ref) / k END FUNCTION WUF_to_WUref !! Calculate the friction velocity (m/s) at non-reference surface conditions !! (surface can be bare ground or above standing biomass canopy) REAL FUNCTION WUS_b (wuf, wzo) REAL, INTENT (IN) :: wuf !! _reference_ friction velocity (m/s) REAL, INTENT (IN) :: wzo !! aerodynamic roughness at surface (mm) ! 0.06676 = 0.4 / ln(Anemht_ref/WZo_ref) (all units in mm) ! 0.4 is Von Karman constant ! "reference" anemometer height is 10m ! "reference" aerodynamic roughness is 25mm WUS_b = wuf * (wzo/WZo_ref)**0.06676 END FUNCTION WUS_b !! Calculate the friction velocity above standing biomass canopy (m/s) REAL FUNCTION WUS_v (wuf, bht, wzo, brcd) REAL, INTENT (IN) :: wuf !! _reference_ friction velocity (m/s) REAL, INTENT (IN) :: bht !! standing biomass height (mm) REAL, INTENT (IN) :: wzo !! aerodynamic roughness at bare ground surface (mm) REAL, INTENT (IN) :: brcd !! effective biomass drag coefficient ! "reference" anemometer height is 10m ! "reference" aerodynamic roughness is 25mm WUS_v = WUS_b(wuf, WZo_v(bht, wzo, brcd)) END FUNCTION WUS_v !! Calculate the friction velocity (m/s) at the surface (below a standing biomass canopy) REAL FUNCTION WUS_bv (wuf, bht, wzo, brcd) REAL, INTENT (IN) :: wuf !! _reference_ friction velocity (m/s) REAL, INTENT (IN) :: bht !! standing biomass height (mm) REAL, INTENT (IN) :: wzo !! aerodynamic roughness at bare ground surface (mm) REAL, INTENT (IN) :: brcd !! effective biomass drag coefficient IF (brcd .eq. 0.0 .or. bht .eq. 0.0) THEN !no cover WUS_bv = WUS_b(wuf, wzo) ELSE IF (brcd .gt. 25.6) THEN !check to avoid underflow WUS_bv = 0.0 ! write(0,*) 'Possible underflow condition (wu_star_bv, brcd)- ', & ! wus_bv, BRcd ELSE IF (brcd .gt. 2.56) THEN !check to avoid underflow WUS_bv = WUS_b(wuf, WZo_v (bht, wzo, brcd)) * (0.025 * exp(-brcd/0.356)) ! write(0,*) 'Possible underflow condition (wu_star_bv, brcd)- ', & ! WUS_bv, brcd ELSE WUS_bv = WUS_b(wuf, WZo_v (bht, wzo, brcd)) * ((0.86 * exp(-brcd/0.0298)) + & (0.025 * exp(-brcd/0.356))) END IF END FUNCTION WUS_bv !! Calculate the friction velocity (m/s) at the simulation surface !! (select minimum value based on surface conditions - between bare surface and below standing biomass canopy) REAL FUNCTION WUS_b_bv (wus_b, wus_bv) REAL, INTENT (IN) :: wus_b !! friction velocity (m/s) at simulation surface (bare surface) REAL, INTENT (IN) :: wus_bv !! friction velocity (m/s) at simulation surface (below standing biomass canopy) WUS_b_bv = min(wus_b, wus_bv) !??????? END FUNCTION WUS_b_bv !! Calculate the _representative_ aerodynamic roughness (mm) !! given a friction velocity (m/s) REAL FUNCTION WZo_rep (wus, wuf) REAL, INTENT (IN) :: wus !! friction velocity (m/s) REAL, INTENT (IN) :: wuf !! _reference_ friction velocity (m/s) WZo_rep = (wus / wuf * (WZo_ref**0.06676))**(1.0/0.06676) END FUNCTION WZo_rep !! Calculate the aerodynamic roughness (mm) due to Allmaras random roughness of surface REAL FUNCTION WZo_rr (slrr) REAL, INTENT (IN) :: slrr !! Allmaras random roughness (mm) ! WZZO_MIN - Minimum wwzo value allowed for grid cell ! bare, non-ridged, surface aerodynamic roughness (mm) ! (currently this corresponds to ~1.67mm RR) ! WZZO_MAX - Maximum wwzo value allowed for grid cell ! bare, non-ridged, surface aerodynamic roughness (mm) ! (currently this corresponds to ~100mm RR) REAL WZZO_MIN, WZZO_MAX PARAMETER (WZZO_MIN = 0.5) !minimum wwzo value (mm) PARAMETER (WZZO_MAX = 30.0) !maximum wwzo value (mm) WZo_rr = 0.3 * slrr ! limit aerodynamic roughness value for RR range of: ~1.67mm <= RR <= ~100.0mm WZo_rr = min(WZZO_MAX, WZo_rr) WZo_rr = max(WZo_rr, WZZO_MIN) END FUNCTION WZo_rr !! Calculate the aerodynamic roughness (mm) due to oriented roughness (ridges) REAL FUNCTION WZo_rg (szrg, sxp_rg) REAL, INTENT (IN) :: szrg !! ridge height (mm) REAL, INTENT (IN) :: sxp_rg !! soil ridge spacing (mm) parallel to wind direction REAL h1 ! tmp variable IF (szrg .gt. 0.0) THEN ! no oriented roughness (ridges) h1 = szrg/sxp_rg ! ridge height/spacing ratio (normal to wind) ELSE WZo_rg = 0.0 RETURN END IF IF (szrg .lt. 5.0) THEN ! write (0,*) 'szrg < 5.0mm ', szrg, 'WZo_rg is ', & ! & szrg/(-64.1 + 135.5*h1 + (20.84/(h1)**0.5)) WZo_rg = 0.0 RETURN END IF IF (h1 > 0.2) THEN !winds are never continually normal to ridges ! write (0,*)'szrg/sxp_rg > 0.2',h1,'szrg is ',szrg,'WZo_rg is ',& ! & szrg/(-64.1 + 135.5*h1 + (20.84/(h1)**0.5)) h1 = 0.2 END IF IF (szrg .ge. 5.0) THEN WZo_rg = szrg/(-64.1 + 135.5*h1 + (20.84/(h1)**0.5)) END IF END FUNCTION WZo_rg !! Calculate the aerodynamic roughness (mm) due to total surface roughness (random and oriented) REAL FUNCTION WZo_rgrr (slrr, szrg, sxp_rg) REAL, INTENT (IN) :: slrr !! Allmaras random roughness (mm) REAL, INTENT (IN) :: szrg !! ridge height (mm) REAL, INTENT (IN) :: sxp_rg !! soil ridge spacing (mm) parallel to wind direction WZo_rgrr = max(WZo_rr(slrr), WZo_rg(szrg,sxp_rg)) END FUNCTION WZo_rgrr !! Calculate the aerodynamic roughness (mm) above canopy due to standing biomass REAL FUNCTION WZo_v (bht, wzo, brcd) REAL, INTENT (IN) :: bht !! standing biomass height (mm) REAL, INTENT (IN) :: wzo !! aerodynamic roughness (mm) of bare surface REAL, INTENT (IN) :: brcd !! effective biomass drag coefficient IF ( (brcd .eq. 0.0) .or. (bht .eq. 0.0) ) THEN ! no canopy WZo_v = wzo ELSE IF (brcd .ge. 0.1) THEN ! based upon average stem diameter of 20mm WZo_v = bht / (17.27 - (1.254*log(brcd)/brcd) - (3.714/brcd)) ELSE WZo_v = bht * ( (wzo/bht) + ((0.11-(wzo/bht))/4.60517)*log(brcd/0.001) ) END IF END FUNCTION WZo_v !! Calculate the _representative_ aerodynamic roughness (mm) at the surface below a standing biomass canopy REAL FUNCTION WZo_bv (slrr, szrg, sxp_rg, bht, brcd) REAL, INTENT (IN) :: slrr !! Allmaras random roughness (mm) REAL, INTENT (IN) :: szrg !! ridge height (mm) REAL, INTENT (IN) :: sxp_rg !! soil ridge spacing (mm) parallel to wind direction REAL, INTENT (IN) :: bht !! standing biomass height (mm) REAL, INTENT (IN) :: brcd !! effective biomass drag coefficient ! Doesn't matter what the actual reference friction velocity is, so we set it to 1.0 REAL, PARAMETER :: WUF = 1.0 !! _reference_ friction velocity (m/s) WZo_bv = WZo_rep(WUS_bv(WUF, bht, WZo_rgrr(slrr, szrg, sxp_rg), brcd), WUF) END FUNCTION WZo_bv !! Calculate the "effective" biomass drag coefficient due to standing biomass REAL FUNCTION BRcd_0 (br_lai, br_sai) REAL, INTENT (IN) :: br_lai !! biomass effective leaf area index (m^2/m^2) REAL, INTENT (IN) :: br_sai !! biomass effective stem area index (m^2/m^2) BRcd_0 = 0.2 * br_lai * (1.0 - exp(-br_lai)) + br_sai END FUNCTION BRcd_0 !! Calculate the "effective" biomass drag coefficient due to standing biomass REAL FUNCTION BRcd_1 (Rcrow, Rrg, Clai, bd_lai, Csai, bd_sai) REAL, INTENT (IN) :: Rcrow !! reduction of _effective_ leaf area when crop rows spaced > 5 crop heights REAL, INTENT (IN) :: Rrg !! reduction of _effective_ leaf and stem area when crop partly sheltered in furrow REAL, INTENT (IN) :: Clai !! live biomass (growing crop) leaf area index (m^2/m^2) REAL, INTENT (IN) :: bd_lai !! decomposing biomass (residue) leaf area index (m^2/m^2) REAL, INTENT (IN) :: Csai !! live biomass (growing crop) stem area index (m^2/m^2) REAL, INTENT (IN) :: bd_sai !! decomposing biomass (residue) stem area index (m^2/m^2) BRcd_1 = 0.2 * BRlai(Rcrow, Rrg, Clai, bd_lai) * (1.0 - exp(-BRlai(Rcrow, Rrg, Clai, bd_lai))) + BRsai(Rrg, Csai, bd_sai) END FUNCTION BRcd_1 !! Calculate the "effective" biomass drag coefficient due to standing biomass REAL FUNCTION BRcd_2 (surf_rough, crop, res) TYPE (WUS_roughness), INTENT (IN) :: surf_rough !! surface roughness parms (RR & OR) TYPE (WUS_crop), INTENT (IN) :: crop !! crop parms TYPE (WUS_residue), INTENT (IN) :: res !! residue parms BRcd_2 = 0.2 * BRlai(Rcrow(crop%ht, surf_rough%r_ht, crop%row_sp, crop%ridge_plant_flag), & Rrg(crop%ht, surf_rough%r_ht, crop%ridge_plant_flag), & crop%lai, res%lai) & * (1.0 - exp(-BRlai(Rcrow(crop%ht, surf_rough%r_ht, crop%row_sp, crop%ridge_plant_flag), & Rrg(crop%ht, surf_rough%r_ht, crop%ridge_plant_flag), & crop%lai, res%lai))) & + BRsai(Rrg(crop%ht, surf_rough%r_ht, crop%ridge_plant_flag), crop%sai, res%sai) END FUNCTION BRcd_2 !! Calculate the "effective" biomass drag coefficient due to standing biomass REAL FUNCTION BRcd_3 (crop_ht, ridge_ht, row_sp, ridge_plant_flag, Clai, bd_lai, Csai, bd_sai) REAL, INTENT (IN) :: crop_ht !! crop height (mm) REAL, INTENT (IN) :: ridge_ht !! ridge height (mm) REAL, INTENT (IN) :: row_sp !! crop row spacing (mm) INTEGER, INTENT (IN) :: ridge_plant_flag !! crop planted on ridges(1) or in furrows(0)? REAL, INTENT (IN) :: Clai !! live biomass (growing crop) leaf area index (m^2/m^2) REAL, INTENT (IN) :: bd_lai !! decomposing biomass (residue) leaf area index (m^2/m^2) REAL, INTENT (IN) :: Csai !! live biomass (growing crop) stem area index (m^2/m^2) REAL, INTENT (IN) :: bd_sai !! decomposing biomass (residue) stem area index (m^2/m^2) BRcd_3 = 0.2 * BRlai(Rcrow(crop_ht, ridge_ht, row_sp, ridge_plant_flag), & Rrg(crop_ht, ridge_ht, ridge_plant_flag), & Clai, bd_lai) & * (1.0 - exp(-BRlai(Rcrow(crop_ht, ridge_ht, row_sp, ridge_plant_flag), & Rrg(crop_ht, ridge_ht, ridge_plant_flag), & Clai, bd_lai))) & + BRsai(Rrg(crop_ht, ridge_ht, ridge_plant_flag), Csai, bd_sai) END FUNCTION BRcd_3 !! Calculate the _effective_ biomass stem area index (m^2/m^2) REAL FUNCTION BRsai (Rrg, Csai, bd_sai) REAL, INTENT (IN) :: Rrg !! reduction of _effective_ leaf and stem area when crop partly sheltered in furrow REAL, INTENT (IN) :: Csai !! live biomass (growing crop) stem area index (m^2/m^2) REAL, INTENT (IN) :: bd_sai !! decomposing biomass (residue) stem area index (m^2/m^2) BRsai = Rrg*Csai + bd_sai END FUNCTION BRsai !! Calculate the _effective_ biomass leaf area index (m^2/m^2) REAL FUNCTION BRlai (Rcrow, Rrg, Clai, bd_lai) REAL, INTENT (IN) :: Rcrow !! reduction of _effective_ leaf area when crop rows spaced > 5 crop heights REAL, INTENT (IN) :: Rrg !! reduction of _effective_ leaf and stem area when crop partly sheltered in furrow REAL, INTENT (IN) :: Clai !! live biomass (growing crop) leaf area index (m^2/m^2) REAL, INTENT (IN) :: bd_lai !! decomposing biomass (residue) leaf area index (m^2/m^2) BRlai = Rcrow*Rrg*Clai + bd_lai END FUNCTION BRlai !! Calculate reduction of _effective_ leaf area when crop rows spaced > 5 crop heights !! Reduction factor is set to zero if crop height is zero or "effective crop height" !! is less than the minimum crop height relative to ridge height when crop is planted !! in the furrow. REAL FUNCTION Rcrow (crop_ht, ridge_ht, row_sp, ridge_plant_flag) REAL, INTENT (IN) :: crop_ht !! crop height (mm) REAL, INTENT (IN) :: ridge_ht !! ridge height (mm) REAL, INTENT (IN) :: row_sp !! crop row spacing (mm) INTEGER, INTENT (IN) :: ridge_plant_flag !! crop planted on ridges(1) or in furrows(0)? REAL, PARAMETER :: c_ht_fract = 0.5 !! minimum crop height fraction relative to ridge height REAL, PARAMETER :: min_crop_ht_dist = 5.0 !! minimum (crop ht) distance threshold to trigger reduction factor REAL :: eff_crop_ht !! effective crop height (when crop planted in furrow) IF (ridge_plant_flag == 0) THEN ! crop planted in furrows eff_crop_ht = crop_ht - (c_ht_fract * ridge_ht) ELSE ! crop planted on ridges eff_crop_ht = crop_ht END IF IF (eff_crop_ht > 0.0) THEN IF (row_sp .gt. (min_crop_ht_dist * eff_crop_ht)) THEN Rcrow = 1.0 / (0.92 + (0.021*row_sp)/eff_crop_ht) ELSE Rcrow = 1.0 END IF ELSE Rcrow = 0.0 ! No effective crop height END IF END FUNCTION Rcrow !! Calculate reduction of _effective_ leaf and stem area when crop partly sheltered in furrow !! The intent is to reduce the lai and sai that can influence the wind flow since some of the !! crop mass is "sheltered" within the furrow. !! There are no adjustments made here with respect to ridge spacing relative to ridge height !! or wind direction relative to ridge direction. REAL FUNCTION Rrg (crop_ht, ridge_ht, ridge_plant_flag) REAL, INTENT (IN) :: crop_ht !! crop height (mm) REAL, INTENT (IN) :: ridge_ht !! ridge height (mm) INTEGER, INTENT (IN) :: ridge_plant_flag !! crop planted on ridges(1) or in furrows(0)? REAL, PARAMETER :: c_ht_fract = 0.5 !! minimum crop height fraction relative to ridge height IF (ridge_plant_flag == 0) THEN ! crop planted in furrows IF (crop_ht .gt. (c_ht_fract*ridge_ht)) THEN Rrg = 1.0 - (c_ht_fract*ridge_ht/crop_ht) ELSE Rrg = 0.0 ! crop height assumed too small (relative to ridge height) to influence wind END IF ELSE ! crop planted on ridges Rrg = 1.0 ! no reduction because crop is planted on ridges END IF END FUNCTION Rrg !! Calculate normalized angle relative to North (0,180] given direction angle !! Normalizes direction angle and returns equivalent angle in I & II quadrants REAL FUNCTION n180_angle (angle) REAL, INTENT (IN) :: angle !! direction angle relative to North (degrees) ! mod(angle, 360.0) --> put in unit circle [-360, 360] ! mod(mod(angle, 360.0)+360.0,360.0) --> make all positive values, e.g. (0,360] ! mod(mod(angle, 360.0)+360.0,180.0) --> put in first 2 quadrants, e.g. (0,180] n180_angle = mod(mod(angle, 360.0)+360.0,180.0) END FUNCTION n180_angle !! Calculate difference (acute angle) between two specified angles (intersecting lines) REAL FUNCTION rel_angle (angle1, angle2) REAL, INTENT (IN) :: angle1 !! direction angle relative to North (degrees) REAL, INTENT (IN) :: angle2 !! direction angle relative to North (degrees) REAL :: n1_angle, n2_angle n1_angle = n180_angle(angle1) n2_angle = n180_angle(angle2) IF (n1_angle > n2_angle) THEN rel_angle = n1_angle-n2_angle IF (rel_angle > 90.0) rel_angle = 180.0 - rel_angle ELSE rel_angle = n2_angle-n1_angle IF (rel_angle > 90.0) rel_angle = 180.0 - rel_angle END IF END FUNCTION rel_angle !! Calculate "_representative_" height when both standing crop and standing residue exist !! The value is determined based upon SAI weighting !! No consideration of planting location (ridge/furrow) or row spacing taken REAL FUNCTION rep_biomass_height (crop_ht, crop_sai, res_ht, res_sai) REAL, INTENT (IN) :: crop_ht !! crop height (mm) REAL, INTENT (IN) :: crop_sai !! live biomass (growing crop) stem area index (m^2/m^2) REAL, INTENT (IN) :: res_ht !! residue height (mm) REAL, INTENT (IN) :: res_sai !! decomposing biomass (residue) stem area index (m^2/m^2) REAL :: bio_sai = 0.0 REAL :: bio_ht = 0.0 bio_sai = crop_sai + res_sai IF (bio_sai .le. 0.0) THEN bio_ht = 0.0 ELSE bio_ht = crop_ht * crop_sai/bio_sai + res_ht * res_sai/bio_sai END IF rep_biomass_height = bio_ht END FUNCTION rep_biomass_height END MODULE Friction_Velocity