!> @file biometeorology_mod.f90
!--------------------------------------------------------------------------------!
! This file is part of PALM-4U.
!
! PALM-4U is free software: you can redistribute it and/or modify it under the
! terms of the GNU General Public License as published by the Free Software
! Foundation, either version 3 of the License, or (at your option) any later
! version.
!
! PALM-4U is distributed in the hope that it will be useful, but WITHOUT ANY
! WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
! A PARTICULAR PURPOSE. See the GNU General Public License for more details.
!
! You should have received a copy of the GNU General Public License along with
! PALM. If not, see .
!
! Copyright 2018-2018 Deutscher Wetterdienst (DWD)
! Copyright 2018-2018 Institute of Computer Science, Academy of Sciences, Prague
! Copyright 2018-2018 Leibniz Universitaet Hannover
!--------------------------------------------------------------------------------!
!
! Current revisions:
! -----------------
!
!
! Former revisions:
! -----------------
! $Id: biometeorology_mod.f90 3479 2018-11-01 16:00:30Z suehring $
! - reworked check for output quantities
! - assign new palm-error numbers
! - set unit according to data standard.
!
! 3475 2018-10-30 21:16:31Z kanani
! Add option for using MRT from RTM instead of simplified global value
!
! 3464 2018-10-30 18:08:55Z kanani
! From branch resler@3462, pavelkrc:
! make use of basic_constants_and_equations_mod
!
! 3448 2018-10-29 18:14:31Z kanani
! Initial revision
!
!
!
! Authors:
! --------
! @author Dominik Froehlich
! @author Jaroslav Resler
!
!
! Description:
! ------------
!> Human thermal comfort module calculating thermal perception of a sample
!> human being under the current meteorological conditions.
!>
!> @todo Alphabetical sorting of "USE ..." lists, "ONLY" list, variable declarations
!> (per subroutine: first all CHARACTERs, then INTEGERs, LOGICALs, REALs, )
!> @todo Comments start with capital letter --> "!-- Include..."
!> @todo Variable and routine names strictly in lowercase letters and in English
!>
!> @note nothing now
!>
!> @bug no known bugs by now
!------------------------------------------------------------------------------!
MODULE biometeorology_mod
USE arrays_3d, &
ONLY: pt, p, u, v, w, q
USE averaging, &
ONLY: pt_av, q_av, u_av, v_av, w_av
USE basic_constants_and_equations_mod, &
ONLY: degc_to_k, magnus, sigma_sb
USE biometeorology_ipt_mod
USE biometeorology_pet_mod
USE biometeorology_pt_mod, &
ONLY: calculate_pt_static
USE biometeorology_utci_mod
USE control_parameters, &
ONLY: average_count_3d, biometeorology, dz, dz_stretch_factor, &
dz_stretch_level, humidity, initializing_actions, nz_do3d, &
simulated_time, surface_pressure
USE grid_variables, &
ONLY: ddx, dx, ddy, dy
USE indices, &
ONLY: nxl, nxr, nys, nyn, nzb, nzt, nys, nyn, nxl, nxr
USE kinds !< Set precision of INTEGER and REAL arrays according to PALM
!-- Import radiation model to obtain input for mean radiant temperature
USE radiation_model_mod, &
ONLY: ix, iy, iz, id, mrt_nlevels, mrt_include_sw, &
mrtinsw, mrtinlw, mrtbl, nmrtbl, radiation, &
radiation_interactions, rad_sw_in, &
rad_sw_out, rad_lw_in, rad_lw_out
USE surface_mod, &
ONLY: get_topography_top_index_ji
IMPLICIT NONE
PRIVATE
!-- Declare all global variables within the module (alphabetical order)
REAL(wp), DIMENSION(:,:), ALLOCATABLE :: tmrt_grid !< tmrt results (°C)
REAL(wp), DIMENSION(:,:), ALLOCATABLE :: pt_grid !< PT results (°C)
REAL(wp), DIMENSION(:,:), ALLOCATABLE :: utci_grid !< UTCI results (°C)
REAL(wp), DIMENSION(:,:), ALLOCATABLE :: pet_grid !< PET results (°C)
!-- Grids for averaged thermal indices
REAL(wp), DIMENSION(:,:), ALLOCATABLE :: pt_av_grid !< PT results (aver. input) (°C)
REAL(wp), DIMENSION(:,:), ALLOCATABLE :: utci_av_grid !< UTCI results (aver. input) (°C)
REAL(wp), DIMENSION(:,:), ALLOCATABLE :: pet_av_grid !< PET results (aver. input) (°C)
!-- Grids for radiation_model
REAL(wp), DIMENSION(:), ALLOCATABLE :: biom_mrt !< biometeorology mean radiant temperature for each MRT box
REAL(wp), DIMENSION(:), ALLOCATABLE :: biom_mrt_av !< time average mean
INTEGER( iwp ) :: biom_cell_level !< cell level biom calculates for
REAL ( wp ) :: biom_output_height !< height output is calculated in m
REAL ( wp ) :: time_biom_results !< the time, the last set of biom results have been calculated for
REAL ( wp ), PARAMETER :: human_absorb = 0.7_wp !< SW absorbtivity of a human body (Fanger 1972)
REAL ( wp ), PARAMETER :: human_emiss = 0.97_wp !< LW emissivity of a human body after (Fanger 1972)
!--
LOGICAL :: aver_pt = .FALSE. !< switch: do pt averaging in this module? (if .FALSE. this is done globally)
LOGICAL :: aver_q = .FALSE. !< switch: do e averaging in this module?
LOGICAL :: aver_u = .FALSE. !< switch: do u averaging in this module?
LOGICAL :: aver_v = .FALSE. !< switch: do v averaging in this module?
LOGICAL :: aver_w = .FALSE. !< switch: do w averaging in this module?
LOGICAL :: average_trigger_pt = .FALSE. !< update averaged input on call to biom_pt?
LOGICAL :: average_trigger_utci = .FALSE. !< update averaged input on call to biom_utci?
LOGICAL :: average_trigger_pet = .FALSE. !< update averaged input on call to biom_pet?
LOGICAL :: mrt_from_rtm = .TRUE. !< switch: use mrt calculated by RTM for calculation of thermal indices
LOGICAL :: biom_pt = .TRUE. !< Turn index PT (instant. input) on or off
LOGICAL :: biom_pt_av = .TRUE. !< Turn index PT (averaged input) on or off
LOGICAL :: biom_pet = .TRUE. !< Turn index PET (instant. input) on or off
LOGICAL :: biom_pet_av = .TRUE. !< Turn index PET (averaged input) on or off
LOGICAL :: biom_utci = .TRUE. !< Turn index UTCI (instant. input) on or off
LOGICAL :: biom_utci_av = .TRUE. !< Turn index UTCI (averaged input) on or off
!-- Add INTERFACES that must be available to other modules (alphabetical order)
PUBLIC biom_3d_data_averaging, biom_check_data_output, &
biom_calculate_static_grid, biom_calc_ipt, &
biom_check_parameters, biom_data_output_3d, biom_data_output_2d, &
biom_define_netcdf_grid, biom_determine_input_at, biom_header, biom_init, &
biom_init_arrays, biom_parin, biom_pt, biom_pt_av, biom_pet, biom_pet_av, &
biom_utci, biom_utci_av, time_biom_results
!
!-- PALM interfaces:
!
!-- 3D averaging for HTCM _INPUT_ variables
INTERFACE biom_3d_data_averaging
MODULE PROCEDURE biom_3d_data_averaging
END INTERFACE biom_3d_data_averaging
!-- Calculate static thermal indices PT, UTCI and/or PET
INTERFACE biom_calculate_static_grid
MODULE PROCEDURE biom_calculate_static_grid
END INTERFACE biom_calculate_static_grid
!-- Calculate the dynamic index iPT (to be caled by the agent model)
INTERFACE biom_calc_ipt
MODULE PROCEDURE biom_calc_ipt
END INTERFACE biom_calc_ipt
!-- Data output checks for 2D/3D data to be done in check_parameters
INTERFACE biom_check_data_output
MODULE PROCEDURE biom_check_data_output
END INTERFACE biom_check_data_output
!-- Input parameter checks to be done in check_parameters
INTERFACE biom_check_parameters
MODULE PROCEDURE biom_check_parameters
END INTERFACE biom_check_parameters
!-- Data output of 2D quantities
INTERFACE biom_data_output_2d
MODULE PROCEDURE biom_data_output_2d
END INTERFACE biom_data_output_2d
!-- no 3D data, thus, no averaging of 3D data, removed
INTERFACE biom_data_output_3d
MODULE PROCEDURE biom_data_output_3d
END INTERFACE biom_data_output_3d
!-- Definition of data output quantities
INTERFACE biom_define_netcdf_grid
MODULE PROCEDURE biom_define_netcdf_grid
END INTERFACE biom_define_netcdf_grid
!-- Obtains all relevant input values to estimate local thermal comfort/stress
INTERFACE biom_determine_input_at
MODULE PROCEDURE biom_determine_input_at
END INTERFACE biom_determine_input_at
!-- Output of information to the header file
INTERFACE biom_header
MODULE PROCEDURE biom_header
END INTERFACE biom_header
!-- Initialization actions
INTERFACE biom_init
MODULE PROCEDURE biom_init
END INTERFACE biom_init
!-- Initialization of arrays
INTERFACE biom_init_arrays
MODULE PROCEDURE biom_init_arrays
END INTERFACE biom_init_arrays
!-- Reading of NAMELIST parameters
INTERFACE biom_parin
MODULE PROCEDURE biom_parin
END INTERFACE biom_parin
CONTAINS
!------------------------------------------------------------------------------!
! Description:
! ------------
!> Sum up and time-average biom input quantities as well as allocate
!> the array necessary for storing the average.
!------------------------------------------------------------------------------!
SUBROUTINE biom_3d_data_averaging( mode, variable )
IMPLICIT NONE
CHARACTER (LEN=*) :: mode !<
CHARACTER (LEN=*) :: variable !<
INTEGER(iwp) :: i !<
INTEGER(iwp) :: j !<
INTEGER(iwp) :: k !<
IF ( mode == 'allocate' ) THEN
SELECT CASE ( TRIM( variable ) )
CASE ( 'biom_mrt' )
IF ( .NOT. ALLOCATED( biom_mrt_av ) ) THEN
ALLOCATE( biom_mrt_av(nmrtbl) )
ENDIF
biom_mrt_av = 0.0_wp
CASE ( 'biom_pt', 'biom_utci', 'biom_pet' )
!-- Indices in unknown order as depending on input file, determine
! first index to average und update only once
IF ( .NOT. average_trigger_pt .AND. .NOT. average_trigger_utci &
.AND. .NOT. average_trigger_pet ) THEN
IF ( TRIM( variable ) == 'biom_pt' ) THEN
average_trigger_pt = .TRUE.
ENDIF
IF ( TRIM( variable ) == 'biom_utci' ) THEN
average_trigger_utci = .TRUE.
ENDIF
IF ( TRIM( variable ) == 'biom_pet' ) THEN
average_trigger_pet = .TRUE.
ENDIF
ENDIF
!-- Only continue if updateindex
IF ( average_trigger_pt .AND. TRIM( variable ) /= 'biom_pt') RETURN
IF ( average_trigger_utci .AND. TRIM( variable ) /= 'biom_utci') RETURN
IF ( average_trigger_pet .AND. TRIM( variable ) /= 'biom_pet') RETURN
!-- Set averaging switch to .TRUE. if not done by other module before!
IF ( .NOT. ALLOCATED( pt_av ) ) THEN
ALLOCATE( pt_av(nzb:nzt+1,nys:nyn,nxl:nxr) )
aver_pt = .TRUE.
ENDIF
pt_av = 0.0_wp
IF ( .NOT. ALLOCATED( q_av ) ) THEN
ALLOCATE( q_av(nzb:nzt+1,nys:nyn,nxl:nxr) )
aver_q = .TRUE.
ENDIF
q_av = 0.0_wp
IF ( .NOT. ALLOCATED( u_av ) ) THEN
ALLOCATE( u_av(nzb:nzt+1,nys:nyn,nxl:nxr) )
aver_u = .TRUE.
ENDIF
u_av = 0.0_wp
IF ( .NOT. ALLOCATED( v_av ) ) THEN
ALLOCATE( v_av(nzb:nzt+1,nys:nyn,nxl:nxr) )
aver_v = .TRUE.
ENDIF
v_av = 0.0_wp
IF ( .NOT. ALLOCATED( w_av ) ) THEN
ALLOCATE( w_av(nzb:nzt+1,nys:nyn,nxl:nxr) )
aver_w = .TRUE.
ENDIF
w_av = 0.0_wp
CASE DEFAULT
CONTINUE
END SELECT
ELSEIF ( mode == 'sum' ) THEN
SELECT CASE ( TRIM( variable ) )
CASE ( 'biom_mrt' )
IF ( ALLOCATED( biom_mrt_av ) ) THEN
IF ( nmrtbl > 0 ) THEN
IF ( mrt_include_sw ) THEN
biom_mrt_av(:) = biom_mrt_av(:) + &
((human_absorb*mrtinsw(:) + human_emiss*mrtinlw(:)) &
/ (human_emiss*sigma_sb)) ** .25_wp - degc_to_k
ELSE
biom_mrt_av(:) = biom_mrt_av(:) + &
(human_emiss * mrtinlw(:) / sigma_sb) ** .25_wp &
- degc_to_k
ENDIF
ENDIF
ENDIF
CASE ( 'biom_pt', 'biom_utci', 'biom_pet' )
!-- Only continue if updateindex
IF ( average_trigger_pt .AND. TRIM( variable ) /= 'biom_pt') &
RETURN
IF ( average_trigger_utci .AND. TRIM( variable ) /= 'biom_utci') &
RETURN
IF ( average_trigger_pet .AND. TRIM( variable ) /= 'biom_pet') &
RETURN
IF ( ALLOCATED( pt_av ) .AND. aver_pt ) THEN
DO i = nxl, nxr
DO j = nys, nyn
DO k = nzb, nzt+1
pt_av(k,j,i) = pt_av(k,j,i) + pt(k,j,i)
ENDDO
ENDDO
ENDDO
ENDIF
IF ( ALLOCATED( q_av ) .AND. aver_q ) THEN
DO i = nxl, nxr
DO j = nys, nyn
DO k = nzb, nzt+1
q_av(k,j,i) = q_av(k,j,i) + q(k,j,i)
ENDDO
ENDDO
ENDDO
ENDIF
IF ( ALLOCATED( u_av ) .AND. aver_u ) THEN
DO i = nxl, nxr
DO j = nys, nyn
DO k = nzb, nzt+1
u_av(k,j,i) = u_av(k,j,i) + u(k,j,i)
ENDDO
ENDDO
ENDDO
ENDIF
IF ( ALLOCATED( v_av ) .AND. aver_v ) THEN
DO i = nxl, nxr
DO j = nys, nyn
DO k = nzb, nzt+1
v_av(k,j,i) = v_av(k,j,i) + v(k,j,i)
ENDDO
ENDDO
ENDDO
ENDIF
IF ( ALLOCATED( w_av ) .AND. aver_w ) THEN
DO i = nxl, nxr
DO j = nys, nyn
DO k = nzb, nzt+1
w_av(k,j,i) = w_av(k,j,i) + w(k,j,i)
ENDDO
ENDDO
ENDDO
ENDIF
CASE DEFAULT
CONTINUE
END SELECT
ELSEIF ( mode == 'average' ) THEN
SELECT CASE ( TRIM( variable ) )
CASE ( 'biom_mrt' )
IF ( ALLOCATED( biom_mrt_av ) ) THEN
biom_mrt_av(:) = biom_mrt_av(:) / REAL( average_count_3d, KIND=wp )
ENDIF
CASE ( 'biom_pt', 'biom_utci', 'biom_pet' )
!-- Only continue if update index
IF ( average_trigger_pt .AND. TRIM( variable ) /= 'biom_pt') &
RETURN
IF ( average_trigger_utci .AND. TRIM( variable ) /= 'biom_utci') &
RETURN
IF ( average_trigger_pet .AND. TRIM( variable ) /= 'biom_pet') &
RETURN
IF ( ALLOCATED( pt_av ) .AND. aver_pt ) THEN
DO i = nxl, nxr
DO j = nys, nyn
DO k = nzb, nzt+1
pt_av(k,j,i) = pt_av(k,j,i) / REAL( average_count_3d, KIND=wp )
ENDDO
ENDDO
ENDDO
ENDIF
IF ( ALLOCATED( q_av ) .AND. aver_q ) THEN
DO i = nxl, nxr
DO j = nys, nyn
DO k = nzb, nzt+1
q_av(k,j,i) = q_av(k,j,i) / REAL( average_count_3d, KIND=wp )
ENDDO
ENDDO
ENDDO
ENDIF
IF ( ALLOCATED( u_av ) .AND. aver_u ) THEN
DO i = nxl, nxr
DO j = nys, nyn
DO k = nzb, nzt+1
u_av(k,j,i) = u_av(k,j,i) / REAL( average_count_3d, KIND=wp )
ENDDO
ENDDO
ENDDO
ENDIF
IF ( ALLOCATED( v_av ) .AND. aver_v ) THEN
DO i = nxl, nxr
DO j = nys, nyn
DO k = nzb, nzt+1
v_av(k,j,i) = v_av(k,j,i) / REAL( average_count_3d, KIND=wp )
ENDDO
ENDDO
ENDDO
ENDIF
IF ( ALLOCATED( w_av ) .AND. aver_w ) THEN
DO i = nxl, nxr
DO j = nys, nyn
DO k = nzb, nzt+1
w_av(k,j,i) = w_av(k,j,i) / REAL( average_count_3d, KIND=wp )
ENDDO
ENDDO
ENDDO
ENDIF
!-- Udate thermal indices with derived averages
CALL biom_calculate_static_grid ( .TRUE. )
END SELECT
ENDIF
END SUBROUTINE biom_3d_data_averaging
!------------------------------------------------------------------------------!
! Description:
! ------------
!> Check data output for biometeorology model
!------------------------------------------------------------------------------!
SUBROUTINE biom_check_data_output( var, unit )
USE control_parameters, &
ONLY: data_output, message_string
IMPLICIT NONE
CHARACTER (LEN=*) :: unit !< The unit for the variable var
CHARACTER (LEN=*) :: var !< The variable in question
SELECT CASE ( TRIM( var ) )
CASE( 'biom_tmrt', 'biom_mrt', 'biom_pet', 'biom_pt', 'biom_utci' )
unit = 'degree_C'
CASE DEFAULT
unit = 'illegal'
END SELECT
IF ( unit /= 'illegal' ) THEN
!
!-- Futher checks if output belongs to biometeorology
IF ( .NOT. biometeorology .OR. .NOT. radiation ) THEN
message_string = 'output of "' // TRIM( var ) // '" req' // &
'uires biometeorology = .TRUE. and radiation = .TRUE. ' &
// 'in initialisation_parameters'
CALL message( 'biom_check_data_output', 'PA0561', 1, 2, 0, 6, 0 )
ENDIF
IF ( mrt_nlevels == 0 ) THEN
message_string = 'output of "' // TRIM( var ) // '" require'&
// 's mrt_nlevels > 0'
CALL message( 'biom_check_data_output', 'PA0562', 1, 2, 0, 6, 0 )
ENDIF
ENDIF
END SUBROUTINE biom_check_data_output
!------------------------------------------------------------------------------!
! Description:
! ------------
!> Check parameters routine for biom module
!> check_parameters 1302
!------------------------------------------------------------------------------!
SUBROUTINE biom_check_parameters
USE control_parameters, &
ONLY: message_string
IMPLICIT NONE
!-- Disabled as long as radiation model not available
IF ( .NOT. radiation ) THEN
message_string = 'The human thermal comfort module does require ' // &
'radiation information in terms of the mean ' // &
'radiant temperature, but radiation is not enabled!'
CALL message( 'check_parameters', 'PAHU01', 0, 0, 0, 6, 0 )
ENDIF
IF ( .NOT. humidity ) THEN
message_string = 'The human thermal comfort module does require ' // &
'air humidity information, but humidity module ' // &
'is disabled!'
CALL message( 'check_parameters', 'PAHU02', 0, 0, 0, 6, 0 )
ENDIF
END SUBROUTINE biom_check_parameters
!------------------------------------------------------------------------------!
!
! Description:
! ------------
!> Subroutine defining 3D output variables (dummy, always 2d!)
!> data_output_3d 709ff
!------------------------------------------------------------------------------!
SUBROUTINE biom_data_output_3d( av, variable, found, local_pf, nzb_do, nzt_do )
USE indices
USE kinds
IMPLICIT NONE
!-- Input variables
CHARACTER (LEN=*), INTENT(IN) :: variable !< Char identifier to select var for output
INTEGER(iwp), INTENT(IN) :: av !< Use averaged data? 0 = no, 1 = yes?
INTEGER(iwp), INTENT(IN) :: nzb_do !< Unused. 2D. nz bottom to nz top
INTEGER(iwp), INTENT(IN) :: nzt_do !< Unused.
!-- Output variables
LOGICAL, INTENT(OUT) :: found !< Output found?
REAL(sp), DIMENSION(nxl:nxr,nys:nyn,nzb_do:nzt_do) :: local_pf !< Temp. result grid to return
!-- Internal variables
INTEGER(iwp) :: l !< Running index, radiation grid
INTEGER(iwp) :: i !< Running index, x-dir
INTEGER(iwp) :: j !< Running index, y-dir
INTEGER(iwp) :: k !< Running index, z-dir
CHARACTER (LEN=:), allocatable :: variable_short !< Trimmed version of char variable
REAL(wp), PARAMETER :: fill_value = -999._wp
REAL(wp) :: mrt !< Buffer for mean radiant temperature
found = .TRUE.
variable_short = TRIM( variable )
IF ( variable_short(1:4) /= 'biom' ) THEN
!-- Nothing to do, set found to FALSE and return immediatelly
found = .FALSE.
RETURN
ENDIF
SELECT CASE ( variable_short )
CASE ( 'biom_mrt' )
local_pf = REAL( fill_value, KIND = wp )
DO l = 1, nmrtbl
i = mrtbl(ix,l)
j = mrtbl(iy,l)
k = mrtbl(iz,l)
IF ( mrt_include_sw ) THEN
mrt = ((human_absorb*mrtinsw(l) + human_emiss*mrtinlw(l)) &
/ (human_emiss*sigma_sb)) ** .25_wp
ELSE
mrt = (human_emiss*mrtinlw(l) / sigma_sb) ** .25_wp
ENDIF
local_pf(i,j,k) = mrt
ENDDO
CASE ( 'biom_tmrt' ) ! 2d-array
DO i = nxl, nxr
DO j = nys, nyn
local_pf(i,j,nzb_do) = REAL( tmrt_grid(j,i), sp )
ENDDO
ENDDO
CASE ( 'biom_pt' ) ! 2d-array
IF ( av == 0 ) THEN
DO i = nxl, nxr
DO j = nys, nyn
local_pf(i,j,nzb_do) = REAL( pt_grid(j,i), sp )
ENDDO
ENDDO
ELSE
DO i = nxl, nxr
DO j = nys, nyn
local_pf(i,j,nzb_do) = REAL( pt_av_grid(j,i), sp )
ENDDO
ENDDO
END IF
CASE ( 'biom_utci' ) ! 2d-array
IF ( av == 0 ) THEN
DO i = nxl, nxr
DO j = nys, nyn
local_pf(i,j,nzb_do) = REAL( utci_grid(j,i), sp )
ENDDO
ENDDO
ELSE
DO i = nxl, nxr
DO j = nys, nyn
local_pf(i,j,nzb_do) = REAL( utci_av_grid(j,i), sp )
ENDDO
ENDDO
END IF
CASE ( 'biom_pet' ) ! 2d-array
IF ( av == 0 ) THEN
DO i = nxl, nxr
DO j = nys, nyn
local_pf(i,j,nzb_do) = REAL( pet_grid(j,i), sp )
ENDDO
ENDDO
ELSE
DO i = nxl, nxr
DO j = nys, nyn
local_pf(i,j,nzb_do) = REAL( pet_av_grid(j,i), sp )
ENDDO
ENDDO
END IF
CASE DEFAULT
found = .FALSE.
END SELECT
END SUBROUTINE biom_data_output_3d
!------------------------------------------------------------------------------!
!
! Description:
! ------------
!> Subroutine defining 2D output variables
!> data_output_2d 1188ff
!------------------------------------------------------------------------------!
SUBROUTINE biom_data_output_2d( av, variable, found, grid, local_pf, &
two_d, nzb_do, nzt_do, fill_value )
USE indices, &
ONLY: nxl, nxl, nxr, nxr, nyn, nyn, nys, nys, nzb, nzt
USE kinds
IMPLICIT NONE
!-- Input variables
CHARACTER (LEN=*), INTENT(IN) :: variable !< Char identifier to select var for output
INTEGER(iwp), INTENT(IN) :: av !< Use averaged data? 0 = no, 1 = yes?
INTEGER(iwp), INTENT(IN) :: nzb_do !< Unused. 2D. nz bottom to nz top
INTEGER(iwp), INTENT(IN) :: nzt_do !< Unused.
REAL(wp), INTENT(in) :: fill_value !< Fill value for unassigned locations
!-- Output variables
CHARACTER (LEN=*), INTENT(OUT) :: grid !< Grid type (always "zu1" for biom)
LOGICAL, INTENT(OUT) :: found !< Output found?
LOGICAL, INTENT(OUT) :: two_d !< Flag parameter that indicates 2D variables, horizontal cross sections, must be .TRUE.
REAL(wp), DIMENSION(nxl:nxr,nys:nyn,nzb_do:nzt_do) :: local_pf !< Temp. result grid to return
!-- Internal variables
CHARACTER (LEN=:), allocatable :: variable_short !< Trimmed version of char variable
INTEGER(iwp) :: i !< Running index, x-dir
INTEGER(iwp) :: j !< Running index, y-dir
INTEGER(iwp) :: k !< Running index, z-dir
found = .TRUE.
variable_short = TRIM( variable )
IF ( variable_short(1:4) == 'biom' ) THEN
two_d = .TRUE.
grid = 'zu1'
ELSE
found = .FALSE.
grid = 'none'
RETURN
ENDIF
local_pf = fill_value
SELECT CASE ( variable_short )
CASE ( 'biom_tmrt_xy' ) ! 2d-array
DO i = nxl, nxr
DO j = nys, nyn
local_pf(i,j,1) = tmrt_grid(j,i)
ENDDO
ENDDO
CASE ( 'biom_pt_xy' ) ! 2d-array
IF ( av == 0 ) THEN
DO i = nxl, nxr
DO j = nys, nyn
local_pf(i,j,nzb+1) = pt_grid(j,i)
ENDDO
ENDDO
ELSE
DO i = nxl, nxr
DO j = nys, nyn
local_pf(i,j,nzb+1) = pt_av_grid(j,i)
ENDDO
ENDDO
END IF
CASE ( 'biom_utci_xy' ) ! 2d-array
IF ( av == 0 ) THEN
DO i = nxl, nxr
DO j = nys, nyn
local_pf(i,j,nzb+1) = utci_grid(j,i)
ENDDO
ENDDO
ELSE
DO i = nxl, nxr
DO j = nys, nyn
local_pf(i,j,nzb+1) = utci_av_grid(j,i)
ENDDO
ENDDO
END IF
CASE ( 'biom_pet_xy' ) ! 2d-array
IF ( av == 0 ) THEN
DO i = nxl, nxr
DO j = nys, nyn
local_pf(i,j,nzb+1) = pet_grid(j,i)
ENDDO
ENDDO
ELSE
DO i = nxl, nxr
DO j = nys, nyn
local_pf(i,j,nzb+1) = pet_av_grid(j,i)
ENDDO
ENDDO
END IF
CASE DEFAULT
found = .FALSE.
grid = 'none'
END SELECT
END SUBROUTINE biom_data_output_2d
!------------------------------------------------------------------------------!
! Description:
! ------------
!> Subroutine defining appropriate grid for netcdf variables.
!> It is called out from subroutine netcdf_interface_mod.
!> netcdf_interface_mod 918ff
!------------------------------------------------------------------------------!
SUBROUTINE biom_define_netcdf_grid( var, found, grid_x, grid_y, grid_z )
IMPLICIT NONE
!-- Input variables
CHARACTER (LEN=*), INTENT(IN) :: var !< Name of output variable
!-- Output variables
CHARACTER (LEN=*), INTENT(OUT) :: grid_x !< x grid of output variable
CHARACTER (LEN=*), INTENT(OUT) :: grid_y !< y grid of output variable
CHARACTER (LEN=*), INTENT(OUT) :: grid_z !< z grid of output variable
LOGICAL, INTENT(OUT) :: found !< Flag if output var is found
!-- Local variables
LOGICAL :: is2d !< Var is 2d?
INTEGER(iwp) :: l !< Length of the var array
found = .FALSE.
grid_x = 'none'
grid_y = 'none'
grid_z = 'none'
l = MAX( 2, LEN_TRIM( var ) )
is2d = ( var(l-1:l) == 'xy' )
IF ( var(1:4) == 'biom' ) THEN
found = .TRUE.
grid_x = 'x'
grid_y = 'y'
grid_z = 'zu'
IF ( is2d ) grid_z = 'zu1'
ENDIF
END SUBROUTINE biom_define_netcdf_grid
!------------------------------------------------------------------------------!
! Description:
! ------------
!> Header output for biom module
!> header 982
!------------------------------------------------------------------------------!
SUBROUTINE biom_header( io )
IMPLICIT NONE
!-- Input variables
INTEGER(iwp), INTENT(IN) :: io !< Unit of the output file
!-- Internal variables
CHARACTER (LEN=86) :: output_height_chr !< String for output height
WRITE( output_height_chr, '(F8.1,7X)' ) biom_output_height
!
!-- Write biom header
WRITE( io, 1 )
WRITE( io, 2 ) TRIM( output_height_chr )
WRITE( io, 3 ) TRIM( ACHAR( biom_cell_level ) )
1 FORMAT (//' Human thermal comfort module information:'/ &
' ------------------------------'/)
2 FORMAT (' --> All indices calculated for a height of (m): ', A )
3 FORMAT (' --> This corresponds to cell level : ', A )
END SUBROUTINE biom_header
!------------------------------------------------------------------------------!
! Description:
! ------------
!> Initialization of the HTCM
!> init_3d_model 1987ff
!------------------------------------------------------------------------------!
SUBROUTINE biom_init
USE control_parameters, &
ONLY: message_string
IMPLICIT NONE
!-- Internal vriables
REAL ( wp ) :: height !< current height in meters
INTEGER ( iwp ) :: i !< iteration index
!-- Determine cell level corresponding to 1.1 m above ground level
! (gravimetric center of sample human)
time_biom_results = 0.0_wp
biom_cell_level = 0_iwp
biom_output_height = 0.5_wp * dz(1)
height = 0.0_wp
biom_cell_level = INT ( 1.099_wp / dz(1) )
biom_output_height = biom_output_height + biom_cell_level * dz(1)
IF ( .NOT. radiation_interactions .AND. mrt_from_rtm ) THEN
message_string = 'The mrt_from_rtm switch require ' // &
'enabled radiation_interactions but it ' // &
'is disabled! Set mrt_from_rtm to .F.'
CALL message( 'check_parameters', 'PAHU03', 0, 0, -1, 6, 0 )
mrt_from_rtm = .FALSE.
ENDIF
END SUBROUTINE biom_init
!------------------------------------------------------------------------------!
! Description:
! ------------
!> Allocate biom arrays and define pointers if required
!> init_3d_model 1047ff
!------------------------------------------------------------------------------!
SUBROUTINE biom_init_arrays
IMPLICIT NONE
!-- Allocate a temporary array with the desired output dimensions.
! Initialization omitted for performance, will be overwritten anyway
IF ( .NOT. ALLOCATED( tmrt_grid ) ) THEN
ALLOCATE( tmrt_grid (nys:nyn,nxl:nxr) )
ENDIF
IF ( biom_pt ) THEN
IF ( .NOT. ALLOCATED( pt_grid ) ) THEN
ALLOCATE( pt_grid (nys:nyn,nxl:nxr) )
ENDIF
ENDIF
IF ( biom_utci ) THEN
IF ( .NOT. ALLOCATED( utci_grid ) ) THEN
ALLOCATE( utci_grid (nys:nyn,nxl:nxr) )
ENDIF
ENDIF
IF ( biom_pet ) THEN
IF ( .NOT. ALLOCATED( pet_grid ) ) THEN
ALLOCATE( pet_grid (nys:nyn,nxl:nxr) )
ENDIF
END IF
IF ( biom_pt_av ) THEN
IF ( .NOT. ALLOCATED( pt_av_grid ) ) THEN
ALLOCATE( pt_av_grid (nys:nyn,nxl:nxr) )
ENDIF
ENDIF
IF ( biom_utci_av ) THEN
IF ( .NOT. ALLOCATED( utci_av_grid ) ) THEN
ALLOCATE( utci_av_grid (nys:nyn,nxl:nxr) )
ENDIF
ENDIF
IF ( biom_pet_av ) THEN
IF ( .NOT. ALLOCATED( pet_av_grid ) ) THEN
ALLOCATE( pet_av_grid (nys:nyn,nxl:nxr) )
ENDIF
END IF
END SUBROUTINE biom_init_arrays
!------------------------------------------------------------------------------!
! Description:
! ------------
!> Parin for &biometeorology_parameters for reading biomet parameters
!------------------------------------------------------------------------------!
SUBROUTINE biom_parin
IMPLICIT NONE
!
!-- Internal variables
CHARACTER (LEN=80) :: line !< Dummy string for current line in parameter file
NAMELIST /biometeorology_parameters/ mrt_from_rtm, &
biom_pet, &
biom_pet_av, &
biom_pt, &
biom_pt_av, &
biom_utci, &
biom_utci_av
!-- Try to find biometeorology_parameters namelist
REWIND ( 11 )
line = ' '
DO WHILE ( INDEX( line, '&biometeorology_parameters' ) == 0 )
READ ( 11, '(A)', END = 20 ) line
ENDDO
BACKSPACE ( 11 )
!
!-- Read biometeorology_parameters namelist
READ ( 11, biometeorology_parameters, ERR = 10, END = 20 )
!
!-- Set flag that indicates that the biomet_module is switched on
biometeorology = .TRUE.
GOTO 20
!
!-- In case of error
10 BACKSPACE( 11 )
READ( 11 , '(A)') line
CALL parin_fail_message( 'biometeorology_parameters', line )
!
!-- Complete
20 CONTINUE
END SUBROUTINE biom_parin
!------------------------------------------------------------------------------!
! Description:
! ------------
!> Calculates the mean radiant temperature (tmrt) based on the Six-directions
!> method according to VDI 3787 2.
!------------------------------------------------------------------------------!
SUBROUTINE calculate_tmrt_6_directions( SW_N, SW_E, SW_S, SW_W, &
SW_U, SW_D, LW_N, LW_E, LW_S, LW_W, LW_U, LW_D, tmrt )
IMPLICIT NONE
!-- Type of input of the argument list
! Short- (SW_) and longwave (LW_) radiation fluxes from the six directions
! North (N), East (E), South (S), West (W), up (U) and down (D)
REAL(wp), INTENT ( IN ) :: SW_N !< Sw radflux density from N (W/m²)
REAL(wp), INTENT ( IN ) :: SW_E !< Sw radflux density from E (W/m²)
REAL(wp), INTENT ( IN ) :: SW_S !< Sw radflux density from S (W/m²)
REAL(wp), INTENT ( IN ) :: SW_W !< Sw radflux density from W (W/m²)
REAL(wp), INTENT ( IN ) :: SW_U !< Sw radflux density from U (W/m²)
REAL(wp), INTENT ( IN ) :: SW_D !< Sw radflux density from D (W/m²)
REAL(wp), INTENT ( IN ) :: LW_N !< Lw radflux density from N (W/m²)
REAL(wp), INTENT ( IN ) :: LW_E !< Lw radflux density from E (W/m²)
REAL(wp), INTENT ( IN ) :: LW_S !< Lw radflux density from S (W/m²)
REAL(wp), INTENT ( IN ) :: LW_W !< Lw radflux density from W (W/m²)
REAL(wp), INTENT ( IN ) :: LW_U !< Lw radflux density from U (W/m²)
REAL(wp), INTENT ( IN ) :: LW_D !< Lw radflux density from D (W/m²)
!-- Type of output of the argument list
REAL(wp), INTENT ( OUT ) :: tmrt !< Mean radiant temperature (°C)
!-- Directional weighting factors
REAL(wp), PARAMETER :: weight_h = 0.22_wp
REAL(wp), PARAMETER :: weight_v = 0.06_wp
REAL(wp) :: nrfd !< Net radiation flux density (W/m²)
!-- Initialise
nrfd = 0._wp
!-- Compute mean radiation flux density absorbed by rotational symetric human
nrfd = ( weight_h * ( human_absorb * SW_N + human_emiss * LW_N ) ) + &
( weight_h * ( human_absorb * SW_E + human_emiss * LW_E ) ) + &
( weight_h * ( human_absorb * SW_S + human_emiss * LW_S ) ) + &
( weight_h * ( human_absorb * SW_W + human_emiss * LW_W ) ) + &
( weight_v * ( human_absorb * SW_U + human_emiss * LW_U ) ) + &
( weight_v * ( human_absorb * SW_D + human_emiss * LW_D ) )
!-- Compute mean radiant temperature
tmrt = ( nrfd / (human_emiss * sigma_sb) )**0.25_wp - degc_to_k
END SUBROUTINE calculate_tmrt_6_directions
!------------------------------------------------------------------------------!
! Description:
! ------------
!> Very crude approximation of mean radiant temperature based on upwards and
!> downwards radiation fluxes only (other directions curently not available,
!> replace as soon as possible!)
!------------------------------------------------------------------------------!
SUBROUTINE calculate_tmrt_2_directions( sw_u, sw_d, lw_u, lw_d, ta, tmrt )
IMPLICIT NONE
!-- Type of input of the argument list
REAL(wp), INTENT ( IN ) :: sw_u !< Shortwave radiation flux density from upper direction (W/m²)
REAL(wp), INTENT ( IN ) :: sw_d !< Shortwave radiation flux density from lower direction (W/m²)
REAL(wp), INTENT ( IN ) :: lw_u !< Longwave radiation flux density from upper direction (W/m²)
REAL(wp), INTENT ( IN ) :: lw_d !< Longwave radiation flux density from lower direction (W/m²)
REAL(wp), INTENT ( IN ) :: ta !< Air temperature (°C)
!-- Type of output of the argument list
REAL(wp), INTENT ( OUT ) :: tmrt !< mean radiant temperature, (°C)
!-- Directional weighting factors and parameters
REAL(wp), PARAMETER :: weight_h = 0.22_wp !< Weight for horizontal radiational gain after Fanger (1972)
REAL(wp), PARAMETER :: weight_v = 0.06_wp !< Weight for vertical radiational gain after Fanger (1972)
!-- Other internal variables
REAL(wp) :: sw_in
REAL(wp) :: sw_out
REAL(wp) :: lw_in
REAL(wp) :: lw_out
REAL(wp) :: nrfd !< Net radiation flux density (W/m²)
REAL(wp) :: lw_air !< Longwave emission by surrounding air volume (W/m²)
REAL(wp) :: sw_side !< Shortwave immission from the sides (W/m²)
INTEGER(iwp) :: no_input !< Count missing input radiation fluxes
!-- initialise
sw_in = sw_u
sw_out = sw_d
lw_in = lw_u
lw_out = lw_d
nrfd = 0._wp
no_input = 0_iwp
!-- test for missing input data
IF ( sw_in <= -998._wp .OR. sw_out <= -998._wp .OR. lw_in <= -998._wp .OR. &
lw_out <= -998._wp .OR. ta <= -998._wp ) THEN
IF ( sw_in <= -998._wp ) THEN
sw_in = 0.
no_input = no_input + 1
ENDIF
IF ( sw_out <= -998._wp ) THEN
sw_out = 0.
no_input = no_input + 1
ENDIF
IF ( lw_in <= -998._wp ) THEN
lw_in = 0.
no_input = no_input + 1
ENDIF
IF ( lw_out <= -998._wp ) THEN
lw_out = 0.
no_input = no_input + 1
ENDIF
!-- Accept two missing radiation flux directions, fail otherwise as error might become too large
IF ( ta <= -998._wp .OR. no_input >= 3 ) THEN
tmrt = -999._wp
RETURN
ENDIF
ENDIF
sw_side = sw_in * 0.125_wp ! distribute half of upper sw_in to the 4 sides
lw_air = ( sigma_sb * 0.95_wp * ( ta + degc_to_k )**4 )
!-- Compute mean radiation flux density absorbed by rotational symetric human
nrfd = ( weight_h * ( human_absorb * sw_side + human_emiss * lw_air ) ) + &
( weight_h * ( human_absorb * sw_side + human_emiss * lw_air ) ) + &
( weight_h * ( human_absorb * sw_side + human_emiss * lw_air ) ) + &
( weight_h * ( human_absorb * sw_side + human_emiss * lw_air ) ) + &
( weight_v * ( human_absorb * (sw_in * 0.5_wp) + human_emiss * lw_in ) ) + &
( weight_v * ( human_absorb * sw_out + human_emiss * lw_out ) )
!-- Compute mean radiant temperature
tmrt = ( nrfd / (human_emiss * sigma_sb) )**0.25_wp - degc_to_k
END SUBROUTINE calculate_tmrt_2_directions
!------------------------------------------------------------------------------!
! Description:
! ------------
!> Calculate static thermal indices for given meteorological conditions
!------------------------------------------------------------------------------!
SUBROUTINE calculate_static_thermal_indices ( ta, vp, ws, pair, tmrt, &
pt_static, utci_static, pet_static )
IMPLICIT NONE
!-- Input parameters
REAL(wp), INTENT ( IN ) :: ta !< Air temperature (°C)
REAL(wp), INTENT ( IN ) :: vp !< Vapour pressure (hPa)
REAL(wp), INTENT ( IN ) :: ws !< Wind speed (local level) (m/s)
REAL(wp), INTENT ( IN ) :: pair !< Air pressure (hPa)
REAL(wp), INTENT ( IN ) :: tmrt !< Mean radiant temperature (°C)
!-- Output parameters
REAL(wp), INTENT ( OUT ) :: pt_static !< Perceived temperature (°C)
REAL(wp), INTENT ( OUT ) :: utci_static !< Universal thermal climate index (°C)
REAL(wp), INTENT ( OUT ) :: pet_static !< Physiologically equivalent temp. (°C)
!-- Temporary field, not used here
REAL(wp) :: clo !< Clothing index (no dim.)
clo = -999._wp
IF ( biom_pt ) THEN
!-- Estimate local perceived temperature
CALL calculate_pt_static( ta, vp, ws, tmrt, pair, clo, pt_static )
ENDIF
IF ( biom_utci ) THEN
!-- Estimate local universal thermal climate index
CALL calculate_utci_static( ta, vp, ws, tmrt, biom_output_height, &
utci_static )
ENDIF
IF ( biom_pet ) THEN
!-- Estimate local physiologically equivalent temperature
CALL calculate_pet_static( ta, vp, ws, tmrt, pair, pet_static )
ENDIF
END SUBROUTINE calculate_static_thermal_indices
!------------------------------------------------------------------------------!
! Description:
! ------------
!> Calculate static thermal indices for 2D grid point i, j
!------------------------------------------------------------------------------!
SUBROUTINE biom_determine_input_at( average_input, i, j, ta, vp, ws, pair, &
tmrt )
IMPLICIT NONE
!-- Input variables
LOGICAL, INTENT ( IN ) :: average_input !< Determine averaged input conditions?
INTEGER(iwp), INTENT ( IN ) :: i !< Running index, x-dir
INTEGER(iwp), INTENT ( IN ) :: j !< Running index, y-dir
!-- Output parameters
REAL(wp), INTENT ( OUT ) :: tmrt !< Mean radiant temperature (°C)
REAL(wp), INTENT ( OUT ) :: ta !< Air temperature (°C)
REAL(wp), INTENT ( OUT ) :: vp !< Vapour pressure (hPa)
REAL(wp), INTENT ( OUT ) :: ws !< Wind speed (local level) (m/s)
REAL(wp), INTENT ( OUT ) :: pair !< Air pressure (hPa)
!-- Internal variables
INTEGER(iwp) :: k !< Running index, z-dir
INTEGER(iwp) :: k_wind !< Running index, z-dir, wind speed only
REAL(wp) :: vp_sat !< Saturation vapor pressure (hPa)
!-- Determine cell level closest to 1.1m above ground
! by making use of truncation due to int cast
k = get_topography_top_index_ji(j, i, 's') + biom_cell_level !< Vertical cell center closest to 1.1m
k_wind = k
IF( k_wind < 1_iwp ) THEN ! Avoid horizontal u and v of 0.0 m/s close to ground
k_wind = 1_iwp
ENDIF
!-- Determine local values:
IF ( average_input ) THEN
!-- Calculate ta from Tp assuming dry adiabatic laps rate
ta = pt_av(k, j, i) - ( 0.0098_wp * dz(1) * ( k + .5_wp ) ) - degc_to_k
vp = 0.034_wp
IF ( humidity .AND. ALLOCATED( q_av ) ) THEN
vp = q_av(k, j, i)
ENDIF
ws = ( 0.5_wp * ABS( u_av(k_wind, j, i) + u_av(k_wind, j, i+1) ) + &
0.5_wp * ABS( v_av(k_wind, j, i) + v_av(k_wind, j+1, i) ) + &
0.5_wp * ABS( w_av(k_wind, j, i) + w_av(k_wind+1, j, i) ) )
ELSE
!-- Calculate ta from Tp assuming dry adiabatic laps rate
ta = pt(k, j, i) - ( 0.0098_wp * dz(1) * ( k + .5_wp ) ) - degc_to_k
vp = q(k, j, i)
ws = ( 0.5_wp * ABS( u(k_wind, j, i) + u(k_wind, j, i+1) ) + &
0.5_wp * ABS( v(k_wind, j, i) + v(k_wind, j+1, i) ) + &
0.5_wp * ABS( w(k_wind, j, i) + w(k_wind+1, j, i) ) )
ENDIF
!-- Local air pressure
pair = surface_pressure
!
!-- Calculate water vapour pressure at saturation and convert to hPa
!-- The magnus formula is limited to temperatures up to 333.15 K to
! avoid negative values of vp_sat
vp_sat = 0.01_wp * magnus( MIN( ta + degc_to_k, 333.15_wp ) )
vp = vp * pair / ( vp + 0.622_wp )
IF ( vp > vp_sat ) vp = vp_sat
tmrt = ta
IF ( radiation ) THEN
IF ( mrt_from_rtm ) THEN
tmrt = tmrt_grid(j, i)
ELSE
CALL calculate_tmrt_2_directions (rad_sw_in(0, j, i), &
rad_sw_out(0, j, i), rad_lw_in(0, j, i), rad_lw_out(0, j, i), ta, &
tmrt )
ENDIF
ENDIF
END SUBROUTINE biom_determine_input_at
!------------------------------------------------------------------------------!
! Description:
! ------------
!> Calculate static thermal indices for any point within a 2D grid
!> time_integration.f90: 1065ff
!------------------------------------------------------------------------------!
SUBROUTINE biom_calculate_static_grid ( average_input )
IMPLICIT NONE
!-- Input attributes
LOGICAL, INTENT ( IN ) :: average_input !< Calculate based on averaged input? conditions?
!-- Internal variables
INTEGER(iwp) :: i, j, k, l !< Running index
REAL(wp) :: ta !< Air temperature (°C)
REAL(wp) :: vp !< Vapour pressure (hPa)
REAL(wp) :: ws !< Wind speed (local level) (m/s)
REAL(wp) :: pair !< Air pressure (hPa)
REAL(wp) :: tmrt_tmp !< Mean radiant temperature
REAL(wp) :: pt_tmp !< Perceived temperature
REAL(wp) :: utci_tmp !< Universal thermal climate index
REAL(wp) :: pet_tmp !< Physiologically equivalent temperature
CALL biom_init_arrays
IF ( mrt_from_rtm ) THEN
tmrt_grid = -999._wp
DO l = 1, nmrtbl
i = mrtbl(ix,l)
j = mrtbl(iy,l)
k = mrtbl(iz,l)
IF ( k - get_topography_top_index_ji( j, i, 's' ) == 1 ) THEN
IF ( mrt_include_sw ) THEN
tmrt_tmp = ((human_absorb*mrtinsw(l) + human_emiss*mrtinlw(l)) &
/ (human_emiss*sigma_sb)) ** .25_wp
ELSE
tmrt_tmp = (human_emiss*mrtinlw(l) / sigma_sb) ** .25_wp
ENDIF
tmrt_grid(j, i) = tmrt_tmp - degc_to_k
ENDIF
ENDDO
ENDIF
DO i = nxl, nxr
DO j = nys, nyn
!-- Determine local input conditions
CALL biom_determine_input_at ( average_input, i, j, ta, vp, ws, &
pair, tmrt_tmp )
tmrt_grid(j, i) = tmrt_tmp
!-- Only proceed if tmrt is available
IF ( tmrt_tmp <= -998._wp ) THEN
pt_tmp = -999._wp
utci_tmp = -999._wp
pet_tmp = -999._wp
CYCLE
END IF
!-- Calculate static thermal indices based on local tmrt
CALL calculate_static_thermal_indices ( ta, vp, ws, &
pair, tmrt_tmp, pt_tmp, utci_tmp, pet_tmp )
IF ( average_input ) THEN
!-- Write results for selected averaged indices
IF ( biom_pt_av ) THEN
pt_av_grid(j, i) = pt_tmp
END IF
IF ( biom_utci_av ) THEN
utci_av_grid(j, i) = utci_tmp
END IF
IF ( biom_pet_av ) THEN
pet_av_grid(j, i) = pet_tmp
END IF
ELSE
!-- Write result for selected indices
IF ( biom_pt ) THEN
pt_grid(j, i) = pt_tmp
END IF
IF ( biom_utci ) THEN
utci_grid(j, i) = utci_tmp
END IF
IF ( biom_pet ) THEN
pet_grid(j, i) = pet_tmp
END IF
END IF
END DO
END DO
END SUBROUTINE biom_calculate_static_grid
!------------------------------------------------------------------------------!
! Description:
! ------------
!> Calculate dynamic thermal indices (currently only iPT, but expandable)
!------------------------------------------------------------------------------!
SUBROUTINE biom_calc_ipt( ta, vp, ws, pair, tmrt, dt, energy_storage, &
t_clo, clo, actlev, age, weight, height, work, sex, ipt )
IMPLICIT NONE
!-- Input parameters
REAL(wp), INTENT ( IN ) :: ta !< Air temperature (°C)
REAL(wp), INTENT ( IN ) :: vp !< Vapour pressure (hPa)
REAL(wp), INTENT ( IN ) :: ws !< Wind speed (local level) (m/s)
REAL(wp), INTENT ( IN ) :: pair !< Air pressure (hPa)
REAL(wp), INTENT ( IN ) :: tmrt !< Mean radiant temperature (°C)
REAL(wp), INTENT ( IN ) :: dt !< Time past since last calculation (s)
REAL(wp), INTENT ( IN ) :: age !< Age of agent (y)
REAL(wp), INTENT ( IN ) :: weight !< Weight of agent (Kg)
REAL(wp), INTENT ( IN ) :: height !< Height of agent (m)
REAL(wp), INTENT ( IN ) :: work !< Mechanical workload of agent
! (without metabolism!) (W)
INTEGER(iwp), INTENT ( IN ) :: sex !< Sex of agent (1 = male, 2 = female)
!-- Both, input and output parameters
Real(wp), INTENT ( INOUT ) :: energy_storage !< Energy storage (W/m²)
Real(wp), INTENT ( INOUT ) :: t_clo !< Clothing temperature (°C)
Real(wp), INTENT ( INOUT ) :: clo !< Current clothing in sulation
Real(wp), INTENT ( INOUT ) :: actlev !< Individuals activity level
! per unit surface area (W/m²)
!-- Output parameters
REAL(wp), INTENT ( OUT ) :: ipt !< Instationary perceived temp. (°C)
!-- If clo equals the initial value, this is the initial call
IF ( clo <= -998._wp ) THEN
!-- Initialize instationary perceived temperature with personalized
! PT as an initial guess, set actlev and clo
CALL ipt_init ( age, weight, height, sex, work, actlev, clo, &
ta, vp, ws, tmrt, pair, dt, energy_storage, t_clo, &
ipt )
ELSE
!-- Estimate local instatinoary perceived temperature
CALL ipt_cycle ( ta, vp, ws, tmrt, pair, dt, energy_storage, t_clo, &
clo, actlev, work, ipt )
ENDIF
END SUBROUTINE biom_calc_ipt
END MODULE biometeorology_mod