dicussion needed on calculation of aerodynamic resistance
If PET is calculated using Penman-Monteith, aerodynamic resistance needs to be calculated. This is happening in subroutine aerodynamic_resistance
in mo_multi_param_reg.f90
. It was implemented by Matthias Zink and cites a FAO publication as a source. Now here is the problem:
The FAO paper uses a calculation based on a paper of Jacobs(2002) that seems to be valid for the grass reference surface. What I found in a brief search in the internet is that it is most commonly used for canopy_heights for crops of up to 2m. However, we also use that function for forest canopies. The canopy height for forests in mHM is controlled by a parameter (param(1)
) and can range from 15 (default) to 40m. The measurement height of wind measurement (zm
) is a parameter and set to 10m.
To make that function work, the canopy_height0
is added to zm
, if the canopy_height0
is higher than zm
. I do not know, if that is really the way the function is supposed to be used.
Here is the relevant code:
! regionalization of canopy height
! substitute with canopy height
canopy_height0 = merge(param(1), canopy_height0, LCover0 == 1) ! forest
canopy_height0 = merge(param(2), canopy_height0, LCover0 == 2) ! impervious
! old implementation used values from LUT statically for all cells (Jan-Dec, 12 values):
! 7 Intensive-orchards: 2.0, 2.0, 2.0, 2.0, 3.0, 3.5, 4.0, 4.0, 4.0, 2.5, 2.0, 2.0
maxLAI = MAXVAL(LAI0, dim=2)
do tt = 1, size(LAI0, 2)
! pervious canopy height is scaled with LAI
canopy_height0 = merge((param(3) * LAI0(:, tt) / maxLAI), canopy_height0, LCover0 == 3) ! pervious
! estimation of the aerodynamic resistance on the lower level
! see FAO Irrigation and Drainage Paper No. 56 (p. 19 ff) for more information
zm = WindMeasHeight
! correction: if wind measurement height is below canopy height loagarithm becomes negative
->!! zm = merge(canopy_height0 + zm, zm, ((abs(zm - nodata_dp) .GT. eps_dp) .AND. (zm .LT. canopy_height0)))
! zh = zm
displace = param(4) * canopy_height0
zm_zero = param(5) * canopy_height0
zh_zero = param(6) * zm_zero
!
! calculate aerodynamic resistance (changes monthly)
aerodyn_resistance0(:, tt) = log((zm - displace) / zm_zero) * log((zm - displace) / zh_zero) / (karman**2.0_dp)
aerodyn_resistance1(:, tt) = upscale_arithmetic_mean(n_subcells1, upper_bound1, lower_bound1, &
left_bound1, right_bound1, Id0, mask0, nodata_dp, aerodyn_resistance0(:, tt))
end do
So my questions are (before I adapt this to the new MPR) @lese, @rakovec, @rkumar, @thober:
- are we okay with using this function for calculation of aerodynamic resistance for forests?
- if yes, is the correction, like Matthias did it, okay, or do we come up with a smoother solution, because now: canopy_height: 9.99 m -> zm = 10.0 m, canopy_height: 10.0 m -> zm = 20.0 m
- if no, what other transfer function do we want?