Home


About us
PALM
PALM-4U
Installation
Release notes
Documentation
Getting Help
Gallery
Developers area

Changeset 3448 for palm/trunk


Ignore:
Timestamp:
Oct 29, 2018 6:14:31 PM (6 years ago)
Author:
kanani
Message:

Implementation of human thermal indices (from branch biomet_p2 at r3444)

Location:
palm/trunk/SOURCE
Files:
4 added
15 edited

Legend:

Unmodified
Added
Removed

  • palm/trunk/SOURCE/Makefile

    r3436 r3448  
    2525# -----------------
    2626# $Id$
     27# Adjustment of biometeorology dependencies
     28#
     29# 3436 2018-10-26 18:35:15Z gronemeier
    2730# Add surface_mod to user_data_output_mask
    2831#
    … …  
    540543        header.f90 \
    541544        biometeorology_mod.f90 \
     545        biometeorology_ipt_mod.f90 \
     546        biometeorology_pt_mod.f90 \
     547        biometeorology_pet_mod.f90 \
     548        biometeorology_utci_mod.f90 \
    542549        inflow_turbulence.f90 \
    543550        init_3d_model.f90 \
    … …  
    755762basic_constants_and_equations_mod.o: \
    756763        mod_kinds.o
     764biometeorology_ipt_mod.o: \
     765        biometeorology_pt_mod.o \
     766        mod_kinds.o
    757767biometeorology_mod.o: \
    758         mod_kinds.o \
    759         radiation_model_mod.o
     768        modules.o \
     769        mod_kinds.o \
     770        radiation_model_mod.o \
     771        biometeorology_ipt_mod.o \
     772        biometeorology_pt_mod.o \
     773        biometeorology_pet_mod.o \
     774        biometeorology_utci_mod.o
     775biometeorology_pt_mod.o: \
     776        mod_kinds.o
     777biometeorology_pet_mod.o: \
     778        mod_kinds.o
     779biometeorology_utci_mod.o: \
     780        mod_kinds.o
    760781boundary_conds.o: \
    761782        basic_constants_and_equations_mod.o \
    … …  
    919940data_output_2d.o: \
    920941        basic_constants_and_equations_mod.o \
     942        biometeorology_mod.o \
    921943        bulk_cloud_model_mod.o \
    922944        chemistry_model_mod.o\
    … …  
    10141036header.o: \
    10151037        basic_constants_and_equations_mod.o \
     1038        biometeorology_mod.o \
    10161039        bulk_cloud_model_mod.o \
    10171040        chemistry_model_mod.o \
    … …  
    10421065        advec_ws.o \
    10431066        basic_constants_and_equations_mod.o \
     1067        biometeorology_mod.o \
    10441068        bulk_cloud_model_mod.o \
    10451069        chem_emissions_mod.o \
    … …  
    12831307multi_agent_system_mod.o: \
    12841308        basic_constants_and_equations_mod.o \
     1309        biometeorology_mod.o \
    12851310        chemistry_model_mod.o \
    12861311        cpulog_mod.o \
    … …  
    13151340        uv_exposure_model_mod.o
    13161341nesting_offl_mod.o: \
    1317         cpulog_mod.o \
    1318         mod_kinds.o \
     1342        cpulog_mod.o \
     1343        mod_kinds.o \
    13191344        modules.o \
    13201345        netcdf_data_input_mod.o \
    … …  
    13571382        write_restart_data_mod.o
    13581383parin.o: \
     1384        biometeorology_mod.o \
    13591385        bulk_cloud_model_mod.o \
    13601386        chemistry_model_mod.o \
    … …  
    16231649time_integration.o: \
    16241650        advec_ws.o \
     1651        biometeorology_mod.o \
    16251652        bulk_cloud_model_mod.o \
    16261653        buoyancy.o \

  • palm/trunk/SOURCE/average_3d_data.f90

    r3421 r3448  
    2525! -----------------
    2626! $Id$
     27! Adjustment of biometeorology calls
     28!
     29! 3421 2018-10-24 18:39:32Z gronemeier
    2730! Renamed output variables
    2831!
    … …  
    153156    USE averaging
    154157
     158    USE biometeorology_mod,                                                    &
     159        ONLY:  biom_3d_data_averaging
     160
    155161    USE bulk_cloud_model_mod,                                                  &
    156162        ONLY:  bulk_cloud_model, bcm_3d_data_averaging
    … …  
    160166
    161167    USE control_parameters,                                                    &
    162         ONLY:  air_chemistry, average_count_3d, doav, doav_n, land_surface,    &
    163                ocean_mode, urban_surface, varnamelength
     168        ONLY:  air_chemistry, average_count_3d, biometeorology, doav, doav_n,  &
     169               land_surface, ocean_mode, urban_surface, varnamelength
    164170
    165171    USE cpulog,                                                                &
    … …  
    182188    USE radiation_model_mod,                                                   &
    183189        ONLY:  radiation, radiation_3d_data_averaging
    184 
    185     USE biometeorology_mod,                                                    &
    186         ONLY:  biometeorology_3d_data_averaging
    187190
    188191    USE turbulence_closure_mod,                                                &
    … …  
    567570             ENDIF
    568571
    569              IF ( radiation  .AND.  trimvar(1:4) == 'bio_' )  THEN
    570                 CALL biometeorology_3d_data_averaging( 'average', doav(ii) )
     572             IF ( biometeorology )  THEN
     573                CALL biom_3d_data_averaging( 'average', doav(ii) )
    571574             ENDIF
    572575

  • palm/trunk/SOURCE/biometeorology_mod.f90

    • Property svn:keywords set to Id
    r3361 r3448  
    1 !> @file biometeorology.f90
    2 !------------------------------------------------------------------------------!
     1!> @file biometeorology_mod.f90
     2!--------------------------------------------------------------------------------!
    33! This file is part of PALM-4U.
    44!
    … …  
    1515! PALM. If not, see <http://www.gnu.org/licenses/>.
    1616!
    17 ! Copyright 2018, Institute of Computer Science,Academy of Sciences, Prague
    18 !
    19 ! Calculation of PET:
    20 ! Copyright 2016, Deutscher Wetterdienst (DWD) /
    21 ! German Meteorological Service (DWD)
    22 !------------------------------------------------------------------------------!
     17! Copyright 2018-2018 Deutscher Wetterdienst (DWD)
     18! Copyright 2018-2018 Institute of Computer Science, Academy of Sciences, Prague
     19! Copyright 2018-2018 Leibniz Universitaet Hannover
     20!--------------------------------------------------------------------------------!
    2321!
    2422! Current revisions:
    … …  
    3533! Authors:
    3634! --------
     35! @author Dominik Froehlich <dominik.froehlich@dwd.de>
    3736! @author Jaroslav Resler <resler@cs.cas.cz>
    38 ! @author Dominik Froehlich <dominik.froehlich@dwd.de>
    39 !
    40 !------------------------------------------------------------------------------!
    41 
    42 MODULE biometeorology_mod
    43 !
    44 !-- Load required variables from existing modules
     37!
     38!
     39! Description:
     40! ------------
     41!> Human thermal comfort module calculating thermal perception of a sample
     42!> human being under the current meteorological conditions.
     43!>
     44!> @todo Alphabetical sorting of "USE ..." lists, "ONLY" list, variable declarations
     45!>       (per subroutine: first all CHARACTERs, then INTEGERs, LOGICALs, REALs, )
     46!> @todo Comments start with capital letter --> "!-- Include..."
     47!> @todo Variable and routine names strictly in lowercase letters and in English
     48!>
     49!> @note nothing now
     50!>
     51!> @bug  no known bugs by now
     52!------------------------------------------------------------------------------!
     53 MODULE biometeorology_mod
     54
    4555    USE arrays_3d,                                                             &
    46         ONLY:  hyp, p, pt, q, u, v, w
     56       ONLY:  pt, p, u, v, w, q
     57
     58    USE averaging,                                                             &
     59       ONLY:  pt_av, q_av, u_av, v_av, w_av
    4760
    4861    USE basic_constants_and_equations_mod,                                     &
    49         ONLY:  rd_d_rv
    50 
    51     USE kinds  !< to set precision of INTEGER and REAL arrays according to PALM
    52 
    53     USE radiation_model_mod,  &
    54         ONLY: ix, iy, iz, id, mrt_nlevels, mrt_include_sw,                     &
    55               mrtinsw, mrtinlw, mrtbl, nmrtbl, radiation
    56 
    57 
    58     IMPLICIT NONE
    59 
    60     REAL(wp), DIMENSION(:), ALLOCATABLE            ::  bio_mrt              !< biometeorology mean radiant temperature for each MRT box
    61     REAL(wp), DIMENSION(:), ALLOCATABLE            ::  bio_mrt_av           !< time average mean radiant temperature for each MRT box
    62     REAL(wp), DIMENSION(:), ALLOCATABLE            ::  bio_pet              !< PhysiologiCALLy Equivalent Temperature (PET) for each MRT box
    63     REAL(wp), DIMENSION(:), ALLOCATABLE            ::  bio_pet_av           !< time average PhysiologiCALLy Equivalent Temperature (PET) for each MRT box
    64 
    65 
    66     REAL(wp), PARAMETER :: sigma_sb       = 5.67037321E-8_wp,       & !< Stefan-Boltzmann constant
    67                            t_zero = -273.15_wp,                     & !< temperature 0K in Celsius
    68                            human_absorb = 0.7_wp,                   & !< SW absorbtivity of human body
    69                            human_emiss = 0.97_wp                      !< emissivity of human body (Lindberg 2008)
    70 
    71 
    72     INTERFACE biometeorology_3d_data_averaging
    73        MODULE PROCEDURE biometeorology_3d_data_averaging
    74     END INTERFACE biometeorology_3d_data_averaging
    75 
    76     INTERFACE biometeorology_check_data_output
    77        MODULE PROCEDURE biometeorology_check_data_output
    78     END INTERFACE biometeorology_check_data_output
    79 
    80     INTERFACE biometeorology_data_output_3d
    81        MODULE PROCEDURE biometeorology_data_output_3d
    82     END INTERFACE biometeorology_data_output_3d
    83 
    84     INTERFACE biometeorology_define_netcdf_grid
    85        MODULE PROCEDURE biometeorology_define_netcdf_grid
    86     END INTERFACE biometeorology_define_netcdf_grid
    87 
    88 
    89     SAVE
     62       ONLY:  magnus
     63
     64    USE biometeorology_ipt_mod
     65
     66    USE biometeorology_pet_mod
     67
     68    USE biometeorology_pt_mod,                                                 &
     69       ONLY: calculate_pt_static
     70
     71    USE biometeorology_utci_mod
     72
     73    USE control_parameters,                                                    &
     74       ONLY:  average_count_3d, biometeorology, dz, dz_stretch_factor,         &
     75              dz_stretch_level, humidity, initializing_actions, nz_do3d,       &
     76              simulated_time, surface_pressure
     77
     78    USE grid_variables,                                                        &
     79       ONLY:  ddx, dx, ddy, dy
     80
     81    USE indices,                                                               &
     82       ONLY:  nxl, nxr, nys, nyn, nzb, nzt, nys, nyn, nxl, nxr
     83
     84    USE kinds  !< Set precision of INTEGER and REAL arrays according to PALM
     85
     86!-- Import radiation model to obtain input for mean radiant temperature
     87    USE radiation_model_mod,                                                   &
     88       ONLY: ix, iy, iz, id, mrt_nlevels, mrt_include_sw,                      &
     89             mrtinsw, mrtinlw, mrtbl, nmrtbl, radiation, rad_sw_in,            &
     90             rad_sw_out, rad_lw_in, rad_lw_out
     91
     92    USE surface_mod,                                                           &
     93       ONLY: get_topography_top_index_ji
     94
     95    IMPLICIT NONE
    9096
    9197    PRIVATE
    9298
    93 !
    94 !-- Public functions
    95     PUBLIC biometeorology_3d_data_averaging, biometeorology_check_data_output, &
    96            biometeorology_data_output_3d, biometeorology_define_netcdf_grid
    97 
    98 !
    99 !-- Public variables and constants / NEEDS SORTING
    100 !   PUBLIC
    101 
    102 
     99!-- Declare all global variables within the module (alphabetical order)
     100    REAL(wp), DIMENSION(:,:), ALLOCATABLE ::  tmrt_grid  !< tmrt results (°C)
     101    REAL(wp), DIMENSION(:,:), ALLOCATABLE ::  pt_grid    !< PT results   (°C)
     102    REAL(wp), DIMENSION(:,:), ALLOCATABLE ::  utci_grid  !< UTCI results (°C)
     103    REAL(wp), DIMENSION(:,:), ALLOCATABLE ::  pet_grid   !< PET results  (°C)
     104!-- Grids for averaged thermal indices       
     105    REAL(wp), DIMENSION(:,:), ALLOCATABLE ::  pt_av_grid    !< PT results (aver. input)   (°C)
     106    REAL(wp), DIMENSION(:,:), ALLOCATABLE ::  utci_av_grid  !< UTCI results (aver. input) (°C)
     107    REAL(wp), DIMENSION(:,:), ALLOCATABLE ::  pet_av_grid   !< PET results (aver. input)  (°C)
     108!-- Grids for radiation_model
     109    REAL(wp), DIMENSION(:), ALLOCATABLE ::  biom_mrt     !< biometeorology mean radiant temperature for each MRT box
     110    REAL(wp), DIMENSION(:), ALLOCATABLE ::  biom_mrt_av  !< time average mean
     111
     112    INTEGER( iwp ) ::  biom_cell_level     !< cell level biom calculates for
     113    REAL ( wp )    ::  biom_output_height  !< height output is calculated in m
     114    REAL ( wp )    ::  time_biom_results   !< the time, the last set of biom results have been calculated for
     115    REAL ( wp ), PARAMETER ::  cels_offs = 273.15_wp  !< Kelvin-Celsius offset (K)
     116    REAL ( wp ), PARAMETER ::  sigma_sb  = 5.67037321E-8_wp  !< Stefan-Boltzmann constant
     117    REAL ( wp ), PARAMETER ::  human_absorb = 0.7_wp  !< SW absorbtivity of a human body (Fanger 1972)
     118    REAL ( wp ), PARAMETER ::  human_emiss = 0.97_wp  !< LW emissivity of a human body after (Fanger 1972)
     119!--
     120
     121    LOGICAL ::  aver_pt  = .FALSE.  !< switch: do pt averaging in this module? (if .FALSE. this is done globally)
     122    LOGICAL ::  aver_q   = .FALSE.  !< switch: do e  averaging in this module?
     123    LOGICAL ::  aver_u   = .FALSE.  !< switch: do u  averaging in this module?
     124    LOGICAL ::  aver_v   = .FALSE.  !< switch: do v  averaging in this module?
     125    LOGICAL ::  aver_w   = .FALSE.  !< switch: do w  averaging in this module?
     126    LOGICAL ::  average_trigger_pt   = .FALSE.  !< update averaged input on call to biom_pt?
     127    LOGICAL ::  average_trigger_utci = .FALSE.  !< update averaged input on call to biom_utci?
     128    LOGICAL ::  average_trigger_pet  = .FALSE.  !< update averaged input on call to biom_pet?
     129
     130    LOGICAL ::  biom_pt        = .TRUE.   !< Turn index PT (instant. input) on or off
     131    LOGICAL ::  biom_pt_av     = .TRUE.   !< Turn index PT (averaged input) on or off
     132    LOGICAL ::  biom_pet       = .TRUE.   !< Turn index PET (instant. input) on or off
     133    LOGICAL ::  biom_pet_av    = .TRUE.   !< Turn index PET (averaged input) on or off
     134    LOGICAL ::  biom_utci      = .TRUE.   !< Turn index UTCI (instant. input) on or off
     135    LOGICAL ::  biom_utci_av   = .TRUE.   !< Turn index UTCI (averaged input) on or off
     136
     137!-- Add INTERFACES that must be available to other modules (alphabetical order)
     138
     139    PUBLIC biom_3d_data_averaging, biom_check_data_output,                     &
     140    biom_calculate_static_grid, biom_calc_ipt,                                 &
     141    biom_check_parameters, biom_data_output_3d, biom_data_output_2d,           &
     142    biom_define_netcdf_grid, biom_determine_input_at, biom_header, biom_init,  &
     143    biom_init_arrays, biom_parin, biom_pt, biom_pt_av, biom_pet, biom_pet_av,  &
     144    biom_utci, biom_utci_av, time_biom_results
     145
     146!
     147!-- PALM interfaces:
     148!
     149!-- 3D averaging for HTCM _INPUT_ variables
     150    INTERFACE biom_3d_data_averaging
     151       MODULE PROCEDURE biom_3d_data_averaging
     152    END INTERFACE biom_3d_data_averaging
     153
     154!-- Calculate static thermal indices PT, UTCI and/or PET
     155    INTERFACE biom_calculate_static_grid
     156       MODULE PROCEDURE biom_calculate_static_grid
     157    END INTERFACE biom_calculate_static_grid
     158
     159!-- Calculate the dynamic index iPT (to be caled by the agent model)
     160    INTERFACE biom_calc_ipt
     161       MODULE PROCEDURE biom_calc_ipt
     162    END INTERFACE biom_calc_ipt
     163
     164!-- Data output checks for 2D/3D data to be done in check_parameters
     165     INTERFACE biom_check_data_output
     166        MODULE PROCEDURE biom_check_data_output
     167     END INTERFACE biom_check_data_output
     168
     169!-- Input parameter checks to be done in check_parameters
     170    INTERFACE biom_check_parameters
     171       MODULE PROCEDURE biom_check_parameters
     172    END INTERFACE biom_check_parameters
     173
     174!-- Data output of 2D quantities
     175    INTERFACE biom_data_output_2d
     176       MODULE PROCEDURE biom_data_output_2d
     177    END INTERFACE biom_data_output_2d
     178
     179!-- no 3D data, thus, no averaging of 3D data, removed
     180    INTERFACE biom_data_output_3d
     181       MODULE PROCEDURE biom_data_output_3d
     182    END INTERFACE biom_data_output_3d
     183
     184!-- Definition of data output quantities
     185    INTERFACE biom_define_netcdf_grid
     186       MODULE PROCEDURE biom_define_netcdf_grid
     187    END INTERFACE biom_define_netcdf_grid
     188
     189!-- Obtains all relevant input values to estimate local thermal comfort/stress
     190    INTERFACE biom_determine_input_at
     191       MODULE PROCEDURE biom_determine_input_at
     192    END INTERFACE biom_determine_input_at
     193
     194!-- Output of information to the header file
     195    INTERFACE biom_header
     196       MODULE PROCEDURE biom_header
     197    END INTERFACE biom_header
     198
     199!-- Initialization actions
     200    INTERFACE biom_init
     201       MODULE PROCEDURE biom_init
     202    END INTERFACE biom_init
     203
     204!-- Initialization of arrays
     205    INTERFACE biom_init_arrays
     206       MODULE PROCEDURE biom_init_arrays
     207    END INTERFACE biom_init_arrays
     208
     209!-- Reading of NAMELIST parameters
     210    INTERFACE biom_parin
     211       MODULE PROCEDURE biom_parin
     212    END INTERFACE biom_parin
     213
     214 
    103215 CONTAINS
    104 
    105 
    106 
    107 
    108 !------------------------------------------------------------------------------!
    109 ! Description:
    110 ! ------------
    111 !> Check data output for biometeorology model
    112 !------------------------------------------------------------------------------!
    113     SUBROUTINE biometeorology_check_data_output( var, unit, i, ilen, k )
    114216 
    115217 
     218!------------------------------------------------------------------------------!
     219! Description:
     220! ------------
     221!> Sum up and time-average biom input quantities as well as allocate
     222!> the array necessary for storing the average.
     223!------------------------------------------------------------------------------!
     224 SUBROUTINE biom_3d_data_averaging( mode, variable )
     225
     226    IMPLICIT NONE
     227
     228    CHARACTER (LEN=*) ::  mode     !<
     229    CHARACTER (LEN=*) ::  variable !<
     230
     231    INTEGER(iwp) ::  i  !<
     232    INTEGER(iwp) ::  j  !<
     233    INTEGER(iwp) ::  k  !<
     234
     235
     236    IF ( mode == 'allocate' )  THEN
     237
     238       SELECT CASE ( TRIM( variable ) )
     239
     240          CASE ( 'biom_mrt' )
     241                IF ( .NOT. ALLOCATED( biom_mrt_av ) )  THEN
     242                   ALLOCATE( biom_mrt_av(nmrtbl) )
     243                ENDIF
     244                biom_mrt_av = 0.0_wp
     245
     246          CASE ( 'biom_pt', 'biom_utci', 'biom_pet' )
     247
     248!--          Indices in unknown order as depending on input file, determine
     249!            first index to average und update only once
     250             IF ( .NOT. average_trigger_pt .AND. .NOT. average_trigger_utci    &
     251                .AND. .NOT. average_trigger_pet ) THEN
     252                IF ( TRIM( variable ) == 'biom_pt' ) THEN
     253                    average_trigger_pt = .TRUE.
     254                ENDIF
     255                IF ( TRIM( variable ) == 'biom_utci' ) THEN
     256                    average_trigger_utci = .TRUE.
     257                ENDIF
     258                IF ( TRIM( variable ) == 'biom_pet' ) THEN
     259                    average_trigger_pet = .TRUE.
     260                ENDIF
     261             ENDIF
     262
     263!--          Only continue if updateindex
     264             IF ( average_trigger_pt   .AND. TRIM( variable ) /= 'biom_pt')   RETURN
     265             IF ( average_trigger_utci .AND. TRIM( variable ) /= 'biom_utci') RETURN
     266             IF ( average_trigger_pet  .AND. TRIM( variable ) /= 'biom_pet')  RETURN
     267
     268!--          Set averaging switch to .TRUE. if not done by other module before!
     269             IF ( .NOT. ALLOCATED( pt_av ) )  THEN
     270                ALLOCATE( pt_av(nzb:nzt+1,nys:nyn,nxl:nxr) )
     271                aver_pt = .TRUE.
     272             ENDIF
     273             pt_av = 0.0_wp
     274
     275             IF ( .NOT. ALLOCATED( q_av ) )  THEN
     276                ALLOCATE( q_av(nzb:nzt+1,nys:nyn,nxl:nxr) )
     277                aver_q = .TRUE.
     278             ENDIF
     279             q_av = 0.0_wp
     280
     281             IF ( .NOT. ALLOCATED( u_av ) )  THEN
     282                ALLOCATE( u_av(nzb:nzt+1,nys:nyn,nxl:nxr) )
     283                aver_u = .TRUE.
     284             ENDIF
     285             u_av = 0.0_wp
     286
     287             IF ( .NOT. ALLOCATED( v_av ) )  THEN
     288                ALLOCATE( v_av(nzb:nzt+1,nys:nyn,nxl:nxr) )
     289                aver_v = .TRUE.
     290             ENDIF
     291             v_av = 0.0_wp
     292
     293             IF ( .NOT. ALLOCATED( w_av ) )  THEN
     294                ALLOCATE( w_av(nzb:nzt+1,nys:nyn,nxl:nxr) )
     295                aver_w = .TRUE.
     296             ENDIF
     297             w_av = 0.0_wp
     298
     299          CASE DEFAULT
     300             CONTINUE
     301
     302       END SELECT
     303
     304    ELSEIF ( mode == 'sum' )  THEN
     305
     306       SELECT CASE ( TRIM( variable ) )
     307
     308          CASE ( 'biom_mrt' )
     309             IF ( ALLOCATED( biom_mrt_av ) )  THEN
     310
     311                 IF ( nmrtbl > 0 )  THEN
     312                    IF ( mrt_include_sw )  THEN
     313                       biom_mrt_av(:) = biom_mrt_av(:) + &
     314                          ((human_absorb*mrtinsw(:) + human_emiss*mrtinlw(:))  &
     315                          / (human_emiss*sigma_sb)) ** .25_wp - cels_offs
     316                    ELSE
     317                       biom_mrt_av(:) = biom_mrt_av(:) + &
     318                          (human_emiss * mrtinlw(:) / sigma_sb) ** .25_wp      &
     319                          - cels_offs
     320                    ENDIF
     321                 ENDIF
     322             ENDIF
     323
     324          CASE ( 'biom_pt', 'biom_utci', 'biom_pet' )
     325
     326!--          Only continue if updateindex
     327             IF ( average_trigger_pt   .AND. TRIM( variable ) /= 'biom_pt')    &
     328                RETURN
     329             IF ( average_trigger_utci .AND. TRIM( variable ) /= 'biom_utci')  &
     330                RETURN
     331             IF ( average_trigger_pet  .AND. TRIM( variable ) /= 'biom_pet')   &
     332                RETURN
     333
     334             IF ( ALLOCATED( pt_av ) .AND. aver_pt ) THEN
     335                DO  i = nxl, nxr
     336                   DO  j = nys, nyn
     337                      DO  k = nzb, nzt+1
     338                         pt_av(k,j,i) = pt_av(k,j,i) + pt(k,j,i)
     339                      ENDDO
     340                   ENDDO
     341                ENDDO
     342             ENDIF
     343
     344             IF ( ALLOCATED( q_av )  .AND. aver_q ) THEN
     345                DO  i = nxl, nxr
     346                   DO  j = nys, nyn
     347                      DO  k = nzb, nzt+1
     348                         q_av(k,j,i) = q_av(k,j,i) + q(k,j,i)
     349                      ENDDO
     350                   ENDDO
     351                ENDDO
     352             ENDIF
     353
     354             IF ( ALLOCATED( u_av )  .AND. aver_u ) THEN
     355                DO  i = nxl, nxr
     356                   DO  j = nys, nyn
     357                      DO  k = nzb, nzt+1
     358                         u_av(k,j,i) = u_av(k,j,i) + u(k,j,i)
     359                      ENDDO
     360                   ENDDO
     361                ENDDO
     362             ENDIF
     363
     364             IF ( ALLOCATED( v_av )  .AND. aver_v ) THEN
     365                DO  i = nxl, nxr
     366                   DO  j = nys, nyn
     367                      DO  k = nzb, nzt+1
     368                         v_av(k,j,i) = v_av(k,j,i) + v(k,j,i)
     369                      ENDDO
     370                   ENDDO
     371                ENDDO
     372             ENDIF
     373
     374             IF ( ALLOCATED( w_av )  .AND. aver_w ) THEN
     375                DO  i = nxl, nxr
     376                   DO  j = nys, nyn
     377                      DO  k = nzb, nzt+1
     378                         w_av(k,j,i) = w_av(k,j,i) + w(k,j,i)
     379                      ENDDO
     380                   ENDDO
     381                ENDDO
     382             ENDIF
     383
     384           CASE DEFAULT
     385             CONTINUE
     386
     387       END SELECT
     388
     389    ELSEIF ( mode == 'average' )  THEN
     390
     391       SELECT CASE ( TRIM( variable ) )
     392
     393          CASE ( 'biom_mrt' )
     394             IF ( ALLOCATED( biom_mrt_av ) )  THEN
     395                biom_mrt_av(:) = biom_mrt_av(:) / REAL( average_count_3d, KIND=wp )
     396             ENDIF
     397
     398          CASE ( 'biom_pt', 'biom_utci', 'biom_pet' )
     399
     400!--          Only continue if update index
     401             IF ( average_trigger_pt   .AND. TRIM( variable ) /= 'biom_pt')    &
     402                RETURN
     403             IF ( average_trigger_utci .AND. TRIM( variable ) /= 'biom_utci')  &
     404                RETURN
     405             IF ( average_trigger_pet  .AND. TRIM( variable ) /= 'biom_pet')   &
     406                RETURN
     407
     408             IF ( ALLOCATED( pt_av ) .AND. aver_pt ) THEN
     409                DO  i = nxl, nxr
     410                   DO  j = nys, nyn
     411                      DO  k = nzb, nzt+1
     412                         pt_av(k,j,i) = pt_av(k,j,i) / REAL( average_count_3d, KIND=wp )
     413                      ENDDO
     414                   ENDDO
     415                ENDDO
     416             ENDIF
     417
     418             IF ( ALLOCATED( q_av ) .AND. aver_q ) THEN
     419                DO  i = nxl, nxr
     420                   DO  j = nys, nyn
     421                      DO  k = nzb, nzt+1
     422                         q_av(k,j,i) = q_av(k,j,i) / REAL( average_count_3d, KIND=wp )
     423                      ENDDO
     424                   ENDDO
     425                ENDDO
     426             ENDIF
     427
     428             IF ( ALLOCATED( u_av ) .AND. aver_u ) THEN
     429                DO  i = nxl, nxr
     430                   DO  j = nys, nyn
     431                      DO  k = nzb, nzt+1
     432                         u_av(k,j,i) = u_av(k,j,i) / REAL( average_count_3d, KIND=wp )
     433                      ENDDO
     434                   ENDDO
     435                ENDDO
     436             ENDIF
     437
     438             IF ( ALLOCATED( v_av ) .AND. aver_v ) THEN
     439                DO  i = nxl, nxr
     440                   DO  j = nys, nyn
     441                      DO  k = nzb, nzt+1
     442                         v_av(k,j,i) = v_av(k,j,i) / REAL( average_count_3d, KIND=wp )
     443                      ENDDO
     444                   ENDDO
     445                ENDDO
     446             ENDIF
     447
     448             IF ( ALLOCATED( w_av ) .AND. aver_w ) THEN
     449                DO  i = nxl, nxr
     450                   DO  j = nys, nyn
     451                      DO  k = nzb, nzt+1
     452                         w_av(k,j,i) = w_av(k,j,i) / REAL( average_count_3d, KIND=wp )
     453                      ENDDO
     454                   ENDDO
     455                ENDDO
     456             ENDIF
     457
     458!--          Udate thermal indices with derived averages
     459             CALL biom_calculate_static_grid ( .TRUE. )
     460
     461        END SELECT
     462
     463    ENDIF
     464
     465
     466 END SUBROUTINE biom_3d_data_averaging
     467
     468
     469
     470!------------------------------------------------------------------------------!
     471! Description:
     472! ------------
     473!> Check data output for biometeorology model
     474!------------------------------------------------------------------------------!
     475 SUBROUTINE biom_check_data_output( var, unit )
     476
    116477       USE control_parameters,                                                 &
    117478           ONLY: data_output, message_string
    … …  
    119480       IMPLICIT NONE
    120481
    121        CHARACTER (LEN=*) ::  unit     !<
    122        CHARACTER (LEN=*) ::  var      !<
    123 
    124        INTEGER(iwp) :: i
    125        INTEGER(iwp) :: ilen
    126        INTEGER(iwp) :: k
    127 
    128        SELECT CASE ( TRIM( var ) )
    129 
    130           CASE ( 'bio_mrt', 'bio_pet'  )
    131              IF ( .NOT.  radiation ) THEN
    132                 message_string = 'output of "' // TRIM( var ) // '" require'&
    133                                  // 's radiation = .TRUE.'
    134                 CALL message( 'check_parameters', 'PA0509', 1, 2, 0, 6, 0 )
    135              ENDIF
    136              IF ( mrt_nlevels == 0 ) THEN
    137                 message_string = 'output of "' // TRIM( var ) // '" require'&
    138                                  // 's mrt_nlevels > 0'
    139                 CALL message( 'check_parameters', 'PA0510', 1, 2, 0, 6, 0 )
    140              ENDIF
    141              IF ( TRIM( var ) == 'bio_mrt' ) THEN
    142                 unit = 'K'
    143              ELSE
    144                 unit = 'W m-2'
    145              ENDIF
    146 
    147           CASE DEFAULT
    148              unit = 'illegal'
    149 
    150        END SELECT
    151 
    152     END SUBROUTINE biometeorology_check_data_output
    153 
    154 
    155 
    156 !------------------------------------------------------------------------------!
    157 !
    158 ! Description:
    159 ! ------------
    160 !> Subroutine for averaging 3D data
    161 !------------------------------------------------------------------------------!
    162 SUBROUTINE biometeorology_3d_data_averaging( mode, variable )
    163 
    164     USE control_parameters
     482       CHARACTER (LEN=*) ::  unit     !< The unit for the variable var
     483       CHARACTER (LEN=*) ::  var      !< The variable in question
     484
     485
     486                unit = '°C'
     487                IF ( .NOT. biometeorology ) THEN
     488                  message_string = 'output of "' // TRIM( var ) // '" req' //  &
     489                        'uires biometeorology = .TRUE. in initialisati' &
     490                        // 'on_parameters'
     491                  CALL message( 'check_parameters', 'PA0114', 1, 2, 0, 6, 0 )
     492                  unit = 'illegal'
     493                ENDIF
     494                IF ( .NOT.  radiation ) THEN
     495                   message_string = 'output of "' // TRIM( var ) // '" require'&
     496                                    // 's radiation = .TRUE.'
     497                   CALL message( 'check_parameters', 'PA0509', 1, 2, 0, 6, 0 )
     498                   unit = 'illegal'
     499                ENDIF
     500                IF ( mrt_nlevels == 0 ) THEN
     501                   message_string = 'output of "' // TRIM( var ) // '" require'&
     502                                    // 's mrt_nlevels > 0'
     503                   CALL message( 'check_parameters', 'PA0510', 1, 2, 0, 6, 0 )
     504                   unit = 'illegal'
     505                ENDIF
     506
     507 END SUBROUTINE biom_check_data_output
     508
     509!------------------------------------------------------------------------------!
     510! Description:
     511! ------------
     512!> Check parameters routine for biom module
     513!> check_parameters 1302
     514!------------------------------------------------------------------------------!
     515 SUBROUTINE biom_check_parameters
     516
     517    USE control_parameters,                                                    &
     518        ONLY:  message_string
     519
     520
     521    IMPLICIT NONE
     522
     523
     524!--    Disabled as long as radiation model not available
     525       IF ( .NOT. radiation )  THEN
     526          message_string = 'The human thermal comfort module does require ' // &
     527                           'radiation information in terms of the mean ' //    &
     528                           'radiant temperature, but radiation is not enabled!'
     529          CALL message( 'check_parameters', 'PAHU01', 0, 0, 0, 6, 0 )
     530       ENDIF
     531
     532       IF ( .NOT. humidity )  THEN
     533          message_string = 'The human thermal comfort module does require ' // &
     534                           'air humidity information, but humidity module ' // &
     535                           'is disabled!'
     536          CALL message( 'check_parameters', 'PAHU02', 0, 0, 0, 6, 0 )
     537       ENDIF
     538
     539
     540 END SUBROUTINE biom_check_parameters
     541
     542
     543!------------------------------------------------------------------------------!
     544!
     545! Description:
     546! ------------
     547!> Subroutine defining 3D output variables (dummy, always 2d!)
     548!> data_output_3d 709ff
     549!------------------------------------------------------------------------------!
     550 SUBROUTINE biom_data_output_3d( av, variable, found, local_pf, nzb_do, nzt_do )
    165551
    166552    USE indices
    … …  
    168554    USE kinds
    169555
    170     IMPLICIT NONE
    171 
    172     CHARACTER (LEN=*) ::  mode    !<
    173     CHARACTER (LEN=*) :: variable !<
    174 
    175     INTEGER(iwp) ::  i !<
    176     INTEGER(iwp) ::  j !<
    177     INTEGER(iwp) ::  k !<
    178     INTEGER(iwp) ::  l, m !< index of current surface element
    179 
    180     REAL(wp)     ::  mrt, pet
    181 
    182     IF ( mode == 'allocate' )  THEN
    183 
    184        SELECT CASE ( TRIM( variable ) )
    185              CASE ( 'bio_mrt' )
    186                 IF ( .NOT. ALLOCATED( bio_mrt_av ) )  THEN
    187                    ALLOCATE( bio_mrt_av(nmrtbl) )
    188                 ENDIF
    189                 bio_mrt_av = 0.0_wp
    190 
    191              CASE ( 'bio_pet' )
    192                 IF ( .NOT. ALLOCATED( bio_pet_av ) )  THEN
    193                    ALLOCATE( bio_pet_av(nmrtbl) )
    194                 ENDIF
    195                 bio_pet_av = 0.0_wp
    196 
    197           CASE DEFAULT
    198              CONTINUE
    199 
    200        END SELECT
    201 
    202     ELSEIF ( mode == 'sum' )  THEN
    203 
    204        SELECT CASE ( TRIM( variable ) )
    205 
    206           CASE ( 'bio_mrt' )
    207              IF ( ALLOCATED( bio_mrt_av ) )  THEN
    208 
    209                  IF ( nmrtbl > 0 )  THEN
    210                     IF ( mrt_include_sw )  THEN
    211                        bio_mrt_av(:) = bio_mrt_av(:) + ((human_absorb*mrtinsw(:) &
    212                            + human_emiss*mrtinlw(:)) / (human_emiss*sigma_sb)) ** .25_wp
    213                     ELSE
    214                        bio_mrt_av(:) = bio_mrt_av(:) + &
    215                                        (human_emiss * mrtinlw(:) / sigma_sb) ** .25_wp
    216                     ENDIF
    217                  ENDIF
    218              ENDIF
    219 
    220           CASE ( 'bio_pet' )
    221              IF ( ALLOCATED( bio_pet_av ) )  THEN
    222                 DO  l = 1, nmrtbl
    223                    i = mrtbl(ix,l)
    224                    j = mrtbl(iy,l)
    225                    k = mrtbl(iz,l)
    226 
    227                    IF ( mrt_include_sw )  THEN
    228                       mrt = ((human_absorb*mrtinsw(l) + human_emiss*mrtinlw(l)) &
    229                              / (human_emiss*sigma_sb)) ** .25_wp
    230                    ELSE
    231                       mrt = mrt + (human_emiss * mrtinlw(l) / sigma_sb) ** .25_wp
    232                    ENDIF
    233 
    234                    CALL calculate_pet_static(                                      &
    235                            pt(k,j,i) * (hyp(k) / 100000.0_wp )**0.286_wp + t_zero, &  !< Air temperature (°C)
    236                            q(k,j,i) * hyp(k) / ( q(k,j,i) + rd_d_rv ) / 100._wp,   &  !< Vapor pressure (hPa)
    237                            SQRT( MAX( ( ( u(k,j,i) + u(k,j,i+1) ) * 0.5_wp )**2 +  &
    238                                       ( ( v(k,j,i) + v(k,j+1,i) ) * 0.5_wp )**2 +  &
    239                                       ( ( w(k,j,i) + w(k-1,j,i) ) * 0.5_wp )**2,   &
    240                                       0.01_wp ) ),                                 &  !< Wind speed (at scalar gridpoint) (m/s)
    241                            mrt + t_zero,                                           &  !< Mean radiant temperature (°C)
    242                            (hyp(k) + p(k,j,i)) * 0.01_wp,                          &  !< Air pressure (hPa)
    243                            pet )                                                      !< output variable of PET
    244 
    245                     bio_pet_av(l) = bio_pet_av(l) + pet
    246                 ENDDO
    247              ENDIF
    248 
    249           CASE DEFAULT
    250              CONTINUE
    251 
    252        END SELECT
    253 
    254     ELSEIF ( mode == 'average' )  THEN
    255 
    256        SELECT CASE ( TRIM( variable ) )
    257 
    258           CASE ( 'bio_mrt' )
    259              IF ( ALLOCATED( bio_mrt_av ) )  THEN
    260                 bio_mrt_av(:) = bio_mrt_av(:) / REAL( average_count_3d, KIND=wp )
    261              ENDIF
    262 
    263           CASE ( 'bio_pet' )
    264              IF ( ALLOCATED( bio_pet_av ) )  THEN
    265                 bio_pet_av(:) = bio_pet_av(:) / REAL( average_count_3d, KIND=wp )
    266              ENDIF
    267 
    268        END SELECT
    269 
    270     ENDIF
    271 
    272 END SUBROUTINE biometeorology_3d_data_averaging
    273 
    274 
    275 !------------------------------------------------------------------------------!
    276 !
    277 ! Description:
    278 ! ------------
    279 !> Subroutine defining appropriate grid for netcdf variables.
    280 !> It is called out from subroutine netcdf.
    281 !------------------------------------------------------------------------------!
    282  SUBROUTINE biometeorology_define_netcdf_grid( var, found, grid_x, grid_y, grid_z )
    283 
    284     IMPLICIT NONE
    285 
    286     CHARACTER (LEN=*), INTENT(IN)  ::  var         !<
    287     LOGICAL, INTENT(OUT)           ::  found       !<
    288     CHARACTER (LEN=*), INTENT(OUT) ::  grid_x      !<
    289     CHARACTER (LEN=*), INTENT(OUT) ::  grid_y      !<
    290     CHARACTER (LEN=*), INTENT(OUT) ::  grid_z      !<
    291 
    292     found  = .TRUE.
    293 
    294 
    295 !
    296 !-- Check for the grid
    297     SELECT CASE ( TRIM( var ) )
    298 
    299        CASE ( 'bio_mrt', 'bio_pet' )
    300           grid_x = 'x'
    301           grid_y = 'y'
    302           grid_z = 'zu'
    303 
    304        CASE DEFAULT
    305           found  = .FALSE.
    306           grid_x = 'none'
    307           grid_y = 'none'
    308           grid_z = 'none'
    309 
    310         END SELECT
    311 
    312  END SUBROUTINE biometeorology_define_netcdf_grid
    313 
    314 
    315 !------------------------------------------------------------------------------!
    316 !
    317 ! Description:
    318 ! ------------
    319 !> Subroutine defining 3D output variables
    320 !------------------------------------------------------------------------------!
    321  SUBROUTINE biometeorology_data_output_3d( av, variable, found, local_pf, nzb_do, nzt_do )
    322  
    323 
    324     USE indices
    325 
    326     USE kinds
    327 
    328 
    329     IMPLICIT NONE
    330 
    331     CHARACTER (LEN=*) ::  variable !<
    332 
    333     INTEGER(iwp) ::  av          !<
    334     INTEGER(iwp) ::  i, j, k, l  !<
    335     INTEGER(iwp) ::  nzb_do      !<
    336     INTEGER(iwp) ::  nzt_do      !<
    337 
    338     LOGICAL      ::  found       !<
    339 
    340     REAL(wp)     ::  fill_value = -999.0_wp    !< value for the _FillValue attribute
    341 
    342     REAL(sp), DIMENSION(nxl:nxr,nys:nyn,nzb_do:nzt_do) ::  local_pf !<
    343    
    344     REAL(wp)     ::  mrt, pet
     556
     557    IMPLICIT NONE
     558
     559!-- Input variables
     560    CHARACTER (LEN=*), INTENT(IN) ::  variable   !< Char identifier to select var for output
     561    INTEGER(iwp), INTENT(IN) ::  av       !< Use averaged data? 0 = no, 1 = yes?
     562    INTEGER(iwp), INTENT(IN) ::  nzb_do   !< Unused. 2D. nz bottom to nz top
     563    INTEGER(iwp), INTENT(IN) ::  nzt_do   !< Unused.
     564
     565!-- Output variables
     566    LOGICAL, INTENT(OUT) ::  found   !< Output found?
     567    REAL(sp), DIMENSION(nxl:nxr,nys:nyn,nzb_do:nzt_do) ::  local_pf   !< Temp. result grid to return
     568
     569!-- Internal variables
     570    INTEGER(iwp) ::  l    !< Running index, radiation grid
     571    INTEGER(iwp) ::  i    !< Running index, x-dir
     572    INTEGER(iwp) ::  j    !< Running index, y-dir
     573    INTEGER(iwp) ::  k    !< Running index, z-dir
     574
     575    CHARACTER (LEN=:), allocatable ::  variable_short  !< Trimmed version of char variable
     576
     577    REAL(wp), PARAMETER ::  fill_value = -999._wp
     578    REAL(wp) ::  mrt  !< Buffer for mean radiant temperature
    345579
    346580    found = .TRUE.
    347 
    348 
    349     SELECT CASE ( TRIM( variable ) )
    350 
    351 
    352         CASE ( 'bio_mrt' )
     581    variable_short = TRIM( variable )
     582
     583    IF ( variable_short(1:4) /= 'biom' ) THEN
     584!--   Nothing to do, set found to FALSE and return immediatelly
     585      found = .FALSE.
     586      RETURN
     587    ENDIF
     588
     589    SELECT CASE ( variable_short )
     590
     591        CASE ( 'biom_mrt' )
    353592
    354593            local_pf = REAL( fill_value, KIND = wp )
    … …  
    366605            ENDDO
    367606
    368         CASE ( 'bio_pet' )
    369             local_pf = REAL( fill_value, KIND = wp )
     607        CASE ( 'biom_tmrt' )        ! 2d-array
     608              DO  i = nxl, nxr
     609                 DO  j = nys, nyn
     610                    local_pf(i,j,nzb_do) = REAL( tmrt_grid(j,i), sp )
     611                 ENDDO
     612              ENDDO
     613
     614        CASE ( 'biom_pt' )        ! 2d-array
    370615            IF ( av == 0 )  THEN
    371                 DO  l = 1, nmrtbl
    372                     i = mrtbl(ix,l)
    373                     j = mrtbl(iy,l)
    374                     k = mrtbl(iz,l)
    375 
    376                     IF ( mrt_include_sw )  THEN
    377                         mrt = ((human_absorb*mrtinsw(l) + human_emiss*mrtinlw(l)) &
    378                                / (human_emiss*sigma_sb)) ** .25_wp
    379                     ELSE
    380                         mrt = mrt + (human_emiss*mrtinlw(l) / sigma_sb) ** .25_wp
    381                     ENDIF
    382 
    383                     CALL calculate_pet_static(                                      &
    384                        pt(k,j,i) * (hyp(k) / 100000.0_wp )**0.286_wp + t_zero, &  !< Air temperature (°C)
    385                        q(k,j,i) * hyp(k) / ( q(k,j,i) + rd_d_rv ) / 100.0_wp,  &  !< Vapor pressure (hPa)
    386                        SQRT( MAX( ( ( u(k,j,i) + u(k,j,i+1) ) * 0.5_wp )**2 +  &
    387                                   ( ( v(k,j,i) + v(k,j+1,i) ) * 0.5_wp )**2 +  &
    388                                   ( ( w(k,j,i) + w(k-1,j,i) ) * 0.5_wp )**2,   &
    389                                   0.01_wp ) ),                                 &  !< Wind speed (at scalar gridpoint) (m/s)
    390                        mrt + t_zero,                                           &  !< Mean radiant temperature (°C)
    391                        (hyp(k) + p(k,j,i)) * 0.01_wp,                          &  !< Air pressure (hPa)
    392                        pet )                                                      !< output variable of PET
    393 
    394                local_pf(i,j,k) = pet
    395 
    396             ENDDO
    397          ELSE
    398             IF ( ALLOCATED( bio_pet_av ) ) THEN
    399                DO  l = 1, nmrtbl
    400                   i = mrtbl(ix,l)
    401                   j = mrtbl(iy,l)
    402                   k = mrtbl(iz,l)
    403                   local_pf(i,j,k) = bio_pet_av(l)
    404                ENDDO
    405             ENDIF
    406          ENDIF
     616              DO  i = nxl, nxr
     617                 DO  j = nys, nyn
     618                    local_pf(i,j,nzb_do) = REAL( pt_grid(j,i), sp )
     619                 ENDDO
     620              ENDDO
     621            ELSE
     622              DO  i = nxl, nxr
     623                 DO  j = nys, nyn
     624                    local_pf(i,j,nzb_do) = REAL( pt_av_grid(j,i), sp )
     625                 ENDDO
     626              ENDDO
     627            END IF
     628
     629        CASE ( 'biom_utci' )        ! 2d-array
     630            IF ( av == 0 )  THEN
     631              DO  i = nxl, nxr
     632                 DO  j = nys, nyn
     633                    local_pf(i,j,nzb_do) = REAL( utci_grid(j,i), sp )
     634                 ENDDO
     635              ENDDO
     636            ELSE
     637              DO  i = nxl, nxr
     638                 DO  j = nys, nyn
     639                    local_pf(i,j,nzb_do) = REAL( utci_av_grid(j,i), sp )
     640                 ENDDO
     641              ENDDO
     642            END IF
     643
     644        CASE ( 'biom_pet' )        ! 2d-array
     645            IF ( av == 0 )  THEN
     646              DO  i = nxl, nxr
     647                 DO  j = nys, nyn
     648                    local_pf(i,j,nzb_do) = REAL( pet_grid(j,i), sp )
     649                 ENDDO
     650              ENDDO
     651            ELSE
     652              DO  i = nxl, nxr
     653                 DO  j = nys, nyn
     654                    local_pf(i,j,nzb_do) = REAL( pet_av_grid(j,i), sp )
     655                 ENDDO
     656              ENDDO
     657            END IF
    407658
    408659       CASE DEFAULT
    … …  
    411662    END SELECT
    412663
    413 
    414  END SUBROUTINE biometeorology_data_output_3d
    415 
    416 
    417 
    418 !------------------------------------------------------------------------------!
    419 !
    420 ! Description:
    421 ! ------------
    422 !> PhysiologiCALLy Equivalent Temperature (PET),
    423 !  stationary (calculated based on MEMI),
    424 !  Subroutine based on PETBER vers. 1.5.1996 by P. Hoeppe
    425 !
    426 !  Input arguments:
    427 ! ----------------
    428 !  - ta         : Air temperature               (°C)   REAL(wp)
    429 !  - tmrt       : Mean radiant temperature      (°C)   REAL(wp)
    430 !  - v          : Wind speed                    (m/s)   REAL(wp)
    431 !  - vpa        : Vapor pressure                (hPa)   REAL(wp)
    432 !  - p          : Air pressure                  (hPa)   REAL(wp)
    433 !
    434 !  Output arguments:
    435 ! ----------------
    436 !  - tx         : PET                           (°C)   REAL(wp)
    437 !------------------------------------------------------------------------------!
    438 
    439  SUBROUTINE calculate_pet_static(                                              &
    440 !-- Meteorological input
    441   ta, vpa, v, tmrt, p,                                                         &
     664 END SUBROUTINE biom_data_output_3d
     665
     666!------------------------------------------------------------------------------!
     667!
     668! Description:
     669! ------------
     670!> Subroutine defining 2D output variables
     671!> data_output_2d 1188ff
     672!------------------------------------------------------------------------------!
     673 SUBROUTINE biom_data_output_2d( av, variable, found, grid, local_pf,          &
     674                                      two_d, nzb_do, nzt_do, fill_value )
     675
     676    USE indices,                                                               &
     677       ONLY: nxl, nxl, nxr, nxr, nyn, nyn, nys, nys, nzb, nzt
     678
     679    USE kinds
     680
     681
     682    IMPLICIT NONE
     683
     684!-- Input variables
     685    CHARACTER (LEN=*), INTENT(IN) ::  variable    !< Char identifier to select var for output
     686    INTEGER(iwp), INTENT(IN)      ::  av          !< Use averaged data? 0 = no, 1 = yes?
     687    INTEGER(iwp), INTENT(IN)      ::  nzb_do      !< Unused. 2D. nz bottom to nz top
     688    INTEGER(iwp), INTENT(IN)      ::  nzt_do      !< Unused.
     689    REAL(wp), INTENT(in)          ::  fill_value  !< Fill value for unassigned locations
     690
    442691!-- Output variables
    443   tx )  !,                                                                     &
    444 !-- Configure sample person (optional)
    445 !  age, mbody, ht, work, eta, icl, fcl, pos, sex )
    446 
    447    IMPLICIT NONE
    448 
    449    REAL(wp), INTENT( IN ) :: ta, tmrt, v, vpa, p
    450 
    451    REAL(wp), INTENT ( OUT ) :: tx
    452 
    453 !  REAL(wp), INTENT ( in ), optional :: age, mbody, ht, work, eta, icl, fcl
    454    REAL(wp) :: age, mbody, ht, work, eta, icl, fcl
    455 !  INTEGER(iwp), INTENT ( in ), optional :: pos, sex
    456    INTEGER(iwp) :: pos, sex
    457 
    458    REAL(wp) :: acl, adu, aeff, cair, cb, emcl, emsk, ere, erel, esw,           &
    459           evap, facl, feff, food, h, po, rdcl, rdsk, rob, rtv,              &
    460           sigm, vpts,                                                          &
    461    !  former intent (out)
    462    !  - tsk        : Skin temperature              (°C)    real
    463    !  - tcl        : Clothing temperature          (°C)    real
    464    !  - ws         :                                       real
    465    !  - wetsk      : Fraction of wet skin                  real
    466           tsk, tcl, wetsk
    467 
    468 
    469 !--      Person data
    470 !   IF ( .NOT. present( age ) )   age   = 35.
    471 !   IF ( .NOT. present( mbody ) ) mbody = 75.
    472 !   IF ( .NOT. present( ht ) )    ht    =  1.75
    473 !   IF ( .NOT. present( work ) )  work  = 80.
    474 !   IF ( .NOT. present( eta ) )   eta   =  0.
    475 !   IF ( .NOT. present( icl ) )   icl   =  0.9
    476 !   IF ( .NOT. present( fcl ) )   fcl   =  1.15
    477 !   IF ( .NOT. present( pos ) )   pos   =  1
    478 !   IF ( .NOT. present( sex ) )   sex   =  1
    479    
    480    age   = 35.
    481    mbody = 75.
    482    ht    =  1.75
    483    work  = 80.
    484    eta   =  0.
    485    icl   =  0.9
    486    fcl   =  1.15
    487    pos   =  1
    488    sex   =  1
    489 
    490 !-- constants
    491    po   = 1013.25  !< preassure at sea level
    492    ! p    = 1013.25   !< local preassure (hPa), now defined as input variable
    493    rob  =    1.06
    494    cb   = 3.64 * 1000.
    495    food =    0.
    496    emsk =    0.99
    497    emcl =    0.95
    498    evap = 2.42 * 10. ** 6.
    499    sigm = 5.67 * 10. **(-8.)
    500 
    501 
    502 !   write(9,*) 'Call calculate_pet_static(ta=', ta, ', vpa=', vpa, &
    503 !              ', v=', v, ', tmrt=', tmrt, ', p=', p
    504 !   flush(9)
    505 !-- call subfunctions
    506    CALL  INKOERP ( age, cair, eta, ere, erel, evap, h, ht, mbody,              &
    507          p, rtv, sex, ta, vpa, work )
    508 
    509 
    510    CALL BERECH ( acl, adu, aeff, cair, cb, emcl, emsk,                 &
    511         ere, erel, esw, evap, facl, fcl, feff, food, h, ht, icl,               &
    512         mbody, p, po, rdcl, rdsk, rob, sex, sigm,                              &
    513         ta, tcl, tmrt, tsk, v, vpa, vpts, wetsk )
    514 
    515 
    516    CALL PET ( acl, adu, aeff, cair, emcl, emsk, esw, evap,                     &
    517         facl, feff, h, p, po, rdcl, rdsk,                                      &
    518         rtv, sigm, ta, tcl, tsk, tx, vpts, wetsk )
    519 
    520 
    521  END SUBROUTINE calculate_pet_static
    522 
    523 
    524 !------------------------------------------------------------------------------!
    525 ! Description:
    526 ! ------------
    527 !> Calculate internal energy ballance
    528 !------------------------------------------------------------------------------!
    529  SUBROUTINE inkoerp( age, cair, eta, ere, erel, evap, h, ht, mbody,            &
    530    &     p, rtv, sex, ta, vpa, work )
    531 
    532 
    533    REAL(wp) :: age, cair, eta, ere, erel, eres, evap, h, ht, mbody,            &
    534    &           met, p, rtv, ta, tex, vpa, vpex, work
    535 
    536    INTEGER(iwp) :: sex
    537 
    538    IF ( sex .EQ. 1 ) THEN
    539      met = 3.45 * mbody ** ( 3. / 4. ) * (1. + 0.004 *          &
    540            ( 30. - age) + 0.010 * ( ( ht * 100. /                     &
    541            ( mbody ** ( 1. / 3. ) ) ) - 43.4 ) )
    542    ELSE IF ( sex .EQ. 2 ) THEN
    543      met = 3.19 * mbody ** ( 3. / 4. ) * ( 1. + 0.004 *         &
    544            ( 30. - age ) + 0.018 * ( ( ht * 100. / ( mbody **         &
    545            ( 1. / 3. ) ) ) - 42.1 ) )
    546    END IF
    547    met = work + met
    548 
    549    h = met * (1. - eta)
    550 
    551 !-- SENSIBLE RESPIRATION ENERGY
    552 
    553    cair = 1.01 * 1000.
    554    tex  = 0.47 * ta + 21.0
    555    rtv  = 1.44 * 10. ** (-6.) * met
    556    eres = cair * (ta - tex) * rtv
    557 
    558 !-- LATENT RESPIRATION ENERGY
    559 
    560    vpex = 6.11 * 10. ** ( 7.45 * tex / ( 235. + tex ) )
    561    erel = 0.623 * evap / p * ( vpa - vpex ) * rtv
    562 
    563 !-- SUM OF RESULTS
    564 
    565    ere = eres + erel
    566    RETURN
    567  END SUBROUTINE inkoerp
    568 
    569 
    570 !------------------------------------------------------------------------------!
    571 ! Description:
    572 ! ------------
    573 !> Calculate heat gain or loss
    574 !------------------------------------------------------------------------------!
    575  SUBROUTINE BERECH( acl, adu, aeff, cair, cb, emcl, emsk,              &
    576         ere, erel, esw, evap, facl, fcl, feff, food, h, ht, icl,               &
    577         mbody, p, po, rdcl, rdsk, rob, sex, sigm,                              &
    578         ta, tcl, tmrt, tsk, v, vpa, vpts, wetsk )
    579 
    580 
    581    REAL(wp) :: acl, adu, aeff, c(0:10), cair, cb, cbare,                       &
    582         cclo, csum, di, ed, emcl, emsk, enbal,                                 &
    583         enbal2, ere, erel, esw, eswdif, eswphy, eswpot,                        &
    584         evap, facl, fcl, fec, feff, food, h, hc, he, ht, htcl, icl,            &
    585         mbody, p, po, r1, r2, rbare, rcl,                                      &
    586         rclo, rclo2, rdcl, rdsk, rob, rsum, sigm, sw, swf, swm,                &
    587         ta, tbody, tcl, tcore(1:7), tmrt, tsk, v, vb, vb1, vb2,                &
    588         vpa, vpts, wetsk, wd, wr, ws, wsum, xx, y
    589 
    590    INTEGER(iwp) :: count1, count3, j, sex
    591    logical :: skipIncreaseCount
    592 
    593    wetsk = 0.
    594    adu = 0.203 * mbody ** 0.425 * ht ** 0.725
    595 
    596    hc = 2.67 + ( 6.5 * v ** 0.67)
    597    hc   = hc * (p /po) ** 0.55
    598    feff = 0.725  !< Posture: 0.725 for stading, 0.696 for sitting
    599    facl = (- 2.36 + 173.51 * icl - 100.76 * icl * icl + 19.28                  &
    600          * (icl ** 3.)) / 100.
    601 
    602    IF ( facl .GT. 1. )   facl = 1.
    603    rcl = ( icl / 6.45) / facl
    604    IF ( icl .GE. 2. )  y  = 1.
    605 
    606    IF ( ( icl .GT. 0.6 ) .AND. ( icl .LT. 2. ) )  y = ( ht - 0.2 ) / ht
    607    IF ( ( icl .LE. 0.6 ) .AND. ( icl .GT. 0.3 ) ) y = 0.5
    608    IF ( ( icl .LE. 0.3 ) .AND. ( icl .GT. 0. ) )  y = 0.1
    609 
    610    r2   = adu * (fcl - 1. + facl) / (2. * 3.14 * ht * y)
    611    r1   = facl * adu / (2. * 3.14 * ht * y)
    612 
    613    di   = r2 - r1
    614 
    615 !-- SKIN TEMPERATURE
    616 
    617    DO j = 1,7
    618 
    619      tsk    = 34.
    620      count1 = 0
    621      tcl    = ( ta + tmrt + tsk ) / 3.
    622      count3 = 1
    623      enbal2 = 0.
    624 
    625      DO
    626        acl   = adu * facl + adu * ( fcl - 1. )
    627        rclo2 = emcl * sigm * ( ( tcl + 273.2 )**4. -                           &
    628          ( tmrt + 273.2 )** 4. ) * feff
    629        htcl  = 6.28 * ht * y * di / ( rcl * LOG( r2 / r1 ) * acl )
    630        tsk   = 1. / htcl * ( hc * ( tcl - ta ) + rclo2 ) + tcl
    631 
    632 !--    RADIATION SALDO
    633 
    634        aeff  = adu * feff
    635        rbare = aeff * ( 1. - facl ) * emsk * sigm *                            &
    636          ( ( tmrt + 273.2 )** 4. - ( tsk + 273.2 )** 4. )
    637        rclo  = feff * acl * emcl * sigm *                                      &
    638          ( ( tmrt + 273.2 )** 4. - ( tcl + 273.2 )** 4. )
    639        rsum  = rbare + rclo
    640 
    641 !--    CONVECTION
    642 
    643        cbare = hc * ( ta - tsk ) * adu * ( 1. - facl )
    644        cclo  = hc * ( ta - tcl ) * acl
    645        csum  = cbare + cclo
    646 
    647  !--   CORE TEMPERATUR
    648 
    649        c(0)  = h + ere
    650        c(1)  = adu * rob * cb
    651        c(2)  = 18. - 0.5 * tsk
    652        c(3)  = 5.28 * adu * c(2)
    653        c(4)  = 0.0208 * c(1)
    654        c(5)  = 0.76075 * c(1)
    655        c(6)  = c(3) - c(5) - tsk * c(4)
    656        c(7)  = - c(0) * c(2) - tsk * c(3) + tsk * c(5)
    657        c(8)  = c(6) * c(6) - 4. * c(4) * c(7)
    658        c(9)  = 5.28 * adu - c(5) - c(4) * tsk
    659        c(10) = c(9) * c(9) - 4. * c(4) *                                       &
    660          ( c(5) * tsk - c(0) - 5.28 * adu * tsk )
    661 
    662        IF ( tsk .EQ. 36. ) tsk = 36.01
    663        tcore(7) = c(0) / ( 5.28 * adu + c(1) * 6.3 / 3600. ) + tsk
    664        tcore(3) = c(0) / ( 5.28 * adu + ( c(1) * 6.3 / 3600. ) /               &
    665          ( 1. + 0.5 * ( 34. - tsk ) ) ) + tsk
    666        IF ( c(10) .GE. 0.) THEN
    667          tcore(6) = ( - c(9) - c(10)**0.5 ) / ( 2. * c(4) )
    668          tcore(1) = ( - c(9) + c(10)**0.5 ) / ( 2. * c(4) )
    669        END IF
    670        ! 22
    671        IF ( c(8) .GE. 0. ) THEN
    672          tcore(2) = ( - c(6) + ABS( c(8) ) ** 0.5 ) / ( 2. * c(4) )
    673          tcore(5) = ( - c(6) - ABS( c(8) ) ** 0.5 ) / ( 2. * c(4) )
    674          tcore(4) = c(0) / ( 5.28 * adu + c(1) * 1. / 40. ) + tsk
    675        END IF
    676        ! 24
    677 
    678 !--    TRANSPIRATION
    679 
    680        tbody = 0.1 * tsk + 0.9 * tcore(j)
    681        swm   = 304.94 * ( tbody - 36.6 ) * adu / 3600000.
    682        vpts  = 6.11 * 10.**( 7.45 * tsk / ( 235. + tsk ) )
    683 
    684        IF ( tbody .LE. 36.6 ) swm = 0.
    685        ! swf = 0.7 * swm
    686 
    687        IF ( sex .EQ. 1 ) sw = swm
    688        IF ( sex .EQ. 2 ) sw = 0.7 * swm  ! swf
    689        eswphy = - sw * evap
    690        he     = 0.633 * hc / ( p * cair )
    691        fec    = 1. / ( 1. + 0.92 * hc * rcl )
    692        eswpot = he * ( vpa - vpts ) * adu * evap * fec
    693        wetsk  = eswphy / eswpot
    694 
    695        IF ( wetsk .GT. 1. ) wetsk = 1.
    696 
    697        eswdif = eswphy - eswpot
    698 
    699        IF ( eswdif .LE. 0. ) esw = eswpot
    700        IF ( eswdif .GT. 0. ) esw = eswphy
    701        IF ( esw  .GT. 0. )   esw = 0.
    702 
    703 !--    DIFFUSION
    704 
    705        rdsk = 0.79 * 10. ** 7.
    706        rdcl = 0.
    707        ed   = evap / ( rdsk + rdcl ) * adu * ( 1. - wetsk ) * ( vpa - vpts )
    708 
    709 !--    MAX VB
    710 
    711        vb1 = 34. - tsk
    712        vb2 = tcore(j) - 36.6
    713 
    714        IF ( vb2 .LT. 0. ) vb2 = 0.
    715        IF ( vb1 .LT. 0. ) vb1 = 0.
    716        vb = ( 6.3 + 75. * vb2 ) / ( 1. + 0.5 * vb1 )
    717 
    718 !--    ENERGY BALLANCE
    719 
    720        enbal = h + ed + ere + esw + csum + rsum + food
    721 
    722 
    723 !--    CLOTHING TEMPERATURE
    724 
    725        xx = 0.001
    726        IF ( count1 .EQ. 0 ) xx = 1.
    727        IF ( count1 .EQ. 1 ) xx = 0.1
    728        IF ( count1 .EQ. 2 ) xx = 0.01
    729        IF ( count1 .EQ. 3 ) xx = 0.001
    730 
    731        IF ( enbal .GT. 0. ) tcl = tcl + xx
    732        IF ( enbal .LT. 0. ) tcl = tcl - xx
    733 
    734        skipIncreaseCount = .FALSE.
    735        IF ( ( (enbal .LE. 0.) .AND. (enbal2 .GT. 0.) ) .OR.                    &
    736           ( ( enbal .GE. 0. ) .AND. ( enbal2 .LT. 0. ) ) ) THEN
    737          skipIncreaseCount = .TRUE.
    738        ELSE
    739          enbal2 = enbal
    740          count3 = count3 + 1
    741        END IF
    742        
    743        IF ( ( count3 .GT. 200 ) .OR. skipIncreaseCount ) THEN
    744          IF ( count1 .LT. 3 ) THEN
    745            count1 = count1 + 1
    746            enbal2 = 0.
    747          ELSE
    748            EXIT
    749          END IF
    750        END IF
    751      END DO
    752 
    753      IF ( count1 .EQ. 3 ) THEN
    754        SELECT CASE ( j )
    755          CASE ( 2, 5)
    756            IF ( .NOT. ( ( tcore(j) .GE. 36.6 ) .AND.                        &
    757              ( tsk .LE. 34.050 ) ) ) CYCLE
    758          CASE ( 6, 1 )
    759            IF ( c(10) .LT. 0. ) CYCLE
    760            IF ( .NOT. ( ( tcore(j) .GE. 36.6 ) .AND.                        &
    761              ( tsk .GT. 33.850 ) ) ) CYCLE
    762          CASE ( 3 )
    763            IF ( .NOT. ( ( tcore(j) .LT. 36.6 ) .AND.                        &
    764              ( tsk .LE. 34.000 ) ) ) CYCLE
    765          CASE ( 7 )
    766            IF ( .NOT. ( ( tcore(j) .LT. 36.6 ) .AND.                        &
    767              ( tsk .GT. 34.000 ) ) ) CYCLE
    768          CASE default  ! = CASE ( 4 ), does actually nothing
    769        END SELECT
    770      END IF
    771 
    772      IF ( ( j .NE. 4 ) .AND. ( vb .GE. 91. ) ) CYCLE
    773      IF ( ( j .EQ. 4 ) .AND. ( vb .LT. 89. ) ) CYCLE
    774      IF ( vb .GT. 90.) vb = 90.
    775 
    776 !--  LOSSES BY WATER
    777 
    778      ws = sw * 3600. * 1000.
    779      IF ( ws .GT. 2000. ) ws = 2000.
    780      wd = ed / evap * 3600. * ( -1000. )
    781      wr = erel / evap * 3600. * ( -1000. )
    782 
    783      wsum = ws + wr + wd
    784 
    785      RETURN
    786    END DO
    787    RETURN
    788  END SUBROUTINE Berech
    789 
    790 
    791  
    792 !------------------------------------------------------------------------------!
    793 ! Description:
    794 ! ------------
    795 !> Calculate PET
    796 !------------------------------------------------------------------------------!
    797  SUBROUTINE PET ( acl, adu, aeff, cair, emcl, emsk, esw, evap,                 &
    798      facl, feff, h, p, po, rdcl, rdsk, rtv, sigm, ta, tcl, tsk, tx, vpts, wetsk)
    799 
    800    REAL ( wp )  :: acl, adu, aeff, cair, cbare, cclo, csum, ed,                &
    801       emcl, emsk, enbal, enbal2, ere, erel, eres, esw, evap,                   &
    802       facl, feff, h, hc, p, po, rbare, rclo, rdcl, rdsk, rsum,                 &
    803       rtv, sigm, ta, tcl, tex, tsk, tx, vpex, vpts, wetsk, xx
    804 
    805    INTEGER ( iwp ) :: count1
    806 
    807    tx = ta
    808    enbal2 = 0.
    809 
    810     DO count1 = 0, 3
    811      DO
    812        hc = 2.67 + 6.5 * 0.1 ** 0.67
    813        hc = hc * ( p / po ) ** 0.55
    814 
    815 !--    Radiation
    816 
    817        aeff  = adu * feff
    818        rbare = aeff * ( 1. - facl ) * emsk * sigm *                            &
    819            ( ( tx + 273.2 ) ** 4. - ( tsk + 273.2 ) ** 4. )
    820        rclo  = feff * acl * emcl * sigm *                                      &
    821            ( ( tx + 273.2 ) ** 4. - ( tcl + 273.2 ) ** 4. )
    822        rsum  = rbare + rclo
    823 
    824 !--    Covection
    825 
    826        cbare = hc * ( tx - tsk ) * adu * ( 1. - facl )
    827        cclo  = hc * ( tx - tcl ) * acl
    828        csum  = cbare + cclo
    829 
    830 !--    Diffusion
    831 
    832        ed = evap / ( rdsk + rdcl ) * adu * ( 1. - wetsk ) * ( 12. - vpts )
    833 
    834 !--    Respiration
    835 
    836        tex  = 0.47 * tx + 21.
    837        eres = cair * ( tx - tex ) * rtv
    838        vpex = 6.11 * 10. ** ( 7.45 * tex / ( 235. + tex ) )
    839        erel = 0.623 * evap / p * ( 12. - vpex ) * rtv
    840        ere  = eres + erel
    841 
    842 !--    Energy ballance
    843 
    844        enbal = h + ed + ere + esw + csum + rsum
    845 
    846 !--    Iteration concerning ta
    847 
    848        IF ( count1 .EQ. 0 )  xx = 1.
    849        IF ( count1 .EQ. 1 )  xx = 0.1
    850        IF ( count1 .EQ. 2 )  xx = 0.01
    851        IF ( count1 .EQ. 3 )  xx = 0.001
    852        IF ( enbal .GT. 0. )  tx = tx - xx
    853        IF ( enbal .LT. 0. )  tx = tx + xx
    854        IF ( ( enbal .LE. 0. ) .AND. ( enbal2 .GT. 0. ) ) EXIT
    855        IF ( ( enbal .GE. 0. ) .AND. ( enbal2 .LT. 0. ) ) EXIT
    856 
    857        enbal2 = enbal
    858      END DO
    859    END DO
    860    RETURN
    861  END SUBROUTINE PET
    862 
    863 END MODULE biometeorology_mod
     692    CHARACTER (LEN=*), INTENT(OUT) ::  grid   !< Grid type (always "zu1" for biom)
     693    LOGICAL, INTENT(OUT)           ::  found  !< Output found?
     694    LOGICAL, INTENT(OUT)           ::  two_d  !< Flag parameter that indicates 2D variables, horizontal cross sections, must be .TRUE.
     695    REAL(wp), DIMENSION(nxl:nxr,nys:nyn,nzb_do:nzt_do) ::  local_pf  !< Temp. result grid to return
     696
     697!-- Internal variables
     698    CHARACTER (LEN=:), allocatable ::  variable_short  !< Trimmed version of char variable
     699    INTEGER(iwp) ::  i        !< Running index, x-dir
     700    INTEGER(iwp) ::  j        !< Running index, y-dir
     701    INTEGER(iwp) ::  k        !< Running index, z-dir
     702
     703
     704    found = .TRUE.
     705    variable_short = TRIM( variable )
     706    IF ( variable_short(1:4) == 'biom' ) THEN
     707       two_d = .TRUE.
     708       grid = 'zu1'
     709    ELSE
     710       found = .FALSE.
     711       grid  = 'none'
     712       RETURN
     713    ENDIF
     714
     715    local_pf = fill_value
     716
     717    SELECT CASE ( variable_short )
     718
     719
     720        CASE ( 'biom_tmrt_xy' )        ! 2d-array
     721              DO  i = nxl, nxr
     722                 DO  j = nys, nyn
     723                    local_pf(i,j,1) = tmrt_grid(j,i)
     724                 ENDDO
     725              ENDDO
     726
     727
     728        CASE ( 'biom_pt_xy' )        ! 2d-array
     729            IF ( av == 0 )  THEN
     730              DO  i = nxl, nxr
     731                 DO  j = nys, nyn
     732                    local_pf(i,j,nzb+1) = pt_grid(j,i)
     733                 ENDDO
     734              ENDDO
     735            ELSE
     736              DO  i = nxl, nxr
     737                 DO  j = nys, nyn
     738                    local_pf(i,j,nzb+1) = pt_av_grid(j,i)
     739                 ENDDO
     740              ENDDO
     741            END IF
     742
     743
     744        CASE ( 'biom_utci_xy' )        ! 2d-array
     745            IF ( av == 0 )  THEN
     746              DO  i = nxl, nxr
     747                 DO  j = nys, nyn
     748                    local_pf(i,j,nzb+1) = utci_grid(j,i)
     749                 ENDDO
     750              ENDDO
     751            ELSE
     752              DO  i = nxl, nxr
     753                 DO  j = nys, nyn
     754                    local_pf(i,j,nzb+1) = utci_av_grid(j,i)
     755                 ENDDO
     756              ENDDO
     757            END IF
     758
     759
     760        CASE ( 'biom_pet_xy' )        ! 2d-array
     761            IF ( av == 0 )  THEN
     762              DO  i = nxl, nxr
     763                 DO  j = nys, nyn
     764                    local_pf(i,j,nzb+1) = pet_grid(j,i)
     765                 ENDDO
     766              ENDDO
     767            ELSE
     768              DO  i = nxl, nxr
     769                 DO  j = nys, nyn
     770                    local_pf(i,j,nzb+1) = pet_av_grid(j,i)
     771                 ENDDO
     772              ENDDO
     773            END IF
     774
     775
     776       CASE DEFAULT
     777          found = .FALSE.
     778          grid  = 'none'
     779
     780    END SELECT
     781
     782
     783 END SUBROUTINE biom_data_output_2d
     784
     785
     786!------------------------------------------------------------------------------!
     787! Description:
     788! ------------
     789!> Subroutine defining appropriate grid for netcdf variables.
     790!> It is called out from subroutine netcdf_interface_mod.
     791!> netcdf_interface_mod 918ff
     792!------------------------------------------------------------------------------!
     793 SUBROUTINE biom_define_netcdf_grid( var, found, grid_x, grid_y, grid_z )
     794
     795    IMPLICIT NONE
     796
     797!-- Input variables
     798    CHARACTER (LEN=*), INTENT(IN)  ::  var      !< Name of output variable
     799
     800!-- Output variables
     801    CHARACTER (LEN=*), INTENT(OUT) ::  grid_x   !< x grid of output variable
     802    CHARACTER (LEN=*), INTENT(OUT) ::  grid_y   !< y grid of output variable
     803    CHARACTER (LEN=*), INTENT(OUT) ::  grid_z   !< z grid of output variable
     804
     805    LOGICAL, INTENT(OUT)           ::  found    !< Flag if output var is found
     806
     807!-- Local variables
     808    LOGICAL      :: is2d  !< Var is 2d?
     809
     810    INTEGER(iwp) :: l     !< Length of the var array
     811
     812
     813    found  = .FALSE.
     814    grid_x = 'none'
     815    grid_y = 'none'
     816    grid_z = 'none'
     817
     818    l = MAX( 2, LEN_TRIM( var ) )
     819    is2d = ( var(l-1:l) == 'xy' )
     820
     821
     822    IF ( var(1:4) == 'biom' ) THEN
     823       found  = .TRUE.
     824       grid_x = 'x'
     825       grid_y = 'y'
     826       grid_z = 'zu'
     827       IF ( is2d ) grid_z = 'zu1'
     828    ENDIF
     829
     830 END SUBROUTINE biom_define_netcdf_grid
     831
     832!------------------------------------------------------------------------------!
     833! Description:
     834! ------------
     835!> Header output for biom module
     836!> header 982
     837!------------------------------------------------------------------------------!
     838 SUBROUTINE biom_header( io )
     839
     840    IMPLICIT NONE
     841
     842!-- Input variables
     843    INTEGER(iwp), INTENT(IN) ::  io           !< Unit of the output file
     844
     845!-- Internal variables
     846    CHARACTER (LEN=86) ::  output_height_chr  !< String for output height
     847
     848    WRITE( output_height_chr, '(F8.1,7X)' )  biom_output_height
     849!
     850!-- Write biom header
     851    WRITE( io, 1 )
     852    WRITE( io, 2 )  TRIM( output_height_chr )
     853    WRITE( io, 3 )  TRIM( ACHAR( biom_cell_level ) )
     854
     8551   FORMAT (//' Human thermal comfort module information:'/                    &
     856              ' ------------------------------'/)
     8572   FORMAT ('    --> All indices calculated for a height of (m): ', A )
     8583   FORMAT ('    --> This corresponds to cell level : ', A )
     859
     860 END SUBROUTINE biom_header
     861
     862
     863!------------------------------------------------------------------------------!
     864! Description:
     865! ------------
     866!> Initialization of the HTCM
     867!> init_3d_model 1987ff
     868!------------------------------------------------------------------------------!
     869 SUBROUTINE biom_init
     870
     871    IMPLICIT NONE
     872
     873!-- Internal vriables
     874    REAL ( wp )  :: height  !< current height in meters
     875
     876    INTEGER ( iwp )  :: i  !< iteration index
     877
     878!-- Determine cell level corresponding to 1.1 m above ground level
     879!   (gravimetric center of sample human)
     880
     881    time_biom_results = 0.0_wp
     882    biom_cell_level = 0_iwp
     883    biom_output_height = 0.5_wp * dz(1)
     884    height = 0.0_wp
     885
     886    biom_cell_level = INT ( 1.099_wp / dz(1) )
     887    biom_output_height = biom_output_height + biom_cell_level * dz(1)
     888
     889 END SUBROUTINE biom_init
     890
     891
     892!------------------------------------------------------------------------------!
     893! Description:
     894! ------------
     895!> Allocate biom arrays and define pointers if required
     896!> init_3d_model 1047ff
     897!------------------------------------------------------------------------------!
     898 SUBROUTINE biom_init_arrays
     899
     900    IMPLICIT NONE
     901
     902!-- Allocate a temporary array with the desired output dimensions.
     903!   Initialization omitted for performance, will be overwritten anyway
     904    IF ( .NOT. ALLOCATED( tmrt_grid ) ) THEN
     905      ALLOCATE( tmrt_grid (nys:nyn,nxl:nxr) )
     906    ENDIF
     907
     908    IF ( biom_pt ) THEN
     909      IF ( .NOT. ALLOCATED( pt_grid ) ) THEN
     910        ALLOCATE( pt_grid (nys:nyn,nxl:nxr) )
     911      ENDIF
     912    ENDIF
     913
     914    IF ( biom_utci ) THEN
     915      IF ( .NOT. ALLOCATED( utci_grid ) ) THEN
     916        ALLOCATE( utci_grid (nys:nyn,nxl:nxr) )
     917      ENDIF
     918    ENDIF
     919
     920    IF ( biom_pet ) THEN
     921      IF ( .NOT. ALLOCATED( pet_grid ) ) THEN
     922        ALLOCATE( pet_grid (nys:nyn,nxl:nxr) )
     923      ENDIF
     924    END IF
     925
     926    IF ( biom_pt_av ) THEN
     927      IF ( .NOT. ALLOCATED( pt_av_grid ) ) THEN
     928        ALLOCATE( pt_av_grid (nys:nyn,nxl:nxr) )
     929      ENDIF
     930    ENDIF
     931
     932    IF ( biom_utci_av ) THEN
     933      IF ( .NOT. ALLOCATED( utci_av_grid ) ) THEN
     934        ALLOCATE( utci_av_grid (nys:nyn,nxl:nxr) )
     935      ENDIF
     936    ENDIF
     937
     938    IF ( biom_pet_av ) THEN
     939      IF ( .NOT. ALLOCATED( pet_av_grid ) ) THEN
     940        ALLOCATE( pet_av_grid (nys:nyn,nxl:nxr) )
     941      ENDIF
     942
     943    END IF
     944
     945 END SUBROUTINE biom_init_arrays
     946
     947
     948!------------------------------------------------------------------------------!
     949! Description:
     950! ------------
     951!> Parin for &biometeorology_parameters for reading biomet parameters
     952!------------------------------------------------------------------------------!
     953 SUBROUTINE biom_parin
     954
     955    IMPLICIT NONE
     956
     957!
     958!-- Internal variables
     959    CHARACTER (LEN=80) ::  line  !< Dummy string for current line in parameter file
     960
     961    NAMELIST /biometeorology_parameters/  biom_pet,                            &
     962                                          biom_pet_av,                         &
     963                                          biom_pt,                             &
     964                                          biom_pt_av,                          &
     965                                          biom_utci,                           &
     966                                          biom_utci_av
     967
     968
     969!-- Try to find biometeorology_parameters namelist
     970    REWIND ( 11 )
     971    line = ' '
     972    DO WHILE ( INDEX( line, '&biometeorology_parameters' ) == 0 )
     973       READ ( 11, '(A)', END = 20 )  line
     974    ENDDO
     975    BACKSPACE ( 11 )
     976
     977!
     978!-- Read biometeorology_parameters namelist
     979    READ ( 11, biometeorology_parameters, ERR = 10, END = 20 )
     980
     981!
     982!-- Set flag that indicates that the biomet_module is switched on
     983    biometeorology = .TRUE.
     984
     985    GOTO 20
     986
     987!
     988!-- In case of error
     989 10 BACKSPACE( 11 )
     990    READ( 11 , '(A)') line
     991    CALL parin_fail_message( 'biometeorology_parameters', line )
     992
     993!
     994!-- Complete
     995 20 CONTINUE
     996
     997
     998 END SUBROUTINE biom_parin
     999
     1000!------------------------------------------------------------------------------!
     1001! Description:
     1002! ------------
     1003!> Calculates the mean radiant temperature (tmrt) based on the Six-directions
     1004!> method according to VDI 3787 2.
     1005!------------------------------------------------------------------------------!
     1006 SUBROUTINE calculate_tmrt_6_directions( SW_N, SW_E, SW_S, SW_W,               &
     1007    SW_U, SW_D, LW_N, LW_E, LW_S, LW_W, LW_U, LW_D, tmrt )
     1008
     1009    IMPLICIT NONE
     1010
     1011!-- Type of input of the argument list
     1012!   Short- (SW_) and longwave (LW_) radiation fluxes from the six directions
     1013!   North (N), East (E), South (S), West (W), up (U) and down (D)
     1014    REAL(wp), INTENT ( IN )  ::  SW_N  !< Sw radflux density from N    (W/m²)
     1015    REAL(wp), INTENT ( IN )  ::  SW_E  !< Sw radflux density from E    (W/m²)
     1016    REAL(wp), INTENT ( IN )  ::  SW_S  !< Sw radflux density from S    (W/m²)
     1017    REAL(wp), INTENT ( IN )  ::  SW_W  !< Sw radflux density from W    (W/m²)
     1018    REAL(wp), INTENT ( IN )  ::  SW_U  !< Sw radflux density from U    (W/m²)
     1019    REAL(wp), INTENT ( IN )  ::  SW_D  !< Sw radflux density from D    (W/m²)
     1020    REAL(wp), INTENT ( IN )  ::  LW_N  !< Lw radflux density from N    (W/m²)
     1021    REAL(wp), INTENT ( IN )  ::  LW_E  !< Lw radflux density from E    (W/m²)
     1022    REAL(wp), INTENT ( IN )  ::  LW_S  !< Lw radflux density from S    (W/m²)
     1023    REAL(wp), INTENT ( IN )  ::  LW_W  !< Lw radflux density from W    (W/m²)
     1024    REAL(wp), INTENT ( IN )  ::  LW_U  !< Lw radflux density from U    (W/m²)
     1025    REAL(wp), INTENT ( IN )  ::  LW_D  !< Lw radflux density from D    (W/m²)
     1026
     1027!-- Type of output of the argument list
     1028    REAL(wp), INTENT ( OUT ) ::  tmrt  !< Mean radiant temperature (°C)
     1029
     1030!-- Directional weighting factors
     1031    REAL(wp), PARAMETER      ::  weight_h  = 0.22_wp
     1032    REAL(wp), PARAMETER      ::  weight_v  = 0.06_wp
     1033
     1034    REAL(wp) ::  nrfd           !<  Net radiation flux density (W/m²)
     1035
     1036!-- Initialise
     1037    nrfd = 0._wp
     1038
     1039!-- Compute mean radiation flux density absorbed by rotational symetric human
     1040    nrfd = ( weight_h * ( human_absorb * SW_N + human_emiss * LW_N ) ) +       &
     1041       ( weight_h * ( human_absorb * SW_E + human_emiss * LW_E ) ) +           &
     1042       ( weight_h * ( human_absorb * SW_S + human_emiss * LW_S ) ) +           &
     1043       ( weight_h * ( human_absorb * SW_W + human_emiss * LW_W ) ) +           &
     1044       ( weight_v * ( human_absorb * SW_U + human_emiss * LW_U ) ) +           &
     1045       ( weight_v * ( human_absorb * SW_D + human_emiss * LW_D ) )
     1046
     1047!-- Compute mean radiant temperature
     1048    tmrt = ( nrfd / (human_emiss * sigma_sb) )**0.25_wp - cels_offs
     1049
     1050 END SUBROUTINE calculate_tmrt_6_directions
     1051
     1052!------------------------------------------------------------------------------!
     1053! Description:
     1054! ------------
     1055!> Very crude approximation of mean radiant temperature based on upwards and
     1056!> downwards radiation fluxes only (other directions curently not available,
     1057!> replace as soon as possible!)
     1058!------------------------------------------------------------------------------!
     1059 SUBROUTINE calculate_tmrt_2_directions( sw_u, sw_d, lw_u, lw_d, ta, tmrt )
     1060
     1061    IMPLICIT NONE
     1062
     1063!-- Type of input of the argument list
     1064    REAL(wp), INTENT ( IN )  ::  sw_u  !< Shortwave radiation flux density from upper direction (W/m²)
     1065    REAL(wp), INTENT ( IN )  ::  sw_d  !< Shortwave radiation flux density from lower direction (W/m²)
     1066    REAL(wp), INTENT ( IN )  ::  lw_u  !< Longwave radiation flux density from upper direction (W/m²)
     1067    REAL(wp), INTENT ( IN )  ::  lw_d  !< Longwave radiation flux density from lower direction (W/m²)
     1068    REAL(wp), INTENT ( IN )  ::  ta    !< Air temperature (°C)
     1069
     1070!-- Type of output of the argument list
     1071    REAL(wp), INTENT ( OUT ) ::  tmrt  !< mean radiant temperature, (°C)
     1072
     1073!-- Directional weighting factors and parameters
     1074    REAL(wp), PARAMETER ::  weight_h  = 0.22_wp     !< Weight for horizontal radiational gain after Fanger (1972)
     1075    REAL(wp), PARAMETER ::  weight_v  = 0.06_wp     !< Weight for vertical radiational gain after Fanger (1972)
     1076
     1077!-- Other internal variables
     1078    REAL(wp) ::  sw_in
     1079    REAL(wp) ::  sw_out
     1080    REAL(wp) ::  lw_in
     1081    REAL(wp) ::  lw_out
     1082    REAL(wp) ::  nrfd     !< Net radiation flux density (W/m²)
     1083    REAL(wp) ::  lw_air   !< Longwave emission by surrounding air volume (W/m²)
     1084    REAL(wp) ::  sw_side  !< Shortwave immission from the sides (W/m²)
     1085
     1086    INTEGER(iwp) ::  no_input  !< Count missing input radiation fluxes
     1087
     1088!-- initialise
     1089    sw_in    = sw_u
     1090    sw_out   = sw_d
     1091    lw_in    = lw_u
     1092    lw_out   = lw_d
     1093    nrfd     = 0._wp
     1094    no_input = 0_iwp
     1095
     1096!-- test for missing input data
     1097    IF ( sw_in <= -998._wp .OR. sw_out <= -998._wp .OR. lw_in <= -998._wp .OR. &
     1098       lw_out <= -998._wp .OR. ta <= -998._wp ) THEN
     1099       IF ( sw_in <= -998._wp ) THEN
     1100          sw_in = 0.
     1101          no_input = no_input + 1
     1102       ENDIF
     1103       IF ( sw_out <= -998._wp ) THEN
     1104          sw_out = 0.
     1105          no_input = no_input + 1
     1106       ENDIF
     1107       IF ( lw_in <= -998._wp ) THEN
     1108          lw_in = 0.
     1109          no_input = no_input + 1
     1110       ENDIF
     1111       IF ( lw_out <= -998._wp ) THEN
     1112          lw_out = 0.
     1113          no_input = no_input + 1
     1114       ENDIF
     1115
     1116!-- Accept two missing radiation flux directions, fail otherwise as error might become too large
     1117       IF ( ta <= -998._wp .OR. no_input >= 3 ) THEN
     1118          tmrt = -999._wp
     1119          RETURN
     1120       ENDIF
     1121    ENDIF
     1122
     1123    sw_side = sw_in * 0.125_wp  ! distribute half of upper sw_in to the 4 sides
     1124    lw_air = ( sigma_sb * 0.95_wp * ( ta + cels_offs )**4 )
     1125
     1126!-- Compute mean radiation flux density absorbed by rotational symetric human
     1127    nrfd = ( weight_h * ( human_absorb * sw_side + human_emiss * lw_air ) ) +         &
     1128       ( weight_h * ( human_absorb * sw_side + human_emiss * lw_air ) ) +             &
     1129       ( weight_h * ( human_absorb * sw_side + human_emiss * lw_air ) ) +             &
     1130       ( weight_h * ( human_absorb * sw_side + human_emiss * lw_air ) ) +             &
     1131       ( weight_v * ( human_absorb * (sw_in * 0.5_wp) + human_emiss * lw_in ) ) +     &
     1132       ( weight_v * ( human_absorb * sw_out + human_emiss * lw_out ) )
     1133
     1134!-- Compute mean radiant temperature
     1135    tmrt = ( nrfd / (human_emiss * sigma_sb) )**0.25_wp - cels_offs
     1136
     1137 END SUBROUTINE calculate_tmrt_2_directions
     1138
     1139!------------------------------------------------------------------------------!
     1140! Description:
     1141! ------------
     1142!> Calculate static thermal indices for given meteorological conditions
     1143!------------------------------------------------------------------------------!
     1144 SUBROUTINE calculate_static_thermal_indices ( ta, vp, ws, pair, tmrt,         &
     1145    pt_static, utci_static, pet_static )
     1146
     1147    IMPLICIT NONE
     1148
     1149!-- Input parameters
     1150    REAL(wp), INTENT ( IN )  ::  ta           !< Air temperature                  (°C)
     1151    REAL(wp), INTENT ( IN )  ::  vp           !< Vapour pressure                  (hPa)
     1152    REAL(wp), INTENT ( IN )  ::  ws           !< Wind speed    (local level)      (m/s)
     1153    REAL(wp), INTENT ( IN )  ::  pair         !< Air pressure                     (hPa)
     1154    REAL(wp), INTENT ( IN )  ::  tmrt         !< Mean radiant temperature         (°C)
     1155!-- Output parameters
     1156    REAL(wp), INTENT ( OUT ) ::  pt_static    !< Perceived temperature            (°C)
     1157    REAL(wp), INTENT ( OUT ) ::  utci_static  !< Universal thermal climate index  (°C)
     1158    REAL(wp), INTENT ( OUT ) ::  pet_static   !< Physiologically equivalent temp. (°C)
     1159!-- Temporary field, not used here
     1160    REAL(wp)                 ::  clo          !< Clothing index (no dim.)
     1161
     1162    clo = -999._wp
     1163
     1164    IF ( biom_pt ) THEN
     1165!-- Estimate local perceived temperature
     1166       CALL calculate_pt_static( ta, vp, ws, tmrt, pair, clo, pt_static )
     1167    ENDIF
     1168
     1169    IF ( biom_utci ) THEN
     1170!-- Estimate local universal thermal climate index
     1171       CALL calculate_utci_static( ta, vp, ws, tmrt, biom_output_height,       &
     1172          utci_static )
     1173    ENDIF
     1174
     1175    IF ( biom_pet ) THEN
     1176!-- Estimate local physiologically equivalent temperature
     1177       CALL calculate_pet_static( ta, vp, ws, tmrt, pair, pet_static )
     1178    ENDIF
     1179
     1180 END SUBROUTINE calculate_static_thermal_indices
     1181
     1182
     1183!------------------------------------------------------------------------------!
     1184! Description:
     1185! ------------
     1186!> Calculate static thermal indices for 2D grid point i, j
     1187!------------------------------------------------------------------------------!
     1188 SUBROUTINE biom_determine_input_at( average_input, i, j, ta, vp, ws, pair,    &
     1189    tmrt )
     1190
     1191    IMPLICIT NONE
     1192
     1193!-- Input variables
     1194    LOGICAL,      INTENT ( IN ) ::  average_input  !< Determine averaged input conditions?
     1195    INTEGER(iwp), INTENT ( IN ) ::  i     !< Running index, x-dir
     1196    INTEGER(iwp), INTENT ( IN ) ::  j     !< Running index, y-dir
     1197
     1198!-- Output parameters
     1199    REAL(wp), INTENT ( OUT )    ::  tmrt  !< Mean radiant temperature        (°C)
     1200    REAL(wp), INTENT ( OUT )    ::  ta    !< Air temperature                 (°C)
     1201    REAL(wp), INTENT ( OUT )    ::  vp    !< Vapour pressure                 (hPa)
     1202    REAL(wp), INTENT ( OUT )    ::  ws    !< Wind speed    (local level)     (m/s)
     1203    REAL(wp), INTENT ( OUT )    ::  pair  !< Air pressure                    (hPa)
     1204
     1205!-- Internal variables
     1206    INTEGER(iwp)                ::  k     !< Running index, z-dir
     1207    INTEGER(iwp)                ::  k_wind  !< Running index, z-dir, wind speed only
     1208
     1209    REAL(wp)                    ::  vp_sat  !< Saturation vapor pressure     (hPa)
     1210
     1211
     1212!-- Determine cell level closest to 1.1m above ground
     1213!   by making use of truncation due to int cast
     1214    k = get_topography_top_index_ji(j, i, 's') + biom_cell_level  !< Vertical cell center closest to 1.1m
     1215    k_wind = k
     1216    IF( k_wind < 1_iwp ) THEN  ! Avoid horizontal u and v of 0.0 m/s close to ground
     1217       k_wind = 1_iwp
     1218    ENDIF
     1219
     1220!-- Determine local values:
     1221    IF ( average_input ) THEN
     1222!--    Calculate ta from Tp assuming dry adiabatic laps rate
     1223       ta = pt_av(k, j, i) - ( 0.0098_wp * dz(1) * (  k + .5_wp ) ) - cels_offs
     1224
     1225       vp = 0.034_wp
     1226       IF ( humidity .AND. ALLOCATED( q_av ) ) THEN
     1227          vp = q_av(k, j, i)
     1228       ENDIF
     1229
     1230       ws = ( 0.5_wp * ABS( u_av(k_wind, j, i) + u_av(k_wind, j, i+1) )  +     &
     1231          0.5_wp * ABS( v_av(k_wind, j, i) + v_av(k_wind, j+1, i) )  +         &
     1232          0.5_wp * ABS( w_av(k_wind, j, i) + w_av(k_wind+1, j, i) ) )
     1233    ELSE
     1234!-- Calculate ta from Tp assuming dry adiabatic laps rate
     1235       ta = pt(k, j, i) - ( 0.0098_wp * dz(1) * (  k + .5_wp ) ) - cels_offs
     1236
     1237       vp = q(k, j, i)
     1238
     1239       ws = ( 0.5_wp * ABS( u(k_wind, j, i) + u(k_wind, j, i+1) )  +           &
     1240          0.5_wp * ABS( v(k_wind, j, i) + v(k_wind, j+1, i) )  +               &
     1241          0.5_wp * ABS( w(k_wind, j, i) + w(k_wind+1, j, i) ) )
     1242
     1243    ENDIF
     1244
     1245!-- Local air pressure
     1246    pair = surface_pressure
     1247!
     1248!-- Calculate water vapour pressure at saturation and convert to hPa
     1249!-- The magnus formula is limited to temperatures up to 333.15 K to
     1250!   avoid negative values of vp_sat
     1251    vp_sat = 0.01_wp * magnus( MIN( ta + cels_offs, 333.15_wp ) )
     1252    vp  = vp * pair / ( vp + 0.622_wp )
     1253    IF ( vp > vp_sat ) vp = vp_sat
     1254
     1255    tmrt = ta
     1256    IF ( radiation ) THEN
     1257       CALL calculate_tmrt_2_directions (rad_sw_in(0, j, i),                   &
     1258          rad_sw_out(0, j, i), rad_lw_in(0, j, i), rad_lw_out(0, j, i), ta,    &
     1259          tmrt )
     1260    ENDIF
     1261
     1262 END SUBROUTINE biom_determine_input_at
     1263
     1264
     1265!------------------------------------------------------------------------------!
     1266! Description:
     1267! ------------
     1268!> Calculate static thermal indices for any point within a 2D grid
     1269!> time_integration.f90: 1065ff
     1270!------------------------------------------------------------------------------!
     1271 SUBROUTINE biom_calculate_static_grid ( average_input )
     1272
     1273    IMPLICIT NONE
     1274
     1275!-- Input attributes
     1276    LOGICAL, INTENT ( IN ) ::  average_input  !< Calculate based on averaged input? conditions?
     1277
     1278!-- Internal variables
     1279    INTEGER(iwp) ::  i, j      !< Running index
     1280
     1281    REAL(wp) ::  ta            !< Air temperature                  (°C)
     1282    REAL(wp) ::  vp            !< Vapour pressure                  (hPa)
     1283    REAL(wp) ::  ws            !< Wind speed    (local level)      (m/s)
     1284    REAL(wp) ::  pair          !< Air pressure                     (hPa)
     1285    REAL(wp) ::  tmrt_tmp      !< Mean radiant temperature
     1286    REAL(wp) ::  pt_tmp        !< Perceived temperature
     1287    REAL(wp) ::  utci_tmp      !< Universal thermal climate index
     1288    REAL(wp) ::  pet_tmp       !< Physiologically equivalent temperature
     1289
     1290
     1291    CALL biom_init_arrays
     1292
     1293    DO i = nxl, nxr
     1294       DO j = nys, nyn
     1295!--       Determine local input conditions
     1296          CALL biom_determine_input_at ( average_input, i, j, ta, vp, ws,      &
     1297             pair, tmrt_tmp )
     1298          tmrt_grid(j, i) = tmrt_tmp
     1299
     1300!--       Only proceed if tmrt is available
     1301          IF ( tmrt_tmp <= -998._wp ) THEN
     1302             pt_tmp   = -999._wp
     1303             utci_tmp = -999._wp
     1304             pet_tmp  = -999._wp
     1305             CYCLE
     1306          END IF
     1307
     1308!--       Calculate static thermal indices based on local tmrt
     1309          CALL calculate_static_thermal_indices ( ta, vp, ws,                  &
     1310             pair, tmrt_tmp, pt_tmp, utci_tmp, pet_tmp )
     1311
     1312          IF ( average_input ) THEN
     1313!--          Write results for selected averaged indices
     1314             IF ( biom_pt_av )  THEN
     1315                pt_av_grid(j, i)   = pt_tmp
     1316             END IF
     1317             IF ( biom_utci_av ) THEN
     1318                utci_av_grid(j, i) = utci_tmp
     1319             END IF
     1320             IF ( biom_pet_av ) THEN
     1321                pet_av_grid(j, i)  = pet_tmp
     1322             END IF
     1323          ELSE
     1324!--          Write result for selected indices
     1325             IF ( biom_pt )  THEN
     1326                pt_grid(j, i)   = pt_tmp
     1327             END IF
     1328             IF ( biom_utci ) THEN
     1329                utci_grid(j, i) = utci_tmp
     1330             END IF
     1331             IF ( biom_pet ) THEN
     1332                pet_grid(j, i)  = pet_tmp
     1333             END IF
     1334          END IF
     1335
     1336       END DO
     1337    END DO
     1338
     1339 END SUBROUTINE biom_calculate_static_grid
     1340
     1341!------------------------------------------------------------------------------!
     1342! Description:
     1343! ------------
     1344!> Calculate dynamic thermal indices (currently only iPT, but expandable)
     1345!------------------------------------------------------------------------------!
     1346 SUBROUTINE biom_calc_ipt( ta, vp, ws, pair, tmrt, dt, energy_storage,         &
     1347    t_clo, clo, actlev, age, weight, height, work, sex, ipt )
     1348
     1349    IMPLICIT NONE
     1350
     1351!-- Input parameters
     1352    REAL(wp), INTENT ( IN )  ::  ta   !< Air temperature                  (°C)
     1353    REAL(wp), INTENT ( IN )  ::  vp   !< Vapour pressure                  (hPa)
     1354    REAL(wp), INTENT ( IN )  ::  ws   !< Wind speed    (local level)      (m/s)
     1355    REAL(wp), INTENT ( IN )  ::  pair !< Air pressure                     (hPa)
     1356    REAL(wp), INTENT ( IN )  ::  tmrt !< Mean radiant temperature         (°C)
     1357    REAL(wp), INTENT ( IN )  ::  dt   !< Time past since last calculation (s)
     1358    REAL(wp), INTENT ( IN )  ::  age  !< Age of agent                     (y)
     1359    REAL(wp), INTENT ( IN )  ::  weight  !< Weight of agent               (Kg)
     1360    REAL(wp), INTENT ( IN )  ::  height  !< Height of agent               (m)
     1361    REAL(wp), INTENT ( IN )  ::  work    !< Mechanical workload of agent
     1362                                         !  (without metabolism!)         (W)
     1363    INTEGER(iwp), INTENT ( IN ) ::  sex  !< Sex of agent (1 = male, 2 = female)
     1364
     1365!-- Both, input and output parameters
     1366    Real(wp), INTENT ( INOUT )  ::  energy_storage    !< Energy storage   (W/m²)
     1367    Real(wp), INTENT ( INOUT )  ::  t_clo   !< Clothing temperature       (°C)
     1368    Real(wp), INTENT ( INOUT )  ::  clo     !< Current clothing in sulation
     1369    Real(wp), INTENT ( INOUT )  ::  actlev  !< Individuals activity level
     1370                                            !  per unit surface area      (W/m²)
     1371!-- Output parameters
     1372    REAL(wp), INTENT ( OUT ) ::  ipt    !< Instationary perceived temp.   (°C)
     1373
     1374!-- If clo equals the initial value, this is the initial call
     1375    IF ( clo <= -998._wp ) THEN
     1376!--    Initialize instationary perceived temperature with personalized
     1377!      PT as an initial guess, set actlev and clo
     1378       CALL ipt_init ( age, weight, height, sex, work, actlev, clo,            &
     1379          ta, vp, ws, tmrt, pair, dt, energy_storage, t_clo,                   &
     1380          ipt )
     1381    ELSE
     1382!--    Estimate local instatinoary perceived temperature
     1383       CALL ipt_cycle ( ta, vp, ws, tmrt, pair, dt, energy_storage, t_clo,     &
     1384          clo, actlev, work, ipt )
     1385    ENDIF
     1386
     1387 END SUBROUTINE biom_calc_ipt
     1388
     1389 END MODULE biometeorology_mod

  • palm/trunk/SOURCE/check_parameters.f90

    r3421 r3448  
    2525! -----------------
    2626! $Id$
     27! Add biometeorology
     28!
     29! 3421 2018-10-24 18:39:32Z gronemeier
    2730! Renamed output variables
    2831! Add checks for surface data output
    … …  
    704707    USE basic_constants_and_equations_mod
    705708
     709    USE biometeorology_mod,                                                    &
     710        ONLY:  biom_check_data_output, biom_check_parameters
     711
    706712    USE bulk_cloud_model_mod,                                                  &
    707713        ONLY:  bulk_cloud_model, bcm_check_parameters, bcm_check_data_output,  &
    … …  
    769775        ONLY:  radiation, radiation_check_data_output,                         &
    770776               radiation_check_data_output_pr, radiation_check_parameters
    771 
    772     USE biometeorology_mod,                                                    &
    773         ONLY:  biometeorology_check_data_output
    774777
    775778    USE spectra_mod,                                                           &
    … …  
    14881491
    14891492!-- Check the module settings
     1493    IF ( biometeorology )       CALL biom_check_parameters
    14901494    IF ( bulk_cloud_model )     CALL bcm_check_parameters
    14911495    IF ( gust_module_enabled )  CALL gust_check_parameters
    … …  
    31823186             ENDIF
    31833187
    3184              IF ( unit == 'illegal' )  THEN
    3185                 CALL biometeorology_check_data_output( var, unit, i, ilen, k )
    3186              ENDIF
    3187 
    31883188             IF ( unit == 'illegal'  .AND.  gust_module_enabled  )  THEN
    31893189                CALL gust_check_data_output ( var, unit )
     3190             ENDIF
     3191
     3192             IF ( unit == 'illegal'  .AND.  biometeorology )  THEN
     3193                CALL biom_check_data_output( var, unit )
    31903194             ENDIF
    31913195

  • palm/trunk/SOURCE/data_output_2d.f90

    r3421 r3448  
    2525! -----------------
    2626! $Id$
     27! Add biometeorology
     28!
     29! 3421 2018-10-24 18:39:32Z gronemeier
    2730! Renamed output variables
    2831!
    … …  
    256259        ONLY:  c_p, lv_d_cp, l_v
    257260
     261    USE biometeorology_mod,                                                    &
     262        ONLY:  biom_data_output_2d
     263
    258264    USE bulk_cloud_model_mod,                                                  &
    259265        ONLY:  bulk_cloud_model, bcm_data_output_2d
    … …  
    263269
    264270    USE control_parameters,                                                    &
    265         ONLY:  air_chemistry, data_output_2d_on_each_pe, data_output_xy,       &
    266                data_output_xz, data_output_yz, do2d,                           &
     271        ONLY:  air_chemistry, biometeorology, data_output_2d_on_each_pe,       &
     272               data_output_xy, data_output_xz, data_output_yz, do2d,           &
    267273               do2d_xy_last_time, do2d_xy_time_count,                          &
    268274               do2d_xz_last_time, do2d_xz_time_count,                          &
    … …  
    13081314                   CALL gust_data_output_2d( av, do2d(av,if), found, grid,     &
    13091315                                             local_pf, two_d, nzb_do, nzt_do )
     1316                ENDIF
     1317
     1318                IF ( .NOT. found  .AND.  biometeorology                        &
     1319                                  .AND.  mode == 'xy' )  THEN
     1320                   CALL biom_data_output_2d( av, do2d(av,if), found, grid,     &
     1321                                             local_pf, two_d, nzb_do, nzt_do,  &
     1322                                             fill_value )
    13101323                ENDIF
    13111324

  • palm/trunk/SOURCE/data_output_3d.f90

    r3421 r3448  
    2525! -----------------
    2626! $Id$
     27! Adjustment of biometeorology calls
     28!
     29! 3421 2018-10-24 18:39:32Z gronemeier
    2730! Renamed output variables
    2831!
    … …  
    229232        ONLY:  lv_d_cp
    230233
     234    USE biometeorology_mod,                                                    &
     235        ONLY:  biom_data_output_3d
     236
    231237    USE bulk_cloud_model_mod,                                                  &
    232238        ONLY:  bulk_cloud_model, bcm_data_output_3d
    … …  
    236242
    237243    USE control_parameters,                                                    &
    238         ONLY:  air_chemistry, do3d, do3d_no, do3d_time_count,                  &
     244        ONLY:  air_chemistry, biometeorology, do3d, do3d_no, do3d_time_count,  &
    239245               io_blocks, io_group, land_surface, message_string,              &
    240246               ntdim_3d, nz_do3d,  ocean_mode, plant_canopy,                   &
    … …  
    285291    USE radiation_model_mod,                                                   &
    286292        ONLY:  nzub, nzut, radiation, radiation_data_output_3d
    287 
    288     USE biometeorology_mod,                                                    &
    289         ONLY: biometeorology_data_output_3d
    290293
    291294    USE turbulence_closure_mod,                                                &
    … …  
    767770             ENDIF
    768771
    769 !
    770 !--          Biometeorology quantity
    771              IF ( .NOT. found  .AND.  radiation )  THEN
    772                 CALL biometeorology_data_output_3d( av, do3d(av,if), found,    &
    773                                                local_pf, nzb_do, nzt_do )
    774                 resorted = .TRUE.
    775              ENDIF
    776 
    777772             IF ( .NOT. found )  THEN
    778773                CALL tcm_data_output_3d( av, do3d(av,if), found, local_pf,     &
    779774                                         nzb_do, nzt_do )
    780775                resorted = .TRUE.
     776             ENDIF
     777
     778             IF ( .NOT. found  .AND.  biometeorology )  THEN
     779                 CALL biom_data_output_3d( av, do3d(av,if), found, local_pf,   &
     780                                           nzb_do, nzt_do )
    781781             ENDIF
    782782

  • palm/trunk/SOURCE/header.f90

    r3355 r3448  
    2525! -----------------
    2626! $Id$
     27! Add biometeorology
     28!
     29! 3355 2018-10-16 14:03:34Z knoop
    2730! Header output concerning offline nesting
    2831!
    … …  
    238241! 1551 2015-03-03 14:18:16Z maronga
    239242! Added informal output for land surface model and radiation model. Removed typo.
     243!
     244! 1496 2014-12-02 17:25:50Z maronga
     245! Renamed: "radiation -> "cloud_top_radiation"
    240246!
    241247! 1484 2014-10-21 10:53:05Z kanani
    … …  
    395401        ONLY:  g, kappa, l_v, r_d
    396402
     403    USE biometeorology_mod,                                                    &
     404        ONLY:  biom_header
     405
    397406    USE bulk_cloud_model_mod,                                                  &
    398407        ONLY:  bulk_cloud_model, bcm_header
    … …  
    420429    USE gust_mod,                                                              &
    421430        ONLY: gust_header, gust_module_enabled
    422        
     431
    423432    USE indices,                                                               &
    424433        ONLY:  mg_loc_ind, nnx, nny, nnz, nx, ny, nxl_mg, nxr_mg, nyn_mg,      &
    … …  
    438447    USE netcdf_interface,                                                      &
    439448        ONLY:  netcdf_data_format, netcdf_data_format_string, netcdf_deflate
    440        
     449
    441450    USE nesting_offl_mod,                                                      &
    442451        ONLY:  nesting_offl_header
    … …  
    906915!-- Large scale forcing and nudging
    907916    IF ( large_scale_forcing )  CALL lsf_nudging_header( io )
    908    
     917
    909918!
    910919!-- Offline nesting
    … …  
    11331142!
    11341143!-- Header information from other modules
     1144    IF ( biometeorology      )  CALL biom_header ( io )
    11351145    IF ( gust_module_enabled )  CALL gust_header( io )
    11361146    IF ( land_surface        )  CALL lsm_header( io )

  • palm/trunk/SOURCE/init_3d_model.f90

    r3421 r3448  
    2525! -----------------
    2626! $Id$
     27! Add biometeorology
     28!
     29! 3421 2018-10-24 18:39:32Z gronemeier
    2730! Initialize surface data output
    2831!
    … …  
    533536               ideal_gas_law_rho, ideal_gas_law_rho_pt, barometric_formula
    534537
     538    USE biometeorology_mod,                                                    &
     539        ONLY:  biom_init, biom_init_arrays
     540
    535541    USE bulk_cloud_model_mod,                                                  &
    536542        ONLY:  bulk_cloud_model, bcm_init, bcm_init_arrays
    … …  
    552558    USE gust_mod,                                                              &
    553559        ONLY:  gust_init, gust_init_arrays, gust_module_enabled
    554    
     560
    555561    USE indices
    556562   
    … …  
    578584        ONLY:  chem_emis, chem_emis_att, init_3d,                              &
    579585               netcdf_data_input_init_3d, netcdf_data_input_interpolate
    580        
     586
    581587    USE nesting_offl_mod,                                                      &
    582588        ONLY:  nesting_offl_init
    … …  
    628634        ONLY :  init_surface_arrays, init_surfaces, surf_def_h, surf_lsm_h,    &
    629635                surf_usm_h, get_topography_top_index_ji, vertical_surfaces_exist
    630                 
     636   
    631637    USE surface_output_mod,                                                    &
    632638        ONLY:  surface_output_init
    … …  
    10501056!
    10511057!-- Allocate arrays for other modules
     1058    IF ( biometeorology      )  CALL biom_init_arrays
    10521059    IF ( bulk_cloud_model    )  CALL bcm_init_arrays
    10531060    IF ( gust_module_enabled )  CALL gust_init_arrays
    … …  
    24132420
    24142421!
     2422!-- If required initialize biometeorology module
     2423    IF ( biometeorology )  THEN
     2424        CALL location_message( 'initializing biometeorology module', .FALSE. )
     2425        CALL biom_init
     2426        CALL location_message( 'finished', .TRUE. )
     2427    ENDIF
     2428
     2429!
    24152430!-- Initialize the ws-scheme.   
    2416     IF ( ws_scheme_sca .OR. ws_scheme_mom )  CALL ws_init       
     2431    IF ( ws_scheme_sca .OR. ws_scheme_mom )  CALL ws_init
    24172432
    24182433!

  • palm/trunk/SOURCE/modules.f90

    r3435 r3448  
    2525! -----------------
    2626! $Id$
     27! Add biometeorology
     28!
     29! 3435 2018-10-26 18:25:44Z gronemeier
    2730! +mask_k_over_surface, mask_surface
    2831!
    … …  
    425428! Increased pr_palm to 120. Increased length of dots_unit and dots_label to 13
    426429! digits. Increased length of domask, do2d, and do3d to 20 digits.
     430!
     431! 1496 2014-12-02 17:25:50Z maronga
     432! Renamed "radiation" -> "cloud_top_radiation"
    427433!
    428434! 1484 2014-10-21 10:53:05Z kanani
    … …  
    13051311    LOGICAL ::  bc_radiation_r = .FALSE.                         !< radiation boundary condition for outflow at right domain boundary
    13061312    LOGICAL ::  bc_radiation_s = .FALSE.                         !< radiation boundary condition for outflow at south domain boundary
     1313    LOGICAL ::  biometeorology = .FALSE.                         !< biometeorology module switch
    13071314    LOGICAL ::  calc_soil_moisture_during_spinup = .FALSE.       !< namelist parameter
    13081315    LOGICAL ::  call_psolver_at_all_substeps = .TRUE.            !< namelist parameter

  • palm/trunk/SOURCE/multi_agent_system_mod.f90

    r3274 r3448  
    2525! -----------------
    2626! $Id$
     27! Adjustment of biometeorology calls,
     28! implement some agent biometeorology
     29!
     30! 3274 2018-09-24 15:42:55Z knoop
    2731! Modularization of all bulk cloud physics code components
    2832!
    … …  
    8488
    8589    USE control_parameters,                                                    &
    86         ONLY:  dt_3d, message_string, time_since_reference_point, dt_write_agent_data
     90        ONLY:  biometeorology, dt_3d, dt_write_agent_data, message_string,     &
     91               time_since_reference_point
    8792
    8893    USE cpulog,                                                                &
    … …  
    206211        REAL(wp)     ::  clo                  !< clothing index
    207212        REAL(wp)     ::  energy_storage       !< energy stored by agent
     213        REAL(wp)     ::  clothing_temp        !< energy stored by agent
     214        REAL(wp)     ::  actlev               !< metabolic + work energy of the person
     215        REAL(wp)     ::  age_years            !< physical age of the person
     216        REAL(wp)     ::  weight               !< total weight of the person (kg)
     217        REAL(wp)     ::  height               !< height of the person (m)
     218        REAL(wp)     ::  work                 !< workload of the agent (W)
     219        INTEGER(iwp) ::  sex                  !< agents gender: 1 = male, 2 = female
    208220        REAL(wp)     ::  force_x              !< force term x-direction
    209221        REAL(wp)     ::  force_y              !< force term y-direction
    … …  
    218230        REAL(wp)     ::  speed_x              !< speed of agent in x
    219231        REAL(wp)     ::  speed_y              !< speed of agent in y
    220         REAL(wp)     ::  thermal_index        !< the dynamic thermal index
     232        REAL(wp)     ::  ipt                  !< instationary thermal index iPT (°C)
    221233        REAL(wp)     ::  windspeed            !< absolute value of windspeed at agent position
    222234        REAL(wp)     ::  x                    !< x-position
    … …  
    329341 SUBROUTINE multi_agent_system
    330342
    331 !     USE htcm_mod,                                                              &
    332 !         ONLY:  htcm_dynamic
     343    USE biometeorology_mod,                                                   &
     344        ONLY:  biom_calc_ipt, biom_determine_input_at
     345
    333346
    334347    IMPLICIT NONE
    … …  
    341354    INTEGER(iwp)       ::  js                 !< counter
    342355    INTEGER(iwp), SAVE ::  mas_count = 0      !< counts the mas-calls
     356    INTEGER(iwp)                :: a     !< agent iterator
     357    !-- local meteorological conditions
     358    REAL(wp)                    :: tmrt  !< mean radiant temperature        (°C)
     359    REAL(wp)                    :: ta    !< air temperature                 (°C)
     360    REAL(wp)                    :: vp    !< vapour pressure                 (hPa)
     361    REAL(wp)                    :: v     !< wind speed    (local level)     (m/s)
     362    REAL(wp)                    :: pair  !< air pressure                    (hPa)
     363
    343364
    344365    LOGICAL       ::  first_loop_stride   !< flag for first loop stride of agent sub-timesteps
    … …  
    502523       deleted_agents = 0
    503524!
    504 !--    to be included here: call of human thermal comfort mod (and UV exposure)
    505 !        DO  i = nxl, nxr
    506 !           DO  j = nys, nyn
    507 !
    508 !              number_of_agents = agt_count(j,i)
    509 ! !
    510 ! !--          If grid cell gets empty, cycle
    511 !              IF ( number_of_agents <= 0 ) CYCLE
    512 !
    513 !              agents => grid_agents(j,i)%agents(1:number_of_agents)
    514 ! !
    515 ! !--          Evaluation of social forces
    516 !              CALL htcm_dynamic(i,j)
    517 !
    518 !           ENDDO
    519 !        ENDDO
     525       IF ( biometeorology )  THEN
     526!
     527!--       Call of human thermal comfort mod (and UV exposure)
     528          DO  i = nxl, nxr
     529             DO  j = nys, nyn
     530
     531                number_of_agents = agt_count(j,i)
     532!
     533!--             If grid cell gets empty, cycle
     534                IF ( number_of_agents <= 0 )  CYCLE
     535
     536                agents => grid_agents(j,i)%agents(1:number_of_agents)
     537!
     538!--             Evaluation of social forces
     539!                CALL biom_dynamic(i,j)
     540!
     541!--             Determine local meteorological conditions
     542                CALL biom_determine_input_at ( .FALSE., i, j, ta, vp, v,      &
     543                                               pair, tmrt )
     544
     545                DO  a = 1, number_of_agents
     546!
     547!--                Calculate instationary thermal indices based on local tmrt
     548
     549                   CALL biom_calc_ipt ( ta, vp, v, pair, tmrt,                &
     550                                        agents(a)%dt_sum,                     &
     551                                        agents(a)%energy_storage,             &
     552                                        agents(a)%clothing_temp,              &
     553                                        agents(a)%clo,                        &
     554                                        agents(a)%actlev,                     &
     555                                        agents(a)%age_years,                  &
     556                                        agents(a)%weight,                     &
     557                                        agents(a)%height,                     &
     558                                        agents(a)%work,                       &
     559                                        agents(a)%sex,                        &
     560                                        agents(a)%ipt )
     561                END DO
     562
     563             ENDDO
     564          ENDDO
     565       ENDIF
    520566
    521567       IF ( dt_3d_reached_mas )  EXIT
    … …  
    799845                         tmp_agent%age_m         = 0.0_wp
    800846                         tmp_agent%dt_sum        = 0.0_wp
    801                          tmp_agent%clo           = 99999.0_wp
     847                         tmp_agent%clo           = -999.0_wp
    802848                         tmp_agent%energy_storage= 0.0_wp
     849                         tmp_agent%ipt           = 99999.0_wp
     850                         tmp_agent%clothing_temp = -999._wp      !< energy stored by agent (W)
     851                         tmp_agent%actlev        = 134.6862_wp   !< metabolic + work energy of the person
     852                         tmp_agent%age_years     = 35._wp        !< physical age of the person
     853                         tmp_agent%weight        = 75._wp        !< total weight of the person (kg)
     854                         tmp_agent%height        = 1.75_wp       !< height of the person (m)
     855                         tmp_agent%work          = 134.6862_wp   !< workload of the agent (W)
     856                         tmp_agent%sex           = 1             !< agents gender: 1 = male, 2 = female
    803857                         tmp_agent%force_x       = 0.0_wp
    804858                         tmp_agent%force_y       = 0.0_wp
    … …  
    812866                         tmp_agent%speed_x       = 0.0_wp
    813867                         tmp_agent%speed_y       = 0.0_wp
    814                          tmp_agent%thermal_index = 99999.0_wp
    815868                         tmp_agent%x             = pos_x
    816869                         tmp_agent%y             = pos_y
    … …  
    31853238          zero_agent%speed_x       = 0.0_wp
    31863239          zero_agent%speed_y       = 0.0_wp
    3187           zero_agent%thermal_index = 0.0_wp
     3240          zero_agent%ipt          = 0.0_wp
    31883241          zero_agent%x             = 0.0_wp
    31893242          zero_agent%y             = 0.0_wp

  • palm/trunk/SOURCE/netcdf_interface_mod.f90

    r3435 r3448  
    2525! -----------------
    2626! $Id$
     27! Adjustment of biometeorology calls
     28!
     29! 3435 2018-10-26 18:25:44Z gronemeier
    2730! Bugfix: corrected order of calls to define_netcdf_grid for masked output
    2831! Add vertical dimensions to masked output in case of terrain-following output
    … …  
    312315
    313316    USE control_parameters,                                                    &
    314         ONLY:  fl_max, max_masks, multi_agent_system_end,                      &
     317        ONLY:  biometeorology, fl_max,                                         &
     318               max_masks, multi_agent_system_end,                              &
    315319               multi_agent_system_start, var_fl_max, varnamelength
    316320    USE kinds
    … …  
    326330          (/ 'ag_id           ', 'ag_x            ', 'ag_y            ',       &
    327331             'ag_wind         ', 'ag_temp         ', 'ag_group        ',       &
    328              'PM10            ', 'PM25            ', 'ag_therm_comf   ',       &
     332             'PM10            ', 'PM25            ', 'ag_iPT          ',       &
    329333             'ag_uv           ', 'not_used        ', 'not_used        ',       &
    330334             'not_used        ' /)
    … …  
    513517                   id_var_zu_mask, id_var_zw_mask, &
    514518                   id_var_zusi_mask, id_var_zwwi_mask
    515                    
     519
    516520    INTEGER(iwp), DIMENSION(1:max_masks,0:1,0:2) ::  id_var_eutm_mask, &
    517521                                                     id_var_nutm_mask
    … …  
    576580    USE arrays_3d,                                                             &
    577581        ONLY:  zu, zw
     582
     583    USE biometeorology_mod,                                                    &
     584        ONLY:  biom_define_netcdf_grid
    578585
    579586    USE chemistry_model_mod,                                                   &
    … …  
    635642    USE radiation_model_mod,                                                   &
    636643        ONLY: radiation, radiation_define_netcdf_grid
    637 
    638     USE biometeorology_mod,                                                    &
    639         ONLY: biometeorology_define_netcdf_grid
    640644
    641645    USE spectra_mod,                                                           &
    … …  
    11321136                   ENDIF
    11331137!
    1134 !--                Check for biometeorology quantities
    1135                    IF ( .NOT. found  .AND.  radiation )  THEN
    1136                       CALL biometeorology_define_netcdf_grid( domask(mid,av,i),&
    1137                                                          found, grid_x, grid_y,&
    1138                                                          grid_z )
    1139                    ENDIF
    1140 !
    11411138!--                Now check for user-defined quantities
    11421139                   IF ( .NOT. found )  THEN
    … …  
    18721869                   ENDIF
    18731870
    1874 !
    1875 !--                Check for biometeorology quantities
    1876                    IF ( .NOT. found  .AND.  radiation )  THEN
    1877                       CALL biometeorology_define_netcdf_grid( do3d(av,i), found,&
    1878                                                          grid_x, grid_y,       &
    1879                                                          grid_z )
    1880                    ENDIF
    1881 
    18821871!--                Check for gust module quantities
    18831872                   IF ( .NOT. found  .AND.  gust_module_enabled )  THEN
    18841873                      CALL gust_define_netcdf_grid( do3d(av,i), found, grid_x, &
    18851874                                                    grid_y, grid_z )
     1875                   ENDIF
     1876
     1877!
     1878!--                Check for biometeorology quantities
     1879                   IF ( .NOT. found  .AND.  biometeorology )  THEN
     1880                      CALL biom_define_netcdf_grid( do3d(av,i), found,         &
     1881                                                    grid_x, grid_y, grid_z )
    18861882                   ENDIF
    18871883
    … …  
    24692465
    24702466          ENDDO
     2467!
     2468!--       Define vars for biometeorology
     2469          CALL netcdf_create_var( id_set_agt, (/ id_dim_agtnum,                &
     2470                                  id_dim_time_agt /), agt_var_names(9),        &
     2471                                  nc_precision(8), id_var_agt(9),              &
     2472                                  TRIM( agt_var_units(9) ),                    &
     2473                                  TRIM( agt_var_names(9) ), 339, 340, 341 )
     2474
    24712475!
    24722476!--       Leave netCDF define mode
    … …  
    28272831
    28282832!
    2829 !--                      Check for biometeorology quantities
    2830                          IF ( .NOT. found  .AND.  radiation )  THEN
    2831                             CALL biometeorology_define_netcdf_grid( do2d(av,i),&
    2832                                                          found, grid_x, grid_y,&
    2833                                                          grid_z )
    2834                          ENDIF
    2835 
    2836 !
    28372833!--                      Check for gust module quantities
    28382834                         IF ( .NOT. found  .AND.  gust_module_enabled )  THEN
    … …  
    28412837                                                          grid_z )
    28422838                         ENDIF
    2843 
     2839!
     2840!--                      Check for human thermal comfort quantities
     2841                         IF ( .NOT. found  .AND.  biometeorology )  THEN
     2842                            CALL biom_define_netcdf_grid( do2d( av, i), found, &
     2843                                                          grid_x, grid_y,      &
     2844                                                          grid_z )
     2845                         ENDIF
    28442846!
    28452847!--                      Check for chemistry quantities
    … …  
    28942896                   ELSEIF ( grid_z == 'zs' )  THEN
    28952897                      id_z = id_dim_zs_xy(av)
     2898                   ELSEIF ( grid_z == 'zu1' )  THEN
     2899                      id_z = id_dim_zu1_xy(av)
    28962900                   ENDIF
    28972901
    … …  
    31883192                                          count = (/ nx+1, ny+1 /) )
    31893193                  CALL netcdf_handle_error( 'netcdf_define_header', 556 )
    3190                
     3194
    31913195                ENDDO
    31923196                DEALLOCATE( netcdf_data_2d )
    … …  
    37243728
    37253729!
    3726 !--                   Check for biometeorology quantities
    3727                       IF ( .NOT. found  .AND.  radiation )  THEN
    3728                          CALL biometeorology_define_netcdf_grid( do2d(av,i),   &
    3729                                                             found,             &
    3730                                                             grid_x, grid_y,    &
    3731                                                             grid_z )
    3732                       ENDIF
    3733 
    3734 !
    37353730!--                   Check for gust module quantities
    37363731                      IF ( .NOT. found  .AND.  gust_module_enabled )  THEN
    … …  
    45824577
    45834578!
    4584 !--                   Check for biometeorology quantities
    4585                       IF ( .NOT. found  .AND.  radiation )  THEN
    4586                          CALL biometeorology_define_netcdf_grid( do2d(av,i),   &
    4587                                                             found,             &
    4588                                                             grid_x, grid_y,    &
    4589                                                             grid_z )
    4590                       ENDIF
    4591 
    4592 !
    45934579!--                   Check for gust module quantities
    45944580                      IF ( .NOT. found  .AND.  gust_module_enabled )  THEN

  • palm/trunk/SOURCE/parin.f90

    r3435 r3448  
    2525! -----------------
    2626! $Id$
     27! Add biometeorology
     28!
     29! 3435 2018-10-26 18:25:44Z gronemeier
    2730! Add mask_k_over_surface
    2831!
    … …  
    304307! 1560 2015-03-06 10:48:54Z keck
    305308! +recycling_yshift
     309!
     310! 1496 2014-12-02 17:25:50Z maronga
     311! Renamed: "radiation -> "cloud_top_radiation"
    306312!
    307313! 1484 2014-10-21 10:53:05Z kanani
    … …  
    442448               ug, u_init, v_init, vg
    443449
     450    USE biometeorology_mod,                                                    &
     451        ONLY:  biom_parin
     452
    444453    USE bulk_cloud_model_mod,                                                  &
    445454        ONLY:  bcm_parin
    … …  
    494503
    495504    USE pegrid
    496                
     505
    497506    USE plant_canopy_model_mod,                                                &
    498507         ONLY: pcm_parin
    … …  
    702711             dt_do2d_xz, dt_do2d_yz, dt_do3d, dt_max, dt_restart,              &
    703712             dt_run_control,end_time, force_print_header, mask_k_over_surface, &
    704              mask_scale_x,        &
     713             mask_scale_x,                                                     &
    705714             mask_scale_y, mask_scale_z, mask_x, mask_y, mask_z, mask_x_loop,  &
    706715             mask_y_loop, mask_z_loop, netcdf_data_format, netcdf_deflate,     &
    … …  
    724733             dt_do2d_xz, dt_do2d_yz, dt_do3d, dt_max, dt_restart,              &
    725734             dt_run_control,end_time, force_print_header, mask_k_over_surface, &
    726              mask_scale_x,        &
     735             mask_scale_x,                                                     &
    727736             mask_scale_y, mask_scale_z, mask_x, mask_y, mask_z, mask_x_loop,  &
    728737             mask_y_loop, mask_z_loop, netcdf_data_format, netcdf_deflate,     &
    … …  
    884893!
    885894!--       Check for module namelists and read them
     895          CALL biom_parin
    886896          CALL lsm_parin
    887 !
    888 !--       Check if bulk cloud model is used and read namelist if required
    889           CALL bcm_parin
    890 !
    891 !--       Check if surface output is enabled and read
    892 !--       &surface_output_parameters if required         
     897          CALL bcm_parin
    893898          CALL surface_output_parin
    894899          CALL usm_parin

  • palm/trunk/SOURCE/sum_up_3d_data.f90

    r3421 r3448  
    2525! -----------------
    2626! $Id$
     27! Adjustment of biometeorology calls
     28!
     29! 3421 2018-10-24 18:39:32Z gronemeier
    2730! Renamed output variables
    2831!
    … …  
    236239        ONLY:  c_p, lv_d_cp, l_v
    237240
     241    USE biometeorology_mod,                                                    &
     242        ONLY:  biom_3d_data_averaging
     243
    238244    USE bulk_cloud_model_mod,                                                  &
    239245        ONLY:  bulk_cloud_model, bcm_3d_data_averaging
    … …  
    243249
    244250    USE control_parameters,                                                    &
    245         ONLY:  air_chemistry, average_count_3d, doav, doav_n, land_surface,    &
    246                ocean_mode, rho_surface, urban_surface, uv_exposure,            &
    247                varnamelength
     251        ONLY:  air_chemistry, average_count_3d, biometeorology, doav, doav_n,  &
     252               land_surface, ocean_mode, rho_surface, urban_surface,           &
     253               uv_exposure, varnamelength
    248254
    249255    USE cpulog,                                                                &
    … …  
    269275    USE radiation_model_mod,                                                   &
    270276        ONLY:  radiation, radiation_3d_data_averaging
    271 
    272     USE biometeorology_mod,                                                    &
    273         ONLY:  biometeorology_3d_data_averaging
    274277
    275278    USE surface_mod,                                                           &
    … …  
    517520                ENDIF
    518521
    519                 IF ( radiation  .AND.  trimvar(1:4) == 'bio_' )  THEN
    520                    CALL biometeorology_3d_data_averaging( 'allocate', doav(ii) )
    521                 ENDIF
    522 
    523522                IF ( gust_module_enabled )  THEN
    524523                   CALL gust_3d_data_averaging( 'allocate', doav(ii) )
     524                ENDIF
     525
     526                IF ( biometeorology  .AND.  trimvar(1:5) == 'biom_')  THEN
     527                   CALL biom_3d_data_averaging( 'allocate', doav(ii) )
    525528                ENDIF
    526529
    … …  
    11651168             ENDIF
    11661169
    1167              IF ( radiation   .AND.  trimvar(1:4) == 'bio_')  THEN
    1168                 CALL biometeorology_3d_data_averaging( 'sum', doav(ii) )
    1169              ENDIF
    1170 
    11711170             IF ( gust_module_enabled )  THEN
    11721171                CALL gust_3d_data_averaging( 'sum', doav(ii) )
     1172             ENDIF
     1173
     1174             IF ( biometeorology  .AND.  trimvar(1:5) == 'biom_' )  THEN
     1175                   CALL biom_3d_data_averaging( 'sum', doav(ii) )
    11731176             ENDIF
    11741177

  • palm/trunk/SOURCE/time_integration.f90

    r3421 r3448  
    2525! -----------------
    2626! $Id$
     27! Add biometeorology
     28!
     29! 3421 2018-10-24 18:39:32Z gronemeier
    2730! Surface data output
    2831!
    … …  
    393396               v_p, w, w_p
    394397
     398    USE biometeorology_mod,                                                    &
     399        ONLY:  biom_calculate_static_grid, time_biom_results
     400
    395401    USE bulk_cloud_model_mod,                                                  &
    396402        ONLY: bulk_cloud_model, calc_liquid_water_content,                     &
    … …  
    413419        ONLY:  advected_distance_x, advected_distance_y, air_chemistry,        &
    414420               average_count_3d, averaging_interval, averaging_interval_pr,    &
    415                bc_lr_cyc, bc_ns_cyc, bc_pt_t_val, bc_q_t_val,                  &
     421               bc_lr_cyc, bc_ns_cyc, bc_pt_t_val, bc_q_t_val, biometeorology,  &
    416422               call_psolver_at_all_substeps,  child_domain, cloud_droplets,    &
    417423               constant_flux_layer, constant_heatflux,                         &
    … …  
    479485    USE multi_agent_system_mod,                                                &
    480486        ONLY:  agents_active, multi_agent_system
    481        
     487
    482488    USE nesting_offl_mod,                                                      &
    483489        ONLY:  nesting_offl_bc, nesting_offl_mass_conservation
    … …  
    523529    USE surface_mod,                                                           &
    524530        ONLY:  surf_def_h, surf_lsm_h, surf_usm_h
    525        
     531
    526532    USE surface_output_mod,                                                    &
    527533        ONLY:  average_count_surf, averaging_interval_surf, dt_dosurf,         &
    … …  
    12851291
    12861292!
     1293!--    Biometeorology calculation of stationary thermal indices
     1294       IF ( biometeorology  .AND.  time_do3d >= dt_do3d )  THEN
     1295          CALL biom_calculate_static_grid ( .FALSE. )
     1296          time_biom_results = time_since_reference_point
     1297       ENDIF
     1298
     1299!
    12871300!--    Execute the gust module actions
    12881301       IF ( gust_module_enabled )  THEN
    … …  
    15481561
    15491562    ENDDO   ! time loop
    1550    
     1563
    15511564!
    15521565!-- Vertical nesting: Deallocate variables initialized for vertical nesting   

  • palm/trunk/SOURCE/uv_exposure_model_mod.f90

    r3274 r3448  
    2525! -----------------
    2626! $Id$
     27! Temporary rename of namelist until uv model moves to biometeorology module
     28!
     29! 3274 2018-09-24 15:42:55Z knoop
    2730! Modularization of all bulk cloud physics code components
    2831!
    … …  
    454457    CHARACTER (LEN=80) ::  line  !< dummy string for current line in parameter file
    455458   
    456     NAMELIST /biometeorology_parameters/  clothing
     459    NAMELIST /biometeorology_parameters_uv/  clothing
    457460   
    458461    line = ' '
    … …  
    462465    REWIND ( 11 )
    463466    line = ' '
    464     DO WHILE ( INDEX( line, '&biometeorology_parameters' ) == 0 )
     467    DO WHILE ( INDEX( line, '&biometeorology_parameters_uv' ) == 0 )
    465468       READ ( 11, '(A)', END=20 )  line
    466469    ENDDO
    … …  
    469472!
    470473!-- Read user-defined namelist
    471     READ ( 11, biometeorology_parameters, ERR = 10, END = 20 )
     474    READ ( 11, biometeorology_parameters_uv, ERR = 10, END = 20 )
    472475
    473476!
    … …  
    477480 10 BACKSPACE( 11 )
    478481    READ( 11 , '(A)') line
    479     CALL parin_fail_message( 'biometeorology_parameters', line )
     482    CALL parin_fail_message( 'biometeorology_parameters_uv', line )
    480483
    481484
Note: See TracChangeset for help on using the changeset viewer.
                                                                                                                                                                                                                                                                                                                                                                               
  | Impressum | ©Leibniz Universität Hannover |