5.13.3-dev0
mHM
The mesoscale Hydrological Model
Loading...
Searching...
No Matches
mo_mpr_pet.f90
Go to the documentation of this file.
1!> \file mo_mpr_pet.f90
2!> \brief \copybrief mo_mpr_pet
3!> \details \copydetails mo_mpr_pet
4
5!> \brief MPR routine for PET.
6!> \details This module sets up pet correction factor at level-1 based on LAI
7!> \authors Mehmet Cuneyd Demirel, Simon Stisen
8!> \date May 2017
9!> \copyright Copyright 2005-\today, the mHM Developers, Luis Samaniego, Sabine Attinger: All rights reserved.
10!! mHM is released under the LGPLv3+ license \license_note
11!> \ingroup f_mpr
12module mo_mpr_pet
13
14 use mo_kind, only : i4, dp
15
16 implicit none
17
18 PUBLIC :: pet_correctbylai ! estimate PET correction factor with distributed LAI
19 PUBLIC :: pet_correctbyasp ! estimate PET correction factor with distributed Aspect
20 PUBLIC :: priestley_taylor_alpha ! factor (alpha) for Presley-Taylor ET estimation
21 !PUBLIC :: aerodynamical_resistance ! aerodynamical resistance (ra) for Penman-Monteith ET estimation
22 PUBLIC :: bulksurface_resistance ! bulk surface (stomatal) resistance (rs) for Penman-Monteith ET estimation
23
24 private
25
26contains
27
28
29 ! ----------------------------------------------------------------------------
30
31 ! NAME
32 ! pet_correctbyLAI
33
34 ! PURPOSE
35 !> \brief estimate PET correction factor based on LAI at L1
36
37 !> \details estimate PET correction factor based on LAI at L1 for a given
38 !> Leaf Area Index field.
39
40 !> Global parameters needed (see mhm_parameter.nml):
41 !> Process Case 5:
42 !> - param(1) = PET_a_forest
43 !> - param(2) = PET_a_impervious
44 !> - param(3) = PET_a_pervious
45 !> - param(4) = PET_b
46 !> - param(5) = PET_c
47
48 !> Example DSF=PET_a+PET_b*(1-exp(PET_c*LAI))
49
50 !> Similar to the crop coefficient concept Kc=a+b*(1-exp(c*LAI)) by Allen, R. G., L. S. Pereira,
51 !> D. Raes, and M. Smith (1998), Crop evapotranspiration - Guidelines for computing crop water requirements.,
52 !> FAO
53 !> Irrigation and drainage paper 56. See Chapter 9, Equation 97 <http://www.fao.org/docrep/X0490E/x0490e0f.htm>
54 !> Date: 17/5/2017
55
56 ! INTENT(IN)
57 !> \param[in] "real(dp), dimension(5) :: param" parameters
58 !> \param[in] "real(dp) :: nodata" - nodata value
59 !> \param[in] "integer(i4), dimension(:) :: LCOVER0" Land cover at level 0
60 !> \param[in] "real(dp), dimension(:, :) :: LAI0" LAI at level-0
61 !> \param[in] "logical, dimension(:, :) :: mask0" mask at L0
62 !> \param[in] "integer(i4), dimension(:) :: cell_id0" Cell ids of hi res field
63 !> \param[in] "integer(i4), dimension(:) :: upp_row_L1" Upper row of hi res block
64 !> \param[in] "integer(i4), dimension(:) :: low_row_L1" Lower row of hi res block
65 !> \param[in] "integer(i4), dimension(:) :: lef_col_L1" Left column of hi res block
66 !> \param[in] "integer(i4), dimension(:) :: rig_col_L1" Right column of hi res block
67 !> \param[in] "integer(i4), dimension(:) :: nL0_in_L1" Number of L0 cells within a L1 cel
68
69 ! INTENT(INOUT)
70 !> \param[inout] "real(dp), dimension(:, :) :: L1_petLAIcorFactor" pet cor factor at level-1
71
72 ! HISTORY
73 !> \authors M. Cuneyd Demirel and Simon Stisen from GEUS.dk
74
75 !> \date May. 2017
76
77 ! Modifications:
78 ! Robert Schweppe Jun 2018 - refactoring and reformatting
79
80 subroutine pet_correctbylai(param, nodata, LCOVER0, LAI0, mask0, cell_id0, upp_row_L1, low_row_L1, lef_col_L1, &
81 rig_col_L1, nL0_in_L1, L1_petLAIcorFactor)
82
83 use mo_upscaling_operators, only : upscale_harmonic_mean
84 !$ use omp_lib
85
86 implicit none
87
88 ! parameters
89 real(dp), dimension(5), intent(in) :: param
90
91 ! - nodata value
92 real(dp), intent(in) :: nodata
93
94 ! Land cover at level 0
95 integer(i4), dimension(:), intent(in) :: lcover0
96
97 ! LAI at level-0
98 real(dp), dimension(:, :), intent(in) :: lai0
99
100 ! mask at L0
101 logical, dimension(:, :), intent(in) :: mask0
102
103 ! Cell ids of hi res field
104 integer(i4), dimension(:), intent(in) :: cell_id0
105
106 ! Upper row of hi res block
107 integer(i4), dimension(:), intent(in) :: upp_row_l1
108
109 ! Lower row of hi res block
110 integer(i4), dimension(:), intent(in) :: low_row_l1
111
112 ! Left column of hi res block
113 integer(i4), dimension(:), intent(in) :: lef_col_l1
114
115 ! Right column of hi res block
116 integer(i4), dimension(:), intent(in) :: rig_col_l1
117
118 ! Number of L0 cells within a L1 cel
119 integer(i4), dimension(:), intent(in) :: nl0_in_l1
120
121 ! pet cor factor at level-1
122 real(dp), dimension(:, :), intent(inout) :: l1_petlaicorfactor
123
124 ! pet cor factor at level-0
125 real(dp), dimension(size(LCOVER0, 1), size(LAI0, 2)) :: petlaicorfactor_0
126
127 ! loop index
128 integer(i4) :: kk, tt
129
130 ! loop index
131 integer(i4) :: ll
132
133
134 ! ------------------------------------------------------------------
135 ! Estimate DSF=PET_a+PET_b*(1-exp(PET_c*LAI)) to correct PET as PET=DSF*PET
136 ! ------------------------------------------------------------------
137 !$OMP PARALLEL
138 !$OMP DO PRIVATE( LL ) SCHEDULE( STATIC )
139 ! need to be done for every landcover to get DSF
140 do kk = 1, size(lcover0, 1)
141
142 ll = lcover0(kk)
143
144 ! TODO: memory order of arrays is not optimal, how to improve?
145 select case(ll)
146 case(1) ! forest
147 petlaicorfactor_0(kk, :) = param(1) + (param(4) * (1.0_dp - exp(param(5) * lai0(kk, :))))
148 case(2) ! impervious
149 petlaicorfactor_0(kk, :) = param(2) + (param(4) * (1.0_dp - exp(param(5) * lai0(kk, :))))
150 case(3) ! permeable
151 petlaicorfactor_0(kk, :) = param(3) + (param(4) * (1.0_dp - exp(param(5) * lai0(kk, :))))
152 end select
153
154 end do
155 !$OMP END DO
156 !$OMP END PARALLEL
157
158 do tt = 1, size(lai0, 2)
159 l1_petlaicorfactor(:, tt) = upscale_harmonic_mean(nl0_in_l1, upp_row_l1, low_row_l1, &
160 lef_col_l1, rig_col_l1, cell_id0, mask0, nodata, petlaicorfactor_0(:, tt))
161 end do
162
163
164 end subroutine pet_correctbylai
165
166 ! ----------------------------------------------------------------------------
167
168 ! NAME
169 ! pet_correctbyASP
170
171 ! PURPOSE
172 !> \brief correction of PET
173
174 !> \details Correction of PET based on L0 aspect data.
175
176
177 !> Global parameters needed (see mhm_parameter.nml):
178 !> - param(1) = minCorrectionFactorPET
179 !> - param(2) = maxCorrectionFactorPET
180 !> - param(3) = aspectTresholdPET
181
182 ! INTENT(IN)
183 !> \param[in] "integer(i4), dimension(:) :: id0" Level 0 cell id
184 !> \param[in] "real(dp), dimension(:) :: latitude_l0" latitude on l0
185 !> \param[in] "real(dp), dimension(:) :: Asp0" [degree] Aspect at Level 0
186 !> \param[in] "real(dp), dimension(3) :: param" process parameters
187 !> \param[in] "real(dp) :: nodata" - no data value
188
189 ! INTENT(OUT)
190 !> \param[out] "real(dp), dimension(:) :: fAsp0" PET correction for Aspect
191
192 ! HISTORY
193 !> \authors Stephan Thober, Rohini Kumar
194
195 !> \date Dec 2012
196
197 ! Modifications:
198 ! Juliane Mai Oct 2013 - OLD parametrization --> param(1) = minCorrectionFactorPET
199 ! --> param(2) = maxCorrectionFactorPET
200 ! --> param(3) = aspectTresholdPET
201 ! -------------------------------
202 ! maxCorrectionFactorPET = minCorrectionFactorPET + delta
203 ! -------------------------------
204 ! NEW parametrization
205 ! --> param(1) = minCorrectionFactorPET
206 ! --> param(2) = delta
207 ! --> param(3) = aspectTresholdPET
208 ! Stephan Thober Dec 2013 - changed intent(inout) to intent(out)
209 ! Stephan Thober Sep 2015 - Mapping L1 to Lo, latitude on L0
210 ! Luis Samaniego Sep 2015 - PET correction on the southern hemisphere
211 ! Matthias Zink Jun 2017 - renamed and moved from mo_multi_scale_param_reg.f90 to mo_mpr_pet.f90
212 ! Robert Schweppe Jun 2018 - refactoring and reformatting
213
214 subroutine pet_correctbyasp(Id0, latitude_l0, Asp0, param, nodata, fAsp0)
215 !$ use omp_lib
216
217 implicit none
218
219 ! Level 0 cell id
220 integer(i4), dimension(:), intent(in) :: id0
221
222 ! latitude on l0
223 real(dp), dimension(:), intent(in) :: latitude_l0
224
225 ! - no data value
226 real(dp), intent(in) :: nodata
227
228 ! process parameters
229 real(dp), dimension(3), intent(in) :: param
230
231 ! [degree] Aspect at Level 0
232 real(dp), dimension(:), intent(in) :: asp0
233
234 ! PET correction for Aspect
235 real(dp), dimension(:), intent(out) :: fasp0
236
237 ! PET correction for Aspect, south
238 real(dp), dimension(size(id0, 1)) :: fasp0s
239
240 logical, dimension(size(id0, 1)) :: mask_north_hemisphere_l0
241
242 real(dp) :: tmp_maxcorrectionfactorpet
243
244
245 mask_north_hemisphere_l0 = merge(.true., .false., latitude_l0 .gt. 0.0_dp)
246
247 tmp_maxcorrectionfactorpet = param(1) + param(2)
248
249 ! for cells on the northern hemisphere
250 !$OMP PARALLEL
251 fasp0 = merge(&
252 param(1) + (tmp_maxcorrectionfactorpet - param(1)) / param(3) * asp0, &
253 param(1) + (tmp_maxcorrectionfactorpet - param(1)) / (360._dp - param(3)) * (360._dp - asp0), &
254 ! ( asp0 < param(3) ) .and. mask_north_hemisphere_l0 )
255 asp0 < param(3))
256 fasp0 = merge(fasp0, nodata, id0 /= int(nodata, i4))
257 !$OMP END PARALLEL
258
259 ! for cells on the southern hemisphere
260 !$OMP PARALLEL
261 fasp0s = merge(&
262 param(1) + (tmp_maxcorrectionfactorpet - param(1)) / (360._dp - param(3)) * (360._dp - asp0), &
263 param(1) + (tmp_maxcorrectionfactorpet - param(1)) / param(3) * asp0, &
264 asp0 < param(3))
265 fasp0s = merge(fasp0s, nodata, id0 /= int(nodata, i4))
266 !$OMP END PARALLEL
267
268 !$OMP PARALLEL
269 fasp0 = merge(fasp0, fasp0s, mask_north_hemisphere_l0)
270 !$OMP END PARALLEL
271
272 end subroutine pet_correctbyasp
273
274 ! ----------------------------------------------------------------------------
275
276 ! NAME
277 ! priestley_taylor_alpha
278
279 ! PURPOSE
280 !> \brief Regionalization of priestley taylor alpha
281
282 !> \details estimation of priestley taylor alpha
283 !> Global parameters needed (see mhm_parameter.nml):
284 !> - param(1) = PriestleyTaylorCoeff
285 !> - param(2) = PriestleyTaylorLAIcorr
286
287 ! INTENT(IN)
288 !> \param[in] "real(dp), dimension(:, :) :: LAI0" LAI at level-0
289 !> \param[in] "real(dp), dimension(:) :: param" input parameter
290 !> \param[in] "logical, dimension(:, :) :: mask0" mask at level 0
291 !> \param[in] "real(dp) :: nodata" - nodata value
292 !> \param[in] "integer(i4), dimension(:) :: cell_id0" Cell ids of hi res field
293 !> \param[in] "integer(i4), dimension(:) :: nL0_in_L1" number of l0 cells within a l1 cell
294 !> \param[in] "integer(i4), dimension(:) :: Upp_row_L1" upper row of a l1 cell in l0 grid
295 !> \param[in] "integer(i4), dimension(:) :: Low_row_L1" lower row of a l1 cell in l0 grid
296 !> \param[in] "integer(i4), dimension(:) :: Lef_col_L1" left col of a l1 cell in l0 grid
297 !> \param[in] "integer(i4), dimension(:) :: Rig_col_L1" right col of a l1 cell in l0 grid
298
299 ! INTENT(OUT)
300 !> \param[out] "real(dp), dimension(:, :) :: priestley_taylor_alpha1" bulk surface resistance
301
302 ! HISTORY
303 !> \authors Matthias Zink
304
305 !> \date Apr 2013
306
307 ! Modifications:
308 ! Matthias Zink Jun 2017 - moved from mo_multi_scale_param_reg.f90 to mo_mpr_pet.f90
309 ! Robert Schweppe Jun 2018 - refactoring and reformatting
310
311 subroutine priestley_taylor_alpha(LAI0, param, mask0, nodata, cell_id0, nL0_in_L1, Upp_row_L1, Low_row_L1, Lef_col_L1, &
312 Rig_col_L1, priestley_taylor_alpha1)
313
314 use mo_upscaling_operators, only : upscale_arithmetic_mean
315
316 implicit none
317
318 ! LAI at level-0
319 real(dp), dimension(:, :), intent(in) :: lai0
320
321 ! input parameter
322 real(dp), dimension(:), intent(in) :: param
323
324 ! mask at level 0
325 logical, dimension(:, :), intent(in) :: mask0
326
327 ! - nodata value
328 real(dp), intent(in) :: nodata
329
330 ! Cell ids of hi res field
331 integer(i4), dimension(:), intent(in) :: cell_id0
332
333 ! number of l0 cells within a l1 cell
334 integer(i4), dimension(:), intent(in) :: nl0_in_l1
335
336 ! upper row of a l1 cell in l0 grid
337 integer(i4), dimension(:), intent(in) :: upp_row_l1
338
339 ! lower row of a l1 cell in l0 grid
340 integer(i4), dimension(:), intent(in) :: low_row_l1
341
342 ! left col of a l1 cell in l0 grid
343 integer(i4), dimension(:), intent(in) :: lef_col_l1
344
345 ! right col of a l1 cell in l0 grid
346 integer(i4), dimension(:), intent(in) :: rig_col_l1
347
348 ! bulk surface resistance
349 real(dp), dimension(:, :), intent(out) :: priestley_taylor_alpha1
350
351 integer(i4) :: tt
352
353 ! dim 1 = number of cells on level 0, time
354 real(dp), dimension(:, :), allocatable :: priestley_taylor_alpha0
355
356
357 ! initialize some things
358 allocate(priestley_taylor_alpha0(size(lai0, 1), size(lai0, 2))) ; priestley_taylor_alpha0 = nodata
359 priestley_taylor_alpha1 = nodata
360 !
361 do tt = 1, size(lai0, 2)
362 priestley_taylor_alpha0(:, tt) = param(1) + param(2) * lai0(:, tt)
363
364 priestley_taylor_alpha1(:, tt) = upscale_arithmetic_mean(nl0_in_l1, upp_row_l1, low_row_l1, &
365 lef_col_l1, rig_col_l1, cell_id0, mask0, nodata, priestley_taylor_alpha0(:, tt))
366
367 end do
368
369 end subroutine priestley_taylor_alpha
370
371 ! ----------------------------------------------------------------------------
372
373 ! NAME
374 ! bulksurface_resistance
375
376 ! PURPOSE
377 !> \brief Regionalization of bulk surface resistance
378
379 !> \details estimation of bulk surface resistance
380 !> Global parameters needed (see mhm_parameter.nml):
381 !> - param(1) = stomatal_resistance
382
383 ! INTENT(IN)
384 !> \param[in] "real(dp), dimension(:, :) :: LAI0" LAI at level-0
385 !> \param[in] "real(dp) :: param" - global parameter
386 !> \param[in] "logical, dimension(:, :) :: mask0" mask at level 0
387 !> \param[in] "real(dp) :: nodata" - nodata value
388 !> \param[in] "integer(i4), dimension(:) :: cell_id0" Cell ids of hi res field
389 !> \param[in] "integer(i4), dimension(:) :: nL0_in_L1" number of l0 cells within a l1 cell
390 !> \param[in] "integer(i4), dimension(:) :: Upp_row_L1" upper row of a l1 cell in l0 grid
391 !> \param[in] "integer(i4), dimension(:) :: Low_row_L1" lower row of a l1 cell in l0 grid
392 !> \param[in] "integer(i4), dimension(:) :: Lef_col_L1" left col of a l1 cell in l0 grid
393 !> \param[in] "integer(i4), dimension(:) :: Rig_col_L1" right col of a l1 cell in l0 grid
394
395 ! INTENT(OUT)
396 !> \param[out] "real(dp), dimension(:, :) :: bulksurface_resistance1" bulk surface resistance
397
398 ! HISTORY
399 !> \authors Matthias Zink
400
401 !> \date Apr 2013
402
403 ! Modifications:
404 ! Matthias Zink Jun 2017 - moved from mo_multi_scale_param_reg.f90 to mo_mpr_pet.f90
405 ! Robert Schweppe Jun 2018 - refactoring and reformatting
406
407 subroutine bulksurface_resistance(LAI0, param, mask0, nodata, cell_id0, nL0_in_L1, Upp_row_L1, Low_row_L1, Lef_col_L1, &
408 Rig_col_L1, bulksurface_resistance1)
409
410 use mo_mpr_constants, only : lai_factor_surfresi, lai_offset_surfresi, max_surfresist
411 use mo_upscaling_operators, only : upscale_arithmetic_mean
412
413 implicit none
414
415 ! LAI at level-0
416 real(dp), dimension(:, :), intent(in) :: lai0
417
418 ! - global parameter
419 real(dp), intent(in) :: param
420
421 ! mask at level 0
422 logical, dimension(:, :), intent(in) :: mask0
423
424 ! - nodata value
425 real(dp), intent(in) :: nodata
426
427 ! Cell ids of hi res field
428 integer(i4), dimension(:), intent(in) :: cell_id0
429
430 ! number of l0 cells within a l1 cell
431 integer(i4), dimension(:), intent(in) :: nl0_in_l1
432
433 ! upper row of a l1 cell in l0 grid
434 integer(i4), dimension(:), intent(in) :: upp_row_l1
435
436 ! lower row of a l1 cell in l0 grid
437 integer(i4), dimension(:), intent(in) :: low_row_l1
438
439 ! left col of a l1 cell in l0 grid
440 integer(i4), dimension(:), intent(in) :: lef_col_l1
441
442 ! right col of a l1 cell in l0 grid
443 integer(i4), dimension(:), intent(in) :: rig_col_l1
444
445 ! bulk surface resistance
446 real(dp), dimension(:, :), intent(out) :: bulksurface_resistance1
447
448 integer(i4) :: tt
449
450 ! dim 1 = number of cells on level 0,
451 ! dim 2 = number of months in year (12)
452 real(dp), dimension(:, :), allocatable :: bulksurface_resistance0
453
454
455 ! initialize some things
456 allocate(bulksurface_resistance0(size(lai0, 1), size(lai0, 2))) ; bulksurface_resistance0 = nodata
457 bulksurface_resistance1 = nodata
458 !
459
460 do tt = 1, size(lai0, 2)
461 bulksurface_resistance0(:, tt) = param / (lai0(:, tt) / &
462 (lai_factor_surfresi * lai0(:, tt) + lai_offset_surfresi))
463 ! efeective LAI from McMahon et al ,2013 , HESS supplements
464
465 ! since LAI may be very low, rs becomes very high
466 ! thus the values are restricted to maximum literaure values (i.e. McMahon et al ,2013 , HESS)
467 bulksurface_resistance0(:, tt) = merge(max_surfresist, bulksurface_resistance0(:, tt), &
468 bulksurface_resistance0(:, tt) .GT. max_surfresist)
469
470 bulksurface_resistance1(:, tt) = upscale_arithmetic_mean(nl0_in_l1, upp_row_l1, low_row_l1, &
471 lef_col_l1, rig_col_l1, cell_id0, mask0, nodata, bulksurface_resistance0(:, tt))
472
473 end do
474
475 end subroutine bulksurface_resistance
476
477end module mo_mpr_pet
mo_mpr_constants
Provides MPR specific constants.
Definition mo_mpr_constants.f90:12
mo_mpr_constants::lai_offset_surfresi
real(dp), parameter, public lai_offset_surfresi
LAI offset for bulk surface resistance formulation.
Definition mo_mpr_constants.f90:72
mo_mpr_constants::lai_factor_surfresi
real(dp), parameter, public lai_factor_surfresi
LAI factor for bulk surface resistance formulation.
Definition mo_mpr_constants.f90:70
mo_mpr_constants::max_surfresist
real(dp), parameter, public max_surfresist
maximum bulk surface resistance
Definition mo_mpr_constants.f90:74
mo_mpr_pet
MPR routine for PET.
Definition mo_mpr_pet.f90:12
mo_mpr_pet::pet_correctbyasp
subroutine, public pet_correctbyasp(id0, latitude_l0, asp0, param, nodata, fasp0)
correction of PET
Definition mo_mpr_pet.f90:215
mo_mpr_pet::bulksurface_resistance
subroutine, public bulksurface_resistance(lai0, param, mask0, nodata, cell_id0, nl0_in_l1, upp_row_l1, low_row_l1, lef_col_l1, rig_col_l1, bulksurface_resistance1)
Regionalization of bulk surface resistance.
Definition mo_mpr_pet.f90:409
mo_mpr_pet::priestley_taylor_alpha
subroutine, public priestley_taylor_alpha(lai0, param, mask0, nodata, cell_id0, nl0_in_l1, upp_row_l1, low_row_l1, lef_col_l1, rig_col_l1, priestley_taylor_alpha1)
Regionalization of priestley taylor alpha.
Definition mo_mpr_pet.f90:313
mo_mpr_pet::pet_correctbylai
subroutine, public pet_correctbylai(param, nodata, lcover0, lai0, mask0, cell_id0, upp_row_l1, low_row_l1, lef_col_l1, rig_col_l1, nl0_in_l1, l1_petlaicorfactor)
estimate PET correction factor based on LAI at L1
Definition mo_mpr_pet.f90:82
mo_upscaling_operators
Module containing upscaling operators.
Definition mo_upscaling_operators.f90:12
mo_upscaling_operators::upscale_arithmetic_mean
real(dp) function, dimension(size(nl0_cells_in_l1_cell, 1)), public upscale_arithmetic_mean(nl0_cells_in_l1_cell, l1_upper_rowid_cell, l1_lower_rowid_cell, l1_left_colonid_cell, l1_right_colonid_cell, l0_cellid, mask0, nodata_value, l0_finescale_data)
aritmetic mean
Definition mo_upscaling_operators.f90:268
mo_upscaling_operators::upscale_harmonic_mean
real(dp) function, dimension(size(nl0_cells_in_l1_cell, 1)), public upscale_harmonic_mean(nl0_cells_in_l1_cell, l1_upper_rowid_cell, l1_lower_rowid_cell, l1_left_colonid_cell, l1_right_colonid_cell, l0_cellid, mask0, nodata_value, l0_finescale_data)
harmonic mean
Definition mo_upscaling_operators.f90:371