Changeset 3569 for palm/trunk
- Timestamp:
- Nov 27, 2018 5:03:40 PM (6 years ago)
- Location:
- palm/trunk/SOURCE
- Files:
-
- 1 deleted
- 11 edited
Legend:
- Unmodified
- Added
- Removed
-
palm/trunk/SOURCE/Makefile
r3525 r3569 25 25 # ----------------- 26 26 # $Id$ 27 # dom_dwd_user, Schrempf: 28 # New dependencies for biometeorology_mod, 29 # remove uv exposure model code, this is now part of biometeorology_mod. 30 # 31 # 3525 2018-11-14 16:06:14Z kanani 27 32 # Changes related to clean-up of biometeorology (dom_dwd_user) 28 33 # … … 701 706 user_statistics.f90 \ 702 707 user_write_restart_data_mod.f90 \ 703 uv_exposure_model_mod.f90 \704 708 vertical_nesting_mod.f90 \ 705 709 virtual_flight_mod.f90 \ … … 791 795 mod_kinds.o 792 796 biometeorology_mod.o: \ 793 modules.o \ 794 mod_kinds.o \ 797 basic_constants_and_equations_mod.o \ 798 date_and_time_mod.o \ 799 modules.o \ 800 mod_kinds.o \ 801 netcdf_data_input_mod.o \ 795 802 radiation_model_mod.o 796 803 boundary_conds.o: \ … … 853 860 turbulence_closure_mod.o \ 854 861 urban_surface_mod.o \ 855 uv_exposure_model_mod.o \856 862 vertical_nesting_mod.o \ 857 863 wind_turbine_model_mod.o … … 974 980 salsa_mod.o \ 975 981 surface_mod.o \ 976 urban_surface_mod.o \ 977 uv_exposure_model_mod.o 982 urban_surface_mod.o 978 983 data_output_3d.o: \ 979 984 basic_constants_and_equations_mod.o \ … … 1128 1133 turbulence_closure_mod.o \ 1129 1134 urban_surface_mod.o \ 1130 uv_exposure_model_mod.o \1131 1135 virtual_flight_mod.o \ 1132 1136 virtual_measurement_mod.o \ … … 1373 1377 spectra_mod.o \ 1374 1378 turbulence_closure_mod.o \ 1375 urban_surface_mod.o \ 1376 uv_exposure_model_mod.o 1379 urban_surface_mod.o 1377 1380 nesting_offl_mod.o: \ 1378 1381 cpulog_mod.o \ … … 1446 1449 synthetic_turbulence_generator_mod.o \ 1447 1450 turbulence_closure_mod.o \ 1448 uv_exposure_model_mod.o \1449 1451 vertical_nesting_mod.o \ 1450 1452 virtual_flight_mod.o \ … … 1658 1660 surface_mod.o \ 1659 1661 turbulence_closure_mod.o \ 1660 urban_surface_mod.o \ 1661 uv_exposure_model_mod.o 1662 urban_surface_mod.o 1662 1663 surface_coupler.o: \ 1663 1664 basic_constants_and_equations_mod.o \ … … 1746 1747 urban_surface_mod.o \ 1747 1748 user_actions.o \ 1748 uv_exposure_model_mod.o \1749 1749 vertical_nesting_mod.o \ 1750 1750 virtual_flight_mod.o \ … … 1948 1948 user_write_restart_data_mod.o: \ 1949 1949 user_module.o 1950 uv_exposure_model_mod.o: \1951 basic_constants_and_equations_mod.o \1952 date_and_time_mod.o \1953 mod_kinds.o \1954 modules.o \1955 netcdf_data_input_mod.o \1956 radiation_model_mod.o1957 1950 vertical_nesting_mod.o: \ 1958 1951 mod_kinds.o \ -
palm/trunk/SOURCE/average_3d_data.f90
r3525 r3569 25 25 ! ----------------- 26 26 ! $Id$ 27 ! dom_dwd_user: 28 ! Clean up biometeorology call 29 ! 30 ! 3525 2018-11-14 16:06:14Z kanani 27 31 ! Changes related to clean-up of biometeorology (dom_dwd_user) 28 32 ! … … 583 587 ENDIF 584 588 585 IF ( biometeorology .AND. trimvar(1:4) == 'bio_') THEN589 IF ( biometeorology ) THEN 586 590 CALL bio_3d_data_averaging( 'average', doav(ii) ) 587 591 ENDIF -
palm/trunk/SOURCE/biometeorology_mod.f90
r3525 r3569 21 21 ! 22 22 ! Current revisions: 23 ! ----------------- 23 ! ------------------ 24 24 ! 25 25 ! … … 27 27 ! ----------------- 28 28 ! $Id$ 29 ! Consistently use bio_fill_value everywhere, 30 ! move allocation and initialization of output variables to bio_check_data_output 31 ! and bio_3d_data_averaging, 32 ! dom_dwd_user, Schrempf: 33 ! - integration of UVEM module part from r3474 (edited) 34 ! - split UTCI regression into 6 parts 35 ! - all data_output_3d is now explicity casted to sp 36 ! 37 ! 3525 2018-11-14 16:06:14Z kanani 29 38 ! Clean up, renaming from "biom" to "bio", summary of thermal index calculation 30 39 ! into one module (dom_dwd_user) … … 49 58 ! Authors: 50 59 ! -------- 51 ! @author Dominik Froehlich <dominik.froehlich@dwd.de> 52 ! @author Jaroslav Resler <resler@cs.cas.cz> 60 ! @author Dominik Froehlich <dominik.froehlich@dwd.de>, thermal indices 61 ! @author Jaroslav Resler <resler@cs.cas.cz>, mean radiant temperature 62 ! @author Michael Schrempf <schrempf@muk.uni-hannover.de>, uv exposure 53 63 ! 54 64 ! 55 65 ! Description: 56 66 ! ------------ 57 !> Human thermal comfort module calculating thermal perception of a sample 67 !> Biometeorology module consisting of two parts: 68 !> 1.: Human thermal comfort module calculating thermal perception of a sample 58 69 !> human being under the current meteorological conditions. 70 !> 2.: Calculation of vitamin-D weighted UV exposure 59 71 !> 60 72 !> @todo Alphabetical sorting of "USE ..." lists, "ONLY" list, variable declarations 61 73 !> (per subroutine: first all CHARACTERs, then INTEGERs, LOGICALs, REALs, ) 62 74 !> @todo Comments start with capital letter --> "!-- Include..." 63 !> @todo Variable and routine names strictly in lowercase letters and in English 75 !> @todo uv_vitd3dose-->new output type necessary (cumulative) 76 !> @todo consider upwelling radiation in UV 64 77 !> 65 78 !> @note nothing now … … 76 89 77 90 USE basic_constants_and_equations_mod, & 78 ONLY: c_p, degc_to_k, l_v, magnus, sigma_sb 91 ONLY: c_p, degc_to_k, l_v, magnus, sigma_sb, pi 79 92 80 93 USE control_parameters, & … … 83 96 simulated_time, surface_pressure 84 97 98 USE date_and_time_mod, & 99 ONLY: calc_date_and_time, day_of_year, time_utc 100 85 101 USE grid_variables, & 86 102 ONLY: ddx, dx, ddy, dy … … 92 108 USE kinds !< Set precision of INTEGER and REAL arrays according to PALM 93 109 110 USE netcdf_data_input_mod, & 111 ONLY: netcdf_data_input_uvem, uvem_projarea_f, uvem_radiance_f, & 112 uvem_irradiance_f, uvem_integration_f, building_obstruction_f 113 ! 94 114 !-- Import radiation model to obtain input for mean radiant temperature 95 115 USE radiation_model_mod, & … … 106 126 PRIVATE 107 127 128 ! 108 129 !-- Declare all global variables within the module (alphabetical order) 109 130 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: tmrt_grid !< tmrt results (°C) 110 131 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: perct !< PT results (°C) 111 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: utci_grid !< UTCI results (°C) 112 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: pet_grid !< PET results (°C) 132 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: utci !< UTCI results (°C) 133 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: pet !< PET results (°C) 134 ! 113 135 !-- Grids for averaged thermal indices 114 136 REAL(wp), DIMENSION(:), ALLOCATABLE :: mrt_av_grid !< time average mean 115 137 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: perct_av !< PT results (aver. input) (°C) 116 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: utci_av _grid!< UTCI results (aver. input) (°C)117 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: pet_av _grid!< PET results (aver. input) (°C)138 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: utci_av !< UTCI results (aver. input) (°C) 139 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: pet_av !< PET results (aver. input) (°C) 118 140 119 141 … … 123 145 REAL ( wp ), PARAMETER :: human_absorb = 0.7_wp !< SW absorbtivity of a human body (Fanger 1972) 124 146 REAL ( wp ), PARAMETER :: human_emiss = 0.97_wp !< LW emissivity of a human body after (Fanger 1972) 147 REAL ( wp ), PARAMETER :: bio_fill_value = -9999._wp !< set module fill value, replace by global fill value as soon as available 148 ! 125 149 !-- 126 127 150 LOGICAL :: aver_perct = .FALSE. !< switch: do perct averaging in this module? (if .FALSE. this is done globally) 128 151 LOGICAL :: aver_q = .FALSE. !< switch: do e averaging in this module? … … 135 158 LOGICAL :: average_trigger_pet = .FALSE. !< update averaged input on call to bio_pet? 136 159 160 LOGICAL :: thermal_comfort = .TRUE. !< Turn all thermal indices on or off 137 161 LOGICAL :: bio_perct = .TRUE. !< Turn index PT (instant. input) on or off 138 162 LOGICAL :: bio_perct_av = .TRUE. !< Turn index PT (averaged input) on or off … … 142 166 LOGICAL :: bio_utci_av = .TRUE. !< Turn index UTCI (averaged input) on or off 143 167 168 ! 169 !-- UVEM parameters from here 170 ! 171 !-- Declare all global variables within the module (alphabetical order) 172 INTEGER(iwp) :: ai = 0 !< loop index in azimuth direction 173 INTEGER(iwp) :: bi = 0 !< loop index of bit location within an 8bit-integer (one Byte) 174 INTEGER(iwp) :: clothing = 1 !< clothing (0=unclothed, 1=Arms,Hands,Face free, 3=Hand,Face free) 175 INTEGER(iwp) :: iq = 0 !< loop index of irradiance quantity 176 INTEGER(iwp) :: pobi = 0 !< loop index of the position of corresponding byte within ibset byte vektor 177 INTEGER(iwp) :: obstruction_direct_beam = 0 !< Obstruction information for direct beam 178 INTEGER(iwp) :: zi = 0 !< loop index in zenith direction 179 180 INTEGER(KIND=1), DIMENSION(0:44) :: obstruction_temp1 = 0 !< temporary obstruction information stored with ibset 181 INTEGER(iwp), DIMENSION(0:359) :: obstruction_temp2 = 0 !< restored temporary obstruction information from ibset file 182 183 INTEGER(iwp), DIMENSION(0:35,0:9) :: obstruction = 1 !< final 2D obstruction information array 184 185 LOGICAL :: consider_obstructions = .TRUE. !< namelist parameter (see documentation) 186 LOGICAL :: sun_in_south = .FALSE. !< namelist parameter (see documentation) 187 LOGICAL :: turn_to_sun = .TRUE. !< namelist parameter (see documentation) 188 LOGICAL :: uv_exposure = .FALSE. !< namelist parameter (see documentation) 189 190 REAL(wp) :: diffuse_exposure = 0.0_wp !< calculated exposure by diffuse radiation 191 REAL(wp) :: direct_exposure = 0.0_wp !< calculated exposure by direct solar beam 192 REAL(wp) :: orientation_angle = 0.0_wp !< orientation of front/face of the human model 193 REAL(wp) :: projection_area_direct_beam = 0.0_wp !< projection area for direct solar beam 194 REAL(wp) :: saa = 180.0_wp !< solar azimuth angle 195 REAL(wp) :: startpos_human = 0.0_wp !< start value for azimuth interpolation of human geometry array 196 REAL(wp) :: startpos_saa_float = 0.0_wp !< start value for azimuth interpolation of radiance array 197 REAL(wp) :: sza = 20.0_wp !< solar zenith angle 198 REAL(wp) :: xfactor = 0.0_wp !< relative x-position used for interpolation 199 REAL(wp) :: yfactor = 0.0_wp !< relative y-position used for interpolation 200 201 REAL(wp), DIMENSION(0:2) :: irradiance = 0.0_wp !< iradiance values extracted from irradiance lookup table 202 203 REAL(wp), DIMENSION(0:2,0:90) :: irradiance_lookup_table = 0.0_wp !< irradiance lookup table 204 REAL(wp), DIMENSION(0:35,0:9) :: integration_array = 0.0_wp !< solid angle factors for hemispherical integration 205 REAL(wp), DIMENSION(0:35,0:9) :: projection_area = 0.0_wp !< projection areas of a human (all directions) 206 REAL(wp), DIMENSION(0:35,0:9) :: projection_area_lookup_table = 0.0_wp !< human geometry lookup table (projection areas) 207 REAL(wp), DIMENSION(0:71,0:9) :: projection_area_direct_temp = 0.0_wp !< temporary projection area for direct solar beam 208 REAL(wp), DIMENSION(0:71,0:9) :: projection_area_temp = 0.0_wp !< temporary projection area for all directions 209 REAL(wp), DIMENSION(0:35,0:9) :: radiance_array = 0.0_wp !< radiance extracted from radiance_lookup_table 210 REAL(wp), DIMENSION(0:71,0:9) :: radiance_array_temp = 0.0_wp !< temporary radiance data 211 212 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: vitd3_exposure !< result variable for instantaneous vitamin-D weighted exposures 213 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: vitd3_exposure_av !< result variable for summation of vitamin-D weighted exposures 214 215 REAL(wp), DIMENSION(0:35,0:9,0:90) :: radiance_lookup_table = 0.0_wp !< radiance lookup table 144 216 145 217 ! … … 150 222 bio_check_parameters, bio_data_output_3d, bio_data_output_2d, & 151 223 bio_define_netcdf_grid, bio_get_thermal_index_input_ij, bio_header, & 152 bio_init, bio_init_arrays, bio_parin, bio_perct, bio_perct_av, bio_pet, & 153 bio_pet_av, bio_utci, bio_utci_av, time_bio_results 224 bio_init, bio_parin, bio_perct, bio_perct_av, bio_pet, & 225 bio_pet_av, bio_utci, bio_utci_av, thermal_comfort, time_bio_results, & 226 ! 227 !-- UVEM PUBLIC variables and methods 228 uvem_calc_exposure, uv_exposure 154 229 155 230 ! … … 160 235 MODULE PROCEDURE bio_3d_data_averaging 161 236 END INTERFACE bio_3d_data_averaging 162 237 ! 163 238 !-- Calculate mtr from rtm fluxes and assign into 2D grid 164 239 INTERFACE bio_calculate_mrt_grid 165 240 MODULE PROCEDURE bio_calculate_mrt_grid 166 241 END INTERFACE bio_calculate_mrt_grid 167 242 ! 168 243 !-- Calculate static thermal indices PT, UTCI and/or PET 169 244 INTERFACE bio_calculate_thermal_index_maps 170 245 MODULE PROCEDURE bio_calculate_thermal_index_maps 171 246 END INTERFACE bio_calculate_thermal_index_maps 172 247 ! 173 248 !-- Calculate the dynamic index iPT (to be caled by the agent model) 174 249 INTERFACE bio_calc_ipt 175 250 MODULE PROCEDURE bio_calc_ipt 176 251 END INTERFACE bio_calc_ipt 177 252 ! 178 253 !-- Data output checks for 2D/3D data to be done in check_parameters 179 180 181 182 254 INTERFACE bio_check_data_output 255 MODULE PROCEDURE bio_check_data_output 256 END INTERFACE bio_check_data_output 257 ! 183 258 !-- Input parameter checks to be done in check_parameters 184 259 INTERFACE bio_check_parameters 185 260 MODULE PROCEDURE bio_check_parameters 186 261 END INTERFACE bio_check_parameters 187 262 ! 188 263 !-- Data output of 2D quantities 189 264 INTERFACE bio_data_output_2d 190 265 MODULE PROCEDURE bio_data_output_2d 191 266 END INTERFACE bio_data_output_2d 192 267 ! 193 268 !-- no 3D data, thus, no averaging of 3D data, removed 194 269 INTERFACE bio_data_output_3d 195 270 MODULE PROCEDURE bio_data_output_3d 196 271 END INTERFACE bio_data_output_3d 197 272 ! 198 273 !-- Definition of data output quantities 199 274 INTERFACE bio_define_netcdf_grid 200 275 MODULE PROCEDURE bio_define_netcdf_grid 201 276 END INTERFACE bio_define_netcdf_grid 202 277 ! 203 278 !-- Obtains all relevant input values to estimate local thermal comfort/stress 204 279 INTERFACE bio_get_thermal_index_input_ij 205 280 MODULE PROCEDURE bio_get_thermal_index_input_ij 206 281 END INTERFACE bio_get_thermal_index_input_ij 207 282 ! 208 283 !-- Output of information to the header file 209 284 INTERFACE bio_header 210 285 MODULE PROCEDURE bio_header 211 286 END INTERFACE bio_header 212 287 ! 213 288 !-- Initialization actions 214 289 INTERFACE bio_init 215 290 MODULE PROCEDURE bio_init 216 291 END INTERFACE bio_init 217 218 !-- Initialization of arrays 219 INTERFACE bio_init_arrays 220 MODULE PROCEDURE bio_init_arrays 221 END INTERFACE bio_init_arrays 222 292 ! 223 293 !-- Reading of NAMELIST parameters 224 294 INTERFACE bio_parin … … 226 296 END INTERFACE bio_parin 227 297 298 ! 299 !-- Calculate UV exposure grid 300 INTERFACE uvem_calc_exposure 301 MODULE PROCEDURE uvem_calc_exposure 302 END INTERFACE uvem_calc_exposure 228 303 229 304 CONTAINS … … 235 310 !> Sum up and time-average biom input quantities as well as allocate 236 311 !> the array necessary for storing the average. 312 !> There is a considerable difference to the 3d_data_averaging subroutines 313 !> used by other modules: 314 !> For the thermal indices, the module needs to average the input conditions 315 !> not the result! 237 316 !------------------------------------------------------------------------------! 238 317 SUBROUTINE bio_3d_data_averaging( mode, variable ) … … 260 339 CASE ( 'bio_perct*', 'bio_utci*', 'bio_pet*' ) 261 340 262 !-- Indices in unknown order as depending on input file, determine 263 ! first index to average und update only once 264 IF ( .NOT. average_trigger_perct .AND. .NOT. average_trigger_utci & 265 .AND. .NOT. average_trigger_pet ) THEN 341 ! 342 !-- Averaging, as well as the allocation of the required grids must be 343 ! done only once, independent from for how many thermal indices 344 ! averaged output is desired. 345 ! Therefore wee need to memorize which index is the one that controls 346 ! the averaging (what must be the first thermal index called). 347 ! Indices are in unknown order as depending on the input file, 348 ! determine first index to average und update only once 349 350 !-- Only proceed here if this was not done for any index before. This 351 ! is done only once during the whole model run. 352 IF ( .NOT. average_trigger_perct .AND. & 353 .NOT. average_trigger_utci .AND. & 354 .NOT. average_trigger_pet ) THEN 355 ! 356 !-- Allocate the required grids 357 IF ( .NOT. ALLOCATED( perct_av ) ) THEN 358 ALLOCATE( perct_av (nys:nyn,nxl:nxr) ) 359 ENDIF 360 perct_av = REAL( bio_fill_value, KIND = wp ) 361 362 IF ( .NOT. ALLOCATED( utci_av ) ) THEN 363 ALLOCATE( utci_av (nys:nyn,nxl:nxr) ) 364 ENDIF 365 utci_av = REAL( bio_fill_value, KIND = wp ) 366 367 IF ( .NOT. ALLOCATED( pet_av ) ) THEN 368 ALLOCATE( pet_av (nys:nyn,nxl:nxr) ) 369 ENDIF 370 pet_av = REAL( bio_fill_value, KIND = wp ) 371 ! 372 !-- Memorize the first index called to control averaging 266 373 IF ( TRIM( variable ) == 'bio_perct*' ) THEN 267 374 average_trigger_perct = .TRUE. … … 274 381 ENDIF 275 382 ENDIF 276 277 !-- Only continue if updateindex 383 ! 384 !-- Only continue if var is the index, that controls averaging. 385 ! Break immediatelly (doing nothing) for the other indices. 278 386 IF ( average_trigger_perct .AND. TRIM( variable ) /= 'bio_perct*') RETURN 279 387 IF ( average_trigger_utci .AND. TRIM( variable ) /= 'bio_utci*') RETURN 280 388 IF ( average_trigger_pet .AND. TRIM( variable ) /= 'bio_pet*') RETURN 281 282 !-- Set averaging switch to .TRUE. if not done by other module before! 389 ! 390 !-- Now memorize which of the input grids are not averaged by other 391 ! modules. Set averaging switch to .TRUE. in that case. 283 392 IF ( .NOT. ALLOCATED( pt_av ) ) THEN 284 393 ALLOCATE( pt_av(nzb:nzt+1,nys:nyn,nxl:nxr) ) 285 394 aver_perct = .TRUE. 395 pt_av = 0.0_wp 286 396 ENDIF 287 pt_av = 0.0_wp288 397 289 398 IF ( .NOT. ALLOCATED( q_av ) ) THEN 290 399 ALLOCATE( q_av(nzb:nzt+1,nys:nyn,nxl:nxr) ) 291 400 aver_q = .TRUE. 401 q_av = 0.0_wp 292 402 ENDIF 293 q_av = 0.0_wp294 403 295 404 IF ( .NOT. ALLOCATED( u_av ) ) THEN 296 405 ALLOCATE( u_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) ) 297 406 aver_u = .TRUE. 407 u_av = 0.0_wp 298 408 ENDIF 299 u_av = 0.0_wp300 409 301 410 IF ( .NOT. ALLOCATED( v_av ) ) THEN 302 411 ALLOCATE( v_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) ) 303 412 aver_v = .TRUE. 413 v_av = 0.0_wp 304 414 ENDIF 305 v_av = 0.0_wp306 415 307 416 IF ( .NOT. ALLOCATED( w_av ) ) THEN 308 417 ALLOCATE( w_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) ) 309 418 aver_w = .TRUE. 419 w_av = 0.0_wp 310 420 ENDIF 311 w_av = 0.0_wp 421 422 CASE ( 'uvem_vitd3dose*' ) 423 IF ( .NOT. ALLOCATED( vitd3_exposure_av ) ) THEN 424 ALLOCATE( vitd3_exposure_av(nysg:nyng,nxlg:nxrg) ) 425 ENDIF 426 vitd3_exposure_av = 0.0_wp 312 427 313 428 CASE DEFAULT … … 335 450 336 451 CASE ( 'bio_perct*', 'bio_utci*', 'bio_pet*' ) 337 338 !-- Only continue if updateindex 452 ! 453 !-- Only continue if updateindex, see above 339 454 IF ( average_trigger_perct .AND. TRIM( variable ) /= 'bio_perct*') & 340 455 RETURN … … 393 508 ENDDO 394 509 ENDIF 395 396 CASE DEFAULT 510 ! 511 !-- This is a cumulated dose. No mode == 'average' for this quantity. 512 CASE ( 'uvem_vitd3dose*' ) 513 IF ( ALLOCATED( vitd3_exposure_av ) ) THEN 514 DO i = nxlg, nxrg 515 DO j = nysg, nyng 516 vitd3_exposure_av(j,i) = vitd3_exposure_av(j,i) + vitd3_exposure(j,i) 517 ENDDO 518 ENDDO 519 ENDIF 520 521 CASE DEFAULT 397 522 CONTINUE 398 523 … … 409 534 410 535 CASE ( 'bio_perct*', 'bio_utci*', 'bio_pet*' ) 411 412 !-- Only continue if update index 413 IF ( average_trigger_perct .AND. TRIM( variable ) /= 'bio_perct*')&414 RETURN415 IF ( average_trigger_utci .AND. TRIM( variable ) /= 'bio_utci*')&416 RETURN417 IF ( average_trigger_pet .AND. TRIM( variable ) /= 'bio_pet*')&418 RETURN536 ! 537 !-- Only continue if update index, see above 538 IF ( average_trigger_perct .AND. & 539 TRIM( variable ) /= 'bio_perct*' ) RETURN 540 IF ( average_trigger_utci .AND. & 541 TRIM( variable ) /= 'bio_utci*' ) RETURN 542 IF ( average_trigger_pet .AND. & 543 TRIM( variable ) /= 'bio_pet*' ) RETURN 419 544 420 545 IF ( ALLOCATED( pt_av ) .AND. aver_perct ) THEN … … 422 547 DO j = nys, nyn 423 548 DO k = nzb, nzt+1 424 pt_av(k,j,i) = pt_av(k,j,i) / REAL( average_count_3d, KIND=wp ) 549 pt_av(k,j,i) = pt_av(k,j,i) / & 550 REAL( average_count_3d, KIND=wp ) 425 551 ENDDO 426 552 ENDDO … … 432 558 DO j = nys, nyn 433 559 DO k = nzb, nzt+1 434 q_av(k,j,i) = q_av(k,j,i) / REAL( average_count_3d, KIND=wp ) 560 q_av(k,j,i) = q_av(k,j,i) / & 561 REAL( average_count_3d, KIND=wp ) 435 562 ENDDO 436 563 ENDDO … … 442 569 DO j = nysg, nyng 443 570 DO k = nzb, nzt+1 444 u_av(k,j,i) = u_av(k,j,i) / REAL( average_count_3d, KIND=wp ) 571 u_av(k,j,i) = u_av(k,j,i) / & 572 REAL( average_count_3d, KIND=wp ) 445 573 ENDDO 446 574 ENDDO … … 452 580 DO j = nysg, nyng 453 581 DO k = nzb, nzt+1 454 v_av(k,j,i) = v_av(k,j,i) / REAL( average_count_3d, KIND=wp ) 582 v_av(k,j,i) = v_av(k,j,i) / & 583 REAL( average_count_3d, KIND=wp ) 455 584 ENDDO 456 585 ENDDO … … 462 591 DO j = nysg, nyng 463 592 DO k = nzb, nzt+1 464 w_av(k,j,i) = w_av(k,j,i) / REAL( average_count_3d, KIND=wp ) 593 w_av(k,j,i) = w_av(k,j,i) / & 594 REAL( average_count_3d, KIND=wp ) 465 595 ENDDO 466 596 ENDDO 467 597 ENDDO 468 598 ENDIF 469 470 !-- Udate thermal indices with derived averages599 ! 600 !-- Udate all thermal index grids with updated averaged input 471 601 CALL bio_calculate_thermal_index_maps ( .TRUE. ) 602 603 ! 604 !-- No averaging for UVEM since we are calculating a dose (only sum is 605 !-- calculated and saved to av.nc file) 472 606 473 607 END SELECT … … 485 619 !> Check data output for biometeorology model 486 620 !------------------------------------------------------------------------------! 487 SUBROUTINE bio_check_data_output( var, unit )621 SUBROUTINE bio_check_data_output( var, unit, i, ilen, k ) 488 622 489 623 USE control_parameters, & … … 495 629 CHARACTER (LEN=*) :: var !< The variable in question 496 630 631 INTEGER(iwp) :: i !< 632 INTEGER(iwp) :: ilen !< 633 INTEGER(iwp) :: k !< 497 634 498 635 SELECT CASE ( TRIM( var ) ) 499 500 CASE( 'bio_mrt', 'bio_pet*', 'bio_perct*', 'bio_utci*' ) 636 ! 637 !-- Allocate a temporary array with the desired output dimensions. 638 CASE ( 'bio_mrt') 501 639 unit = 'degree_C' 640 IF ( .NOT. ALLOCATED( tmrt_grid ) ) THEN 641 ALLOCATE( tmrt_grid (nys:nyn,nxl:nxr) ) 642 ENDIF 643 tmrt_grid = REAL( bio_fill_value, KIND = wp ) 644 645 CASE ( 'bio_perct*' ) 646 unit = 'degree_C' 647 IF ( .NOT. ALLOCATED( perct ) ) THEN 648 ALLOCATE( perct (nys:nyn,nxl:nxr) ) 649 ENDIF 650 perct = REAL( bio_fill_value, KIND = wp ) 651 652 CASE ( 'bio_utci*' ) 653 unit = 'degree_C' 654 IF ( .NOT. ALLOCATED( utci ) ) THEN 655 ALLOCATE( utci (nys:nyn,nxl:nxr) ) 656 ENDIF 657 utci = REAL( bio_fill_value, KIND = wp ) 658 659 CASE ( 'bio_pet*' ) 660 unit = 'degree_C' 661 IF ( .NOT. ALLOCATED( pet ) ) THEN 662 ALLOCATE( pet (nys:nyn,nxl:nxr) ) 663 ENDIF 664 pet = REAL( bio_fill_value, KIND = wp ) 665 666 CASE ( 'uvem_vitd3*' ) 667 IF ( .NOT. uv_exposure ) THEN 668 message_string = 'output of "' // TRIM( var ) // '" requi' // & 669 'res a namelist &uvexposure_par' 670 CALL message( 'uvem_check_data_output', 'UV0001', 1, 2, 0, 6, 0 ) 671 ENDIF 672 IF ( k == 0 .OR. data_output(i)(ilen-2:ilen) /= '_xy' ) THEN 673 message_string = 'illegal value for data_output: "' // & 674 TRIM( var ) // '" & only 2d-horizontal ' // & 675 'cross sections are allowed for this value' 676 CALL message( 'check_parameters', 'PA0111', 1, 2, 0, 6, 0 ) 677 ENDIF 678 unit = 'IU/s' 679 IF ( .NOT. ALLOCATED( vitd3_exposure ) ) THEN 680 ALLOCATE( vitd3_exposure(nysg:nyng,nxlg:nxrg) ) 681 ENDIF 682 vitd3_exposure = 0.0_wp 683 684 CASE ( 'uvem_vitd3dose*' ) 685 IF ( .NOT. uv_exposure ) THEN 686 message_string = 'output of "' // TRIM( var ) // '" requi' // & 687 'res a namelist &uvexposure_par' 688 CALL message( 'uvem_check_data_output', 'UV0001', 1, 2, 0, 6, 0 ) 689 ENDIF 690 IF ( k == 0 .OR. data_output(i)(ilen-2:ilen) /= '_xy' ) THEN 691 message_string = 'illegal value for data_output: "' // & 692 TRIM( var ) // '" & only 2d-horizontal ' // & 693 'cross sections are allowed for this value' 694 CALL message( 'check_parameters', 'PA0111', 1, 2, 0, 6, 0 ) 695 ENDIF 696 unit = 'IU/av-h' 697 IF ( .NOT. ALLOCATED( vitd3_exposure_av ) ) THEN 698 ALLOCATE( vitd3_exposure_av(nysg:nyng,nxlg:nxrg) ) 699 ENDIF 700 vitd3_exposure_av = 0.0_wp 502 701 503 702 CASE DEFAULT … … 506 705 END SELECT 507 706 508 IF ( unit /= 'illegal' ) THEN707 IF ( thermal_comfort .AND. unit == 'degree_C' ) THEN 509 708 ! 510 709 !-- Futher checks if output belongs to biometeorology 710 ! Break if required modules "radiation" and "humidity" are not running. 511 711 IF ( .NOT. radiation ) THEN 512 712 message_string = 'output of "' // TRIM( var ) // '" require' & … … 541 741 IMPLICIT NONE 542 742 543 743 ! 544 744 !-- Disabled as long as radiation model not available 545 745 546 IF ( .NOT.humidity ) THEN746 IF ( thermal_comfort .AND. .NOT. humidity ) THEN 547 747 message_string = 'The estimation of thermal comfort requires ' // & 548 748 'air humidity information, but humidity module ' // & … … 562 762 !------------------------------------------------------------------------------! 563 763 SUBROUTINE bio_data_output_2d( av, variable, found, grid, local_pf, & 564 two_d, nzb_do, nzt_do, fill_value ) 565 566 USE indices, & 567 ONLY: nxl, nxl, nxr, nxr, nyn, nyn, nys, nys, nzb, nzt 764 two_d, nzb_do, nzt_do ) 765 568 766 569 767 USE kinds … … 571 769 572 770 IMPLICIT NONE 573 771 ! 574 772 !-- Input variables 575 773 CHARACTER (LEN=*), INTENT(IN) :: variable !< Char identifier to select var for output … … 577 775 INTEGER(iwp), INTENT(IN) :: nzb_do !< Unused. 2D. nz bottom to nz top 578 776 INTEGER(iwp), INTENT(IN) :: nzt_do !< Unused. 579 REAL(wp), INTENT(in) :: fill_value !< Fill value for unassigned locations 580 777 ! 581 778 !-- Output variables 582 779 CHARACTER (LEN=*), INTENT(OUT) :: grid !< Grid type (always "zu1" for biom) … … 584 781 LOGICAL, INTENT(OUT) :: two_d !< Flag parameter that indicates 2D variables, horizontal cross sections, must be .TRUE. 585 782 REAL(wp), DIMENSION(nxl:nxr,nys:nyn,nzb_do:nzt_do) :: local_pf !< Temp. result grid to return 586 783 ! 587 784 !-- Internal variables 588 CHARACTER (LEN=:), allocatable :: variable_short !< Trimmed version of char variable589 785 INTEGER(iwp) :: i !< Running index, x-dir 590 786 INTEGER(iwp) :: j !< Running index, y-dir … … 593 789 594 790 595 variable_short = TRIM( variable )596 IF ( variable_short(1:4) /= 'bio_' ) THEN597 found = .FALSE.598 grid = 'none'599 ENDIF600 601 791 found = .TRUE. 602 local_pf = fill_value603 604 SELECT CASE ( variable_short)792 local_pf = bio_fill_value 793 794 SELECT CASE ( TRIM( variable ) ) 605 795 606 796 607 797 CASE ( 'bio_mrt_xy' ) 608 grid = 'zu1' 609 two_d = .FALSE. !< can be calculated for several levels 610 local_pf = REAL( fill_value, KIND = wp ) 798 grid = 'zu1' 799 two_d = .FALSE. !< can be calculated for several levels 800 local_pf = REAL( bio_fill_value, KIND = wp ) 801 DO l = 1, nmrtbl 802 i = mrtbl(ix,l) 803 j = mrtbl(iy,l) 804 k = mrtbl(iz,l) 805 IF ( k < nzb_do .OR. k > nzt_do .OR. j < nys .OR. j > nyn .OR. & 806 i < nxl .OR. i > nxr ) CYCLE 807 IF ( av == 0 ) THEN 808 IF ( mrt_include_sw ) THEN 809 local_pf(i,j,k) = ((human_absorb * mrtinsw(l) + & 810 human_emiss * mrtinlw(l)) / & 811 (human_emiss * sigma_sb)) ** .25_wp - degc_to_k 812 ELSE 813 local_pf(i,j,k) = (human_emiss * mrtinlw(l) / & 814 sigma_sb) ** .25_wp - degc_to_k 815 ENDIF 816 ELSE 817 local_pf(i,j,k) = mrt_av_grid(l) 818 ENDIF 819 ENDDO 820 821 822 CASE ( 'bio_perct*_xy' ) ! 2d-array 823 grid = 'zu1' 824 two_d = .TRUE. 825 IF ( av == 0 ) THEN 826 DO i = nxl, nxr 827 DO j = nys, nyn 828 local_pf(i,j,nzb+1) = perct(j,i) 829 ENDDO 830 ENDDO 831 ELSE 832 DO i = nxl, nxr 833 DO j = nys, nyn 834 local_pf(i,j,nzb+1) = perct_av(j,i) 835 ENDDO 836 ENDDO 837 END IF 838 839 840 CASE ( 'bio_utci*_xy' ) ! 2d-array 841 grid = 'zu1' 842 two_d = .TRUE. 843 IF ( av == 0 ) THEN 844 DO i = nxl, nxr 845 DO j = nys, nyn 846 local_pf(i,j,nzb+1) = utci(j,i) 847 ENDDO 848 ENDDO 849 ELSE 850 DO i = nxl, nxr 851 DO j = nys, nyn 852 local_pf(i,j,nzb+1) = utci_av(j,i) 853 ENDDO 854 ENDDO 855 END IF 856 857 858 CASE ( 'bio_pet*_xy' ) ! 2d-array 859 grid = 'zu1' 860 two_d = .TRUE. 861 IF ( av == 0 ) THEN 862 DO i = nxl, nxr 863 DO j = nys, nyn 864 local_pf(i,j,nzb+1) = pet(j,i) 865 ENDDO 866 ENDDO 867 ELSE 868 DO i = nxl, nxr 869 DO j = nys, nyn 870 local_pf(i,j,nzb+1) = pet_av(j,i) 871 ENDDO 872 ENDDO 873 END IF 874 875 ! 876 !-- Before data is transfered to local_pf, transfer is it 2D dummy variable and exchange ghost points therein. 877 !-- However, at this point this is only required for instantaneous arrays, time-averaged quantities are already exchanged. 878 CASE ( 'uvem_vitd3*_xy' ) ! 2d-array 879 IF ( av == 0 ) THEN 880 DO i = nxl, nxr 881 DO j = nys, nyn 882 local_pf(i,j,nzb+1) = vitd3_exposure(j,i) 883 ENDDO 884 ENDDO 885 ENDIF 886 887 two_d = .TRUE. 888 grid = 'zu1' 889 890 CASE ( 'uvem_vitd3dose*_xy' ) ! 2d-array 891 IF ( av == 1 ) THEN 892 DO i = nxl, nxr 893 DO j = nys, nyn 894 local_pf(i,j,nzb+1) = vitd3_exposure_av(j,i) 895 ENDDO 896 ENDDO 897 ENDIF 898 899 two_d = .TRUE. 900 grid = 'zu1' 901 902 903 CASE DEFAULT 904 found = .FALSE. 905 grid = 'none' 906 907 END SELECT 908 909 910 END SUBROUTINE bio_data_output_2d 911 912 913 !------------------------------------------------------------------------------! 914 ! 915 ! Description: 916 ! ------------ 917 !> Subroutine defining 3D output variables (dummy, always 2d!) 918 !> data_output_3d 709ff 919 !------------------------------------------------------------------------------! 920 SUBROUTINE bio_data_output_3d( av, variable, found, local_pf, nzb_do, nzt_do ) 921 922 USE indices 923 924 USE kinds 925 926 927 IMPLICIT NONE 928 ! 929 !-- Input variables 930 CHARACTER (LEN=*), INTENT(IN) :: variable !< Char identifier to select var for output 931 INTEGER(iwp), INTENT(IN) :: av !< Use averaged data? 0 = no, 1 = yes? 932 INTEGER(iwp), INTENT(IN) :: nzb_do !< Unused. 2D. nz bottom to nz top 933 INTEGER(iwp), INTENT(IN) :: nzt_do !< Unused. 934 ! 935 !-- Output variables 936 LOGICAL, INTENT(OUT) :: found !< Output found? 937 REAL(sp), DIMENSION(nxl:nxr,nys:nyn,nzb_do:nzt_do) :: local_pf !< Temp. result grid to return 938 ! 939 !-- Internal variables 940 INTEGER(iwp) :: l !< Running index, radiation grid 941 INTEGER(iwp) :: i !< Running index, x-dir 942 INTEGER(iwp) :: j !< Running index, y-dir 943 INTEGER(iwp) :: k !< Running index, z-dir 944 945 ! REAL(wp) :: mrt !< Buffer for mean radiant temperature 946 947 found = .TRUE. 948 949 SELECT CASE ( TRIM( variable ) ) 950 951 CASE ( 'bio_mrt' ) 952 local_pf = REAL( bio_fill_value, KIND = sp ) 611 953 DO l = 1, nmrtbl 612 954 i = mrtbl(ix,l) … … 617 959 IF ( av == 0 ) THEN 618 960 IF ( mrt_include_sw ) THEN 619 local_pf(i,j,k) = ((human_absorb * mrtinsw(l) + & 620 human_emiss * mrtinlw(l)) / & 621 (human_emiss * sigma_sb)) ** .25_wp - degc_to_k 961 local_pf(i,j,k) = REAL(((human_absorb * mrtinsw(l) + & 962 human_emiss * mrtinlw(l)) / & 963 (human_emiss * sigma_sb)) ** .25_wp - & 964 degc_to_k, kind=sp ) 622 965 ELSE 623 local_pf(i,j,k) = (human_emiss * mrtinlw(l) /&624 sigma_sb) ** .25_wp - degc_to_k966 local_pf(i,j,k) = REAL((human_emiss * mrtinlw(l) / & 967 sigma_sb) ** .25_wp - degc_to_k, kind=sp ) !< why not (human_emiss * sigma_sb) as above? 625 968 ENDIF 626 969 ELSE 627 local_pf(i,j,k) = mrt_av_grid(l) 628 ENDIF 629 ENDDO 630 631 632 CASE ( 'bio_perct*_xy' ) ! 2d-array 633 grid = 'zu1' 634 two_d = .TRUE. 635 IF ( av == 0 ) THEN 636 DO i = nxl, nxr 637 DO j = nys, nyn 638 local_pf(i,j,nzb+1) = perct(j,i) 639 ENDDO 640 ENDDO 641 ELSE 642 DO i = nxl, nxr 643 DO j = nys, nyn 644 local_pf(i,j,nzb+1) = perct_av(j,i) 645 ENDDO 646 ENDDO 647 END IF 648 649 650 CASE ( 'bio_utci*_xy' ) ! 2d-array 651 grid = 'zu1' 652 two_d = .TRUE. 653 IF ( av == 0 ) THEN 654 DO i = nxl, nxr 655 DO j = nys, nyn 656 local_pf(i,j,nzb+1) = utci_grid(j,i) 657 ENDDO 658 ENDDO 659 ELSE 660 DO i = nxl, nxr 661 DO j = nys, nyn 662 local_pf(i,j,nzb+1) = utci_av_grid(j,i) 663 ENDDO 664 ENDDO 665 END IF 666 667 668 CASE ( 'bio_pet*_xy' ) ! 2d-array 669 grid = 'zu1' 670 two_d = .TRUE. 671 IF ( av == 0 ) THEN 672 DO i = nxl, nxr 673 DO j = nys, nyn 674 local_pf(i,j,nzb+1) = pet_grid(j,i) 675 ENDDO 676 ENDDO 677 ELSE 678 DO i = nxl, nxr 679 DO j = nys, nyn 680 local_pf(i,j,nzb+1) = pet_av_grid(j,i) 681 ENDDO 682 ENDDO 683 END IF 684 685 686 CASE DEFAULT 687 found = .FALSE. 688 grid = 'none' 689 690 END SELECT 691 692 693 END SUBROUTINE bio_data_output_2d 694 695 696 !------------------------------------------------------------------------------! 697 ! 698 ! Description: 699 ! ------------ 700 !> Subroutine defining 3D output variables (dummy, always 2d!) 701 !> data_output_3d 709ff 702 !------------------------------------------------------------------------------! 703 SUBROUTINE bio_data_output_3d( av, variable, found, local_pf, nzb_do, nzt_do ) 704 705 USE indices 706 707 USE kinds 708 709 710 IMPLICIT NONE 711 712 !-- Input variables 713 CHARACTER (LEN=*), INTENT(IN) :: variable !< Char identifier to select var for output 714 INTEGER(iwp), INTENT(IN) :: av !< Use averaged data? 0 = no, 1 = yes? 715 INTEGER(iwp), INTENT(IN) :: nzb_do !< Unused. 2D. nz bottom to nz top 716 INTEGER(iwp), INTENT(IN) :: nzt_do !< Unused. 717 718 !-- Output variables 719 LOGICAL, INTENT(OUT) :: found !< Output found? 720 REAL(sp), DIMENSION(nxl:nxr,nys:nyn,nzb_do:nzt_do) :: local_pf !< Temp. result grid to return 721 722 !-- Internal variables 723 INTEGER(iwp) :: l !< Running index, radiation grid 724 INTEGER(iwp) :: i !< Running index, x-dir 725 INTEGER(iwp) :: j !< Running index, y-dir 726 INTEGER(iwp) :: k !< Running index, z-dir 727 728 CHARACTER (LEN=:), allocatable :: variable_short !< Trimmed version of char variable 729 730 REAL(wp), PARAMETER :: fill_value = -999._wp 731 REAL(wp) :: mrt !< Buffer for mean radiant temperature 732 733 found = .TRUE. 734 variable_short = TRIM( variable ) 735 736 IF ( variable_short(1:4) /= 'bio_' ) THEN 737 !-- Nothing to do, set found to FALSE and return immediatelly 738 found = .FALSE. 739 RETURN 740 ENDIF 741 742 SELECT CASE ( variable_short ) 743 744 CASE ( 'bio_mrt' ) 745 local_pf = REAL( fill_value, KIND = wp ) 746 DO l = 1, nmrtbl 747 i = mrtbl(ix,l) 748 j = mrtbl(iy,l) 749 k = mrtbl(iz,l) 750 IF ( k < nzb_do .OR. k > nzt_do .OR. j < nys .OR. j > nyn .OR. & 751 i < nxl .OR. i > nxr ) CYCLE 752 IF ( av == 0 ) THEN 753 IF ( mrt_include_sw ) THEN 754 local_pf(i,j,k) = ((human_absorb * mrtinsw(l) + human_emiss * mrtinlw(l)) / & 755 (human_emiss * sigma_sb)) ** .25_wp - degc_to_k 756 ELSE 757 local_pf(i,j,k) = (human_emiss * mrtinlw(l) / & 758 sigma_sb) ** .25_wp - degc_to_k 759 ENDIF 760 ELSE 761 local_pf(i,j,k) = mrt_av_grid(l) 970 local_pf(i,j,k) = REAL(mrt_av_grid(l), kind=sp) 762 971 ENDIF 763 972 ENDDO … … 780 989 781 990 IMPLICIT NONE 782 991 ! 783 992 !-- Input variables 784 993 CHARACTER (LEN=*), INTENT(IN) :: var !< Name of output variable 785 994 ! 786 995 !-- Output variables 787 996 CHARACTER (LEN=*), INTENT(OUT) :: grid_x !< x grid of output variable … … 790 999 791 1000 LOGICAL, INTENT(OUT) :: found !< Flag if output var is found 792 1001 ! 793 1002 !-- Local variables 794 1003 LOGICAL :: is2d !< Var is 2d? … … 810 1019 grid_y = 'y' 811 1020 grid_z = 'zu' 812 IF ( is2d ) grid_z = 'zu1' 1021 IF ( is2d .AND. var(1:7) /= 'bio_mrt' ) grid_z = 'zu1' 1022 ENDIF 1023 1024 IF ( is2d .AND. var(1:4) == 'uvem' ) THEN 1025 grid_x = 'x' 1026 grid_y = 'y' 1027 grid_z = 'zu1' 813 1028 ENDIF 814 1029 … … 824 1039 825 1040 IMPLICIT NONE 826 1041 ! 827 1042 !-- Input variables 828 1043 INTEGER(iwp), INTENT(IN) :: io !< Unit of the output file 829 1044 ! 830 1045 !-- Internal variables 831 1046 CHARACTER (LEN=86) :: output_height_chr !< String for output height … … 857 1072 ONLY: message_string 858 1073 1074 USE netcdf_data_input_mod, & 1075 ONLY: netcdf_data_input_uvem, uvem_projarea_f, uvem_radiance_f, & 1076 uvem_irradiance_f, uvem_integration_f, building_obstruction_f 1077 859 1078 IMPLICIT NONE 860 1079 ! 861 1080 !-- Internal vriables 862 1081 REAL ( wp ) :: height !< current height in meters 863 864 INTEGER ( iwp ) :: i !< iteration index 865 1082 ! 866 1083 !-- Determine cell level corresponding to 1.1 m above ground level 867 1084 ! (gravimetric center of sample human) … … 876 1093 877 1094 IF ( .NOT. radiation_interactions ) THEN 878 message_string = 'The mrt calculation requires ' // &879 'enabled radiation_interactions but it ' // &1095 message_string = 'The mrt calculation requires ' // & 1096 'enabled radiation_interactions but it ' // & 880 1097 'is disabled!' 881 1098 CALL message( 'check_parameters', 'PAHU03', 0, 0, -1, 6, 0 ) 882 1099 ENDIF 883 1100 1101 ! 1102 !-- Init UVEM and load lookup tables 1103 CALL netcdf_data_input_uvem 1104 884 1105 END SUBROUTINE bio_init 885 886 887 !------------------------------------------------------------------------------!888 ! Description:889 ! ------------890 !> Allocate biom arrays and define pointers if required891 !> init_3d_model 1047ff892 !------------------------------------------------------------------------------!893 SUBROUTINE bio_init_arrays894 895 IMPLICIT NONE896 897 !-- Allocate a temporary array with the desired output dimensions.898 ! Initialization omitted for performance, will be overwritten anyway899 IF ( .NOT. ALLOCATED( tmrt_grid ) ) THEN900 ALLOCATE( tmrt_grid (nys:nyn,nxl:nxr) )901 ENDIF902 903 IF ( bio_perct ) THEN904 IF ( .NOT. ALLOCATED( perct ) ) THEN905 ALLOCATE( perct (nys:nyn,nxl:nxr) )906 ENDIF907 ENDIF908 909 IF ( bio_utci ) THEN910 IF ( .NOT. ALLOCATED( utci_grid ) ) THEN911 ALLOCATE( utci_grid (nys:nyn,nxl:nxr) )912 ENDIF913 ENDIF914 915 IF ( bio_pet ) THEN916 IF ( .NOT. ALLOCATED( pet_grid ) ) THEN917 ALLOCATE( pet_grid (nys:nyn,nxl:nxr) )918 ENDIF919 END IF920 921 IF ( bio_perct_av ) THEN922 IF ( .NOT. ALLOCATED( perct_av ) ) THEN923 ALLOCATE( perct_av (nys:nyn,nxl:nxr) )924 ENDIF925 ENDIF926 927 IF ( bio_utci_av ) THEN928 IF ( .NOT. ALLOCATED( utci_av_grid ) ) THEN929 ALLOCATE( utci_av_grid (nys:nyn,nxl:nxr) )930 ENDIF931 ENDIF932 933 IF ( bio_pet_av ) THEN934 IF ( .NOT. ALLOCATED( pet_av_grid ) ) THEN935 ALLOCATE( pet_av_grid (nys:nyn,nxl:nxr) )936 ENDIF937 938 END IF939 940 END SUBROUTINE bio_init_arrays941 1106 942 1107 … … 959 1124 bio_perct_av, & 960 1125 bio_utci, & 961 bio_utci_av 1126 bio_utci_av, & 1127 thermal_comfort, & 1128 ! 1129 !-- UVEM namelist parameters 1130 clothing, & 1131 consider_obstructions, & 1132 orientation_angle, & 1133 sun_in_south, & 1134 turn_to_sun, & 1135 uv_exposure 962 1136 963 1137 … … 1003 1177 1004 1178 LOGICAL, INTENT(IN) :: av !< use averaged input? 1179 ! 1005 1180 !-- Internal variables 1006 1181 INTEGER(iwp) :: i !< Running index, x-dir, radiation coordinates … … 1009 1184 INTEGER(iwp) :: l !< Running index, radiation coordinates 1010 1185 1011 1186 ! 1012 1187 !-- Calculate biometeorology MRT from local radiation 1013 ! --fluxes calculated by RTM and assign into 2D grid1014 tmrt_grid = -999.1188 ! fluxes calculated by RTM and assign into 2D grid 1189 tmrt_grid = bio_fill_value 1015 1190 DO l = 1, nmrtbl 1016 1191 i = mrtbl(ix,l) … … 1042 1217 1043 1218 IMPLICIT NONE 1044 1219 ! 1045 1220 !-- Input variables 1046 1221 LOGICAL, INTENT ( IN ) :: average_input !< Determine averaged input conditions? 1047 1222 INTEGER(iwp), INTENT ( IN ) :: i !< Running index, x-dir 1048 1223 INTEGER(iwp), INTENT ( IN ) :: j !< Running index, y-dir 1049 1224 ! 1050 1225 !-- Output parameters 1051 1226 REAL(wp), INTENT ( OUT ) :: tmrt !< Mean radiant temperature (°C) … … 1054 1229 REAL(wp), INTENT ( OUT ) :: ws !< Wind speed (local level) (m/s) 1055 1230 REAL(wp), INTENT ( OUT ) :: pair !< Air pressure (hPa) 1056 1231 ! 1057 1232 !-- Internal variables 1058 1233 INTEGER(iwp) :: k !< Running index, z-dir 1059 INTEGER(iwp) :: ir !< Running index, x-dir, radiation coordinates1060 INTEGER(iwp) :: jr !< Running index, y-dir, radiation coordinates1061 INTEGER(iwp) :: kr !< Running index, y-dir, radiation coordinates1234 ! INTEGER(iwp) :: ir !< Running index, x-dir, radiation coordinates 1235 ! INTEGER(iwp) :: jr !< Running index, y-dir, radiation coordinates 1236 ! INTEGER(iwp) :: kr !< Running index, y-dir, radiation coordinates 1062 1237 INTEGER(iwp) :: k_wind !< Running index, z-dir, wind speed only 1063 INTEGER(iwp) :: l !< Running index, radiation coordinates1238 ! INTEGER(iwp) :: l !< Running index, radiation coordinates 1064 1239 1065 1240 REAL(wp) :: vp_sat !< Saturation vapor pressure (hPa) 1066 1241 1067 1242 ! 1068 1243 !-- Determine cell level closest to 1.1m above ground 1069 1244 ! by making use of truncation due to int cast … … 1076 1251 k_wind = k + 1_iwp 1077 1252 ENDIF 1078 1253 ! 1079 1254 !-- Determine local values: 1080 1255 IF ( average_input ) THEN 1256 ! 1081 1257 !-- Calculate ta from Tp assuming dry adiabatic laps rate 1082 1258 ta = pt_av(k, j, i) - ( 0.0098_wp * dz(1) * ( k + .5_wp ) ) - degc_to_k 1083 1259 1084 vp = -999._wp1260 vp = bio_fill_value 1085 1261 IF ( humidity .AND. ALLOCATED( q_av ) ) THEN 1086 1262 vp = q_av(k, j, i) … … 1091 1267 0.5_wp * ABS( w_av(k_wind, j, i) + w_av(k_wind+1, j, i) ) ) 1092 1268 ELSE 1269 ! 1093 1270 !-- Calculate ta from Tp assuming dry adiabatic laps rate 1094 1271 ta = pt(k, j, i) - ( 0.0098_wp * dz(1) * ( k + .5_wp ) ) - degc_to_k 1095 1272 1096 vp = -999._wp1273 vp = bio_fill_value 1097 1274 IF ( humidity ) THEN 1098 1275 vp = q(k, j, i) … … 1104 1281 1105 1282 ENDIF 1106 1283 ! 1107 1284 !-- Local air pressure 1108 1285 pair = surface_pressure … … 1116 1293 IF ( vp > vp_sat ) vp = vp_sat 1117 1294 ENDIF 1118 1295 ! 1119 1296 !-- local mtr value at [i,j] 1120 tmrt = -999.!< this can be a valid result (e.g. for inside some ostacle)1297 tmrt = bio_fill_value !< this can be a valid result (e.g. for inside some ostacle) 1121 1298 IF ( radiation ) THEN 1299 ! 1122 1300 !-- Use MRT from RTM precalculated in tmrt_grid 1123 1301 tmrt = tmrt_grid(j,i) … … 1136 1314 1137 1315 IMPLICIT NONE 1138 1316 ! 1139 1317 !-- Input attributes 1140 1318 LOGICAL, INTENT ( IN ) :: av !< Calculate based on averaged input conditions? 1141 1319 ! 1142 1320 !-- Internal variables 1143 INTEGER(iwp) :: i, j !< Running index 1144 1145 REAL(wp) :: clo !< Clothing index (no dimension) 1146 REAL(wp) :: ta !< Air temperature (°C) 1147 REAL(wp) :: vp !< Vapour pressure (hPa) 1148 REAL(wp) :: ws !< Wind speed (local level) (m/s) 1149 REAL(wp) :: pair !< Air pressure (hPa) 1150 REAL(wp) :: perct_tmp !< Perceived temperature (°C) 1151 REAL(wp) :: utci_tmp !< Universal thermal climate index (°C) 1152 REAL(wp) :: pet_tmp !< Physiologically equivalent temperature (°C) 1153 REAL(wp) :: tmrt_tmp !< Mean radiant temperature (°C) 1154 1155 CALL bio_init_arrays 1156 1321 INTEGER(iwp) :: i, j !< Running index 1322 1323 REAL(wp) :: clo !< Clothing index (no dimension) 1324 REAL(wp) :: ta !< Air temperature (°C) 1325 REAL(wp) :: vp !< Vapour pressure (hPa) 1326 REAL(wp) :: ws !< Wind speed (local level) (m/s) 1327 REAL(wp) :: pair !< Air pressure (hPa) 1328 REAL(wp) :: perct_ij !< Perceived temperature (°C) 1329 REAL(wp) :: utci_ij !< Universal thermal climate index (°C) 1330 REAL(wp) :: pet_ij !< Physiologically equivalent temperature (°C) 1331 REAL(wp) :: tmrt_ij !< Mean radiant temperature (°C) 1332 1333 ! 1157 1334 !-- fill out the MRT 2D grid from appropriate source (RTM, RRTMG,...) 1158 CALL bio_calculate_mrt_grid ( av ) 1159 1335 IF ( simulated_time > 0.0_wp ) THEN 1336 CALL bio_calculate_mrt_grid ( av ) 1337 ENDIF 1160 1338 1161 1339 DO i = nxl, nxr 1162 1340 DO j = nys, nyn 1341 ! 1163 1342 !-- Determine local input conditions 1164 CALL bio_get_thermal_index_input_ij ( av, i, j, ta, vp, & 1165 ws, pair, tmrt_tmp ) 1166 ! tmrt_grid(j, i) = tmrt_tmp !< already set by bio_calculate_mrt_grid! 1167 1343 tmrt_ij = bio_fill_value 1344 vp = bio_fill_value 1345 ! 1346 !-- Determine input only if 1347 IF ( simulated_time > 0.0_wp ) THEN 1348 CALL bio_get_thermal_index_input_ij ( av, i, j, ta, vp, & 1349 ws, pair, tmrt_ij ) 1350 END IF 1351 ! 1168 1352 !-- Only proceed if input is available 1169 IF ( tmrt_tmp <= -998._wp .OR. vp <= -998._wp ) THEN1170 pet_tmp = -999._wp !< set fail value, e.g. valid for within1171 perct_tmp = -999._wp !< some obstacle1172 utci_tmp = -999._wp1173 ELSE 1353 pet_ij = bio_fill_value !< set fail value, e.g. valid for within 1354 perct_ij = bio_fill_value !< some obstacle 1355 utci_ij = bio_fill_value 1356 IF ( .NOT. ( tmrt_ij <= -998._wp .OR. vp <= -998._wp ) ) THEN 1357 ! 1174 1358 !-- Calculate static thermal indices based on local tmrt 1175 clo = -999._wp1176 1359 clo = bio_fill_value 1360 1177 1361 IF ( bio_perct ) THEN 1362 ! 1178 1363 !-- Estimate local perceived temperature 1179 CALL calculate_perct_static( ta, vp, ws, tmrt_ tmp, pair, clo,&1180 perct_ tmp)1364 CALL calculate_perct_static( ta, vp, ws, tmrt_ij, pair, clo, & 1365 perct_ij ) 1181 1366 ENDIF 1182 1367 1183 1368 IF ( bio_utci ) THEN 1369 ! 1184 1370 !-- Estimate local universal thermal climate index 1185 CALL calculate_utci_static( ta, vp, ws, tmrt_ tmp,&1186 bio_output_height, utci_ tmp)1371 CALL calculate_utci_static( ta, vp, ws, tmrt_ij, & 1372 bio_output_height, utci_ij ) 1187 1373 ENDIF 1188 1374 1189 1375 IF ( bio_pet ) THEN 1376 ! 1190 1377 !-- Estimate local physiologically equivalent temperature 1191 CALL calculate_pet_static( ta, vp, ws, tmrt_ tmp, pair, pet_tmp)1378 CALL calculate_pet_static( ta, vp, ws, tmrt_ij, pair, pet_ij ) 1192 1379 ENDIF 1193 1380 END IF … … 1195 1382 1196 1383 IF ( av ) THEN 1384 ! 1197 1385 !-- Write results for selected averaged indices 1198 1386 IF ( bio_perct_av ) THEN 1199 perct_av(j, i) = perct_ tmp1387 perct_av(j, i) = perct_ij 1200 1388 END IF 1201 1389 IF ( bio_utci_av ) THEN 1202 utci_av _grid(j, i) = utci_tmp1390 utci_av(j, i) = utci_ij 1203 1391 END IF 1204 1392 IF ( bio_pet_av ) THEN 1205 pet_av _grid(j, i) = pet_tmp1393 pet_av(j, i) = pet_ij 1206 1394 END IF 1207 1395 ELSE 1396 ! 1208 1397 !-- Write result for selected indices 1209 1398 IF ( bio_perct ) THEN 1210 perct(j, i) = perct_ tmp1399 perct(j, i) = perct_ij 1211 1400 END IF 1212 1401 IF ( bio_utci ) THEN 1213 utci _grid(j, i) = utci_tmp1402 utci(j, i) = utci_ij 1214 1403 END IF 1215 1404 IF ( bio_pet ) THEN 1216 pet _grid(j, i) = pet_tmp1405 pet(j, i) = pet_ij 1217 1406 END IF 1218 1407 END IF … … 1232 1421 1233 1422 IMPLICIT NONE 1234 1423 ! 1235 1424 !-- Input parameters 1236 1425 REAL(wp), INTENT ( IN ) :: ta !< Air temperature (°C) … … 1246 1435 ! (without metabolism!) (W) 1247 1436 INTEGER(iwp), INTENT ( IN ) :: sex !< Sex of agent (1 = male, 2 = female) 1248 1437 ! 1249 1438 !-- Both, input and output parameters 1250 1439 Real(wp), INTENT ( INOUT ) :: energy_storage !< Energy storage (W/m²) … … 1253 1442 Real(wp), INTENT ( INOUT ) :: actlev !< Individuals activity level 1254 1443 ! per unit surface area (W/m²) 1444 ! 1255 1445 !-- Output parameters 1256 1446 REAL(wp), INTENT ( OUT ) :: ipt !< Instationary perceived temp. (°C) 1257 1447 ! 1258 1448 !-- If clo equals the initial value, this is the initial call 1259 1449 IF ( clo <= -998._wp ) THEN 1450 ! 1260 1451 !-- Initialize instationary perceived temperature with personalized 1261 1452 ! PT as an initial guess, set actlev and clo … … 1264 1455 ipt ) 1265 1456 ELSE 1457 ! 1266 1458 !-- Estimate local instatinoary perceived temperature 1267 1459 CALL ipt_cycle ( ta, vp, ws, tmrt, pair, dt, energy_storage, t_clo, & … … 1287 1479 !> www.utci.org 1288 1480 !------------------------------------------------------------------------------! 1289 SUBROUTINE calculate_utci_static( ta_in, vp, ws_hag, tmrt, hag, utci )1481 SUBROUTINE calculate_utci_static( ta_in, vp, ws_hag, tmrt, hag, utci_ij ) 1290 1482 1291 1483 IMPLICIT NONE 1292 1484 ! 1293 1485 !-- Type of input of the argument list 1294 REAL(WP), INTENT ( IN ) :: ta_in !< Local air temperature (°C) 1295 REAL(WP), INTENT ( IN ) :: vp !< Loacl vapour pressure (hPa) 1296 REAL(WP), INTENT ( IN ) :: ws_hag !< Incident wind speed (m/s) 1297 REAL(WP), INTENT ( IN ) :: tmrt !< Local mean radiant temperature (°C) 1298 REAL(WP), INTENT ( IN ) :: hag !< Height of wind speed input (m) 1486 REAL(WP), INTENT ( IN ) :: ta_in !< Local air temperature (°C) 1487 REAL(WP), INTENT ( IN ) :: vp !< Loacl vapour pressure (hPa) 1488 REAL(WP), INTENT ( IN ) :: ws_hag !< Incident wind speed (m/s) 1489 REAL(WP), INTENT ( IN ) :: tmrt !< Local mean radiant temperature (°C) 1490 REAL(WP), INTENT ( IN ) :: hag !< Height of wind speed input (m) 1491 ! 1299 1492 !-- Type of output of the argument list 1300 REAL(wp), INTENT ( OUT ) :: utci !< Universal Thermal Climate Index (°C) 1301 1302 !-- Make sure precission is sufficient for regression equation 1303 REAL(WP) :: ta !< internal air temperature (°C) 1304 REAL(WP) :: pa !< air pressure in kPa (kPa) 1305 REAL(WP) :: d_tmrt !< delta-tmrt (°C) 1306 REAL(WP) :: va !< wind speed at 10 m above ground level (m/s) 1307 REAL(WP) :: offset !< utci deviation by ta cond. exceeded (°C) 1308 1493 REAL(wp), INTENT ( OUT ) :: utci_ij !< Universal Thermal Climate Index (°C) 1494 1495 REAL(WP) :: ta !< air temperature modified by offset (°C) 1496 REAL(WP) :: pa !< air pressure in kPa (kPa) 1497 REAL(WP) :: d_tmrt !< delta-tmrt (°C) 1498 REAL(WP) :: va !< wind speed at 10 m above ground level (m/s) 1499 REAL(WP) :: offset !< utci deviation by ta cond. exceeded (°C) 1500 REAL(WP) :: part_ta !< Air temperature related part of the regression 1501 REAL(WP) :: ta2 !< 2 times ta 1502 REAL(WP) :: ta3 !< 3 times ta 1503 REAL(WP) :: ta4 !< 4 times ta 1504 REAL(WP) :: ta5 !< 5 times ta 1505 REAL(WP) :: ta6 !< 6 times ta 1506 REAL(WP) :: part_va !< Vapour pressure related part of the regression 1507 REAL(WP) :: va2 !< 2 times va 1508 REAL(WP) :: va3 !< 3 times va 1509 REAL(WP) :: va4 !< 4 times va 1510 REAL(WP) :: va5 !< 5 times va 1511 REAL(WP) :: va6 !< 6 times va 1512 REAL(WP) :: part_d_tmrt !< Mean radiant temp. related part of the reg. 1513 REAL(WP) :: d_tmrt2 !< 2 times d_tmrt 1514 REAL(WP) :: d_tmrt3 !< 3 times d_tmrt 1515 REAL(WP) :: d_tmrt4 !< 4 times d_tmrt 1516 REAL(WP) :: d_tmrt5 !< 5 times d_tmrt 1517 REAL(WP) :: d_tmrt6 !< 6 times d_tmrt 1518 REAL(WP) :: part_pa !< Air pressure related part of the regression 1519 REAL(WP) :: pa2 !< 2 times pa 1520 REAL(WP) :: pa3 !< 3 times pa 1521 REAL(WP) :: pa4 !< 4 times pa 1522 REAL(WP) :: pa5 !< 5 times pa 1523 REAL(WP) :: pa6 !< 6 times pa 1524 REAL(WP) :: part_pa2 !< Air pressure^2 related part of the regression 1525 REAL(WP) :: part_pa3 !< Air pressure^3 related part of the regression 1526 REAL(WP) :: part_pa46 !< Air pressure^4-6 related part of the regression 1527 1528 ! 1309 1529 !-- Initialize 1310 1530 offset = 0._wp 1311 1531 ta = ta_in 1312 1532 d_tmrt = tmrt - ta_in 1313 1533 ! 1314 1534 !-- Use vapour pressure in kpa 1315 1535 pa = vp / 10.0_wp 1316 1536 ! 1317 1537 !-- Wind altitude correction from hag to 10m after Broede et al. (2012), eq.3 1318 1538 ! z(0) is set to 0.01 according to UTCI profile definition 1319 1539 va = ws_hag * log ( 10.0_wp / 0.01_wp ) / log ( hag / 0.01_wp ) 1320 1540 ! 1321 1541 !-- Check if input values in range after Broede et al. (2012) 1322 1542 IF ( ( d_tmrt > 70._wp ) .OR. ( d_tmrt < -30._wp ) .OR. & 1323 1543 ( vp >= 50._wp ) ) THEN 1324 utci = -999._wp1544 utci_ij = bio_fill_value 1325 1545 RETURN 1326 1546 ENDIF 1327 1547 ! 1328 1548 !-- Apply eq. 2 in Broede et al. (2012) for ta out of bounds 1329 1549 IF ( ta > 50._wp ) THEN … … 1335 1555 ta = -50._wp 1336 1556 ENDIF 1337 1557 ! 1338 1558 !-- For routine application. For wind speeds and relative 1339 1559 ! humidity values below 0.5 m/s or 5%, respectively, the … … 1342 1562 IF ( va > 17._wp ) va = 17._wp 1343 1563 1564 ! 1565 !-- Pre-calculate multiples of input parameters to save time later 1566 1567 ta2 = ta * ta 1568 ta3 = ta2 * ta 1569 ta4 = ta3 * ta 1570 ta5 = ta4 * ta 1571 ta6 = ta5 * ta 1572 1573 va2 = va * va 1574 va3 = va2 * va 1575 va4 = va3 * va 1576 va5 = va4 * va 1577 va6 = va5 * va 1578 1579 d_tmrt2 = d_tmrt * d_tmrt 1580 d_tmrt3 = d_tmrt2 * d_tmrt 1581 d_tmrt4 = d_tmrt3 * d_tmrt 1582 d_tmrt5 = d_tmrt4 * d_tmrt 1583 d_tmrt6 = d_tmrt5 * d_tmrt 1584 1585 pa2 = pa * pa 1586 pa3 = pa2 * pa 1587 pa4 = pa3 * pa 1588 pa5 = pa4 * pa 1589 pa6 = pa5 * pa 1590 1591 ! 1592 !-- Pre-calculate parts of the regression equation 1593 part_ta = ( 6.07562052e-01_wp ) + & 1594 ( -2.27712343e-02_wp ) * ta + & 1595 ( 8.06470249e-04_wp ) * ta2 + & 1596 ( -1.54271372e-04_wp ) * ta3 + & 1597 ( -3.24651735e-06_wp ) * ta4 + & 1598 ( 7.32602852e-08_wp ) * ta5 + & 1599 ( 1.35959073e-09_wp ) * ta6 1600 1601 part_va = ( -2.25836520e+00_wp ) * va + & 1602 ( 8.80326035e-02_wp ) * ta * va + & 1603 ( 2.16844454e-03_wp ) * ta2 * va + & 1604 ( -1.53347087e-05_wp ) * ta3 * va + & 1605 ( -5.72983704e-07_wp ) * ta4 * va + & 1606 ( -2.55090145e-09_wp ) * ta5 * va + & 1607 ( -7.51269505e-01_wp ) * va2 + & 1608 ( -4.08350271e-03_wp ) * ta * va2 + & 1609 ( -5.21670675e-05_wp ) * ta2 * va2 + & 1610 ( 1.94544667e-06_wp ) * ta3 * va2 + & 1611 ( 1.14099531e-08_wp ) * ta4 * va2 + & 1612 ( 1.58137256e-01_wp ) * va3 + & 1613 ( -6.57263143e-05_wp ) * ta * va3 + & 1614 ( 2.22697524e-07_wp ) * ta2 * va3 + & 1615 ( -4.16117031e-08_wp ) * ta3 * va3 + & 1616 ( -1.27762753e-02_wp ) * va4 + & 1617 ( 9.66891875e-06_wp ) * ta * va4 + & 1618 ( 2.52785852e-09_wp ) * ta2 * va4 + & 1619 ( 4.56306672e-04_wp ) * va5 + & 1620 ( -1.74202546e-07_wp ) * ta * va5 + & 1621 ( -5.91491269e-06_wp ) * va6 1622 1623 part_d_tmrt = ( 3.98374029e-01_wp ) * d_tmrt + & 1624 ( 1.83945314e-04_wp ) * ta * d_tmrt + & 1625 ( -1.73754510e-04_wp ) * ta2 * d_tmrt + & 1626 ( -7.60781159e-07_wp ) * ta3 * d_tmrt + & 1627 ( 3.77830287e-08_wp ) * ta4 * d_tmrt + & 1628 ( 5.43079673e-10_wp ) * ta5 * d_tmrt + & 1629 ( -2.00518269e-02_wp ) * va * d_tmrt + & 1630 ( 8.92859837e-04_wp ) * ta * va * d_tmrt + & 1631 ( 3.45433048e-06_wp ) * ta2 * va * d_tmrt + & 1632 ( -3.77925774e-07_wp ) * ta3 * va * d_tmrt + & 1633 ( -1.69699377e-09_wp ) * ta4 * va * d_tmrt + & 1634 ( 1.69992415e-04_wp ) * va2 * d_tmrt + & 1635 ( -4.99204314e-05_wp ) * ta * va2 * d_tmrt + & 1636 ( 2.47417178e-07_wp ) * ta2 * va2 * d_tmrt + & 1637 ( 1.07596466e-08_wp ) * ta3 * va2 * d_tmrt + & 1638 ( 8.49242932e-05_wp ) * va3 * d_tmrt + & 1639 ( 1.35191328e-06_wp ) * ta * va3 * d_tmrt + & 1640 ( -6.21531254e-09_wp ) * ta2 * va3 * d_tmrt + & 1641 ( -4.99410301e-06_wp ) * va4 * d_tmrt + & 1642 ( -1.89489258e-08_wp ) * ta * va4 * d_tmrt + & 1643 ( 8.15300114e-08_wp ) * va5 * d_tmrt + & 1644 ( 7.55043090e-04_wp ) * d_tmrt2 + & 1645 ( -5.65095215e-05_wp ) * ta * d_tmrt2 + & 1646 ( -4.52166564e-07_wp ) * ta2 * d_tmrt2 + & 1647 ( 2.46688878e-08_wp ) * ta3 * d_tmrt2 + & 1648 ( 2.42674348e-10_wp ) * ta4 * d_tmrt2 + & 1649 ( 1.54547250e-04_wp ) * va * d_tmrt2 + & 1650 ( 5.24110970e-06_wp ) * ta * va * d_tmrt2 + & 1651 ( -8.75874982e-08_wp ) * ta2 * va * d_tmrt2 + & 1652 ( -1.50743064e-09_wp ) * ta3 * va * d_tmrt2 + & 1653 ( -1.56236307e-05_wp ) * va2 * d_tmrt2 + & 1654 ( -1.33895614e-07_wp ) * ta * va2 * d_tmrt2 + & 1655 ( 2.49709824e-09_wp ) * ta2 * va2 * d_tmrt2 + & 1656 ( 6.51711721e-07_wp ) * va3 * d_tmrt2 + & 1657 ( 1.94960053e-09_wp ) * ta * va3 * d_tmrt2 + & 1658 ( -1.00361113e-08_wp ) * va4 * d_tmrt2 + & 1659 ( -1.21206673e-05_wp ) * d_tmrt3 + & 1660 ( -2.18203660e-07_wp ) * ta * d_tmrt3 + & 1661 ( 7.51269482e-09_wp ) * ta2 * d_tmrt3 + & 1662 ( 9.79063848e-11_wp ) * ta3 * d_tmrt3 + & 1663 ( 1.25006734e-06_wp ) * va * d_tmrt3 + & 1664 ( -1.81584736e-09_wp ) * ta * va * d_tmrt3 + & 1665 ( -3.52197671e-10_wp ) * ta2 * va * d_tmrt3 + & 1666 ( -3.36514630e-08_wp ) * va2 * d_tmrt3 + & 1667 ( 1.35908359e-10_wp ) * ta * va2 * d_tmrt3 + & 1668 ( 4.17032620e-10_wp ) * va3 * d_tmrt3 + & 1669 ( -1.30369025e-09_wp ) * d_tmrt4 + & 1670 ( 4.13908461e-10_wp ) * ta * d_tmrt4 + & 1671 ( 9.22652254e-12_wp ) * ta2 * d_tmrt4 + & 1672 ( -5.08220384e-09_wp ) * va * d_tmrt4 + & 1673 ( -2.24730961e-11_wp ) * ta * va * d_tmrt4 + & 1674 ( 1.17139133e-10_wp ) * va2 * d_tmrt4 + & 1675 ( 6.62154879e-10_wp ) * d_tmrt5 + & 1676 ( 4.03863260e-13_wp ) * ta * d_tmrt5 + & 1677 ( 1.95087203e-12_wp ) * va * d_tmrt5 + & 1678 ( -4.73602469e-12_wp ) * d_tmrt6 1679 1680 part_pa = ( 5.12733497e+00_wp ) * pa + & 1681 ( -3.12788561e-01_wp ) * ta * pa + & 1682 ( -1.96701861e-02_wp ) * ta2 * pa + & 1683 ( 9.99690870e-04_wp ) * ta3 * pa + & 1684 ( 9.51738512e-06_wp ) * ta4 * pa + & 1685 ( -4.66426341e-07_wp ) * ta5 * pa + & 1686 ( 5.48050612e-01_wp ) * va * pa + & 1687 ( -3.30552823e-03_wp ) * ta * va * pa + & 1688 ( -1.64119440e-03_wp ) * ta2 * va * pa + & 1689 ( -5.16670694e-06_wp ) * ta3 * va * pa + & 1690 ( 9.52692432e-07_wp ) * ta4 * va * pa + & 1691 ( -4.29223622e-02_wp ) * va2 * pa + & 1692 ( 5.00845667e-03_wp ) * ta * va2 * pa + & 1693 ( 1.00601257e-06_wp ) * ta2 * va2 * pa + & 1694 ( -1.81748644e-06_wp ) * ta3 * va2 * pa + & 1695 ( -1.25813502e-03_wp ) * va3 * pa + & 1696 ( -1.79330391e-04_wp ) * ta * va3 * pa + & 1697 ( 2.34994441e-06_wp ) * ta2 * va3 * pa + & 1698 ( 1.29735808e-04_wp ) * va4 * pa + & 1699 ( 1.29064870e-06_wp ) * ta * va4 * pa + & 1700 ( -2.28558686e-06_wp ) * va5 * pa + & 1701 ( -3.69476348e-02_wp ) * d_tmrt * pa + & 1702 ( 1.62325322e-03_wp ) * ta * d_tmrt * pa + & 1703 ( -3.14279680e-05_wp ) * ta2 * d_tmrt * pa + & 1704 ( 2.59835559e-06_wp ) * ta3 * d_tmrt * pa + & 1705 ( -4.77136523e-08_wp ) * ta4 * d_tmrt * pa + & 1706 ( 8.64203390e-03_wp ) * va * d_tmrt * pa + & 1707 ( -6.87405181e-04_wp ) * ta * va * d_tmrt * pa + & 1708 ( -9.13863872e-06_wp ) * ta2 * va * d_tmrt * pa + & 1709 ( 5.15916806e-07_wp ) * ta3 * va * d_tmrt * pa + & 1710 ( -3.59217476e-05_wp ) * va2 * d_tmrt * pa + & 1711 ( 3.28696511e-05_wp ) * ta * va2 * d_tmrt * pa + & 1712 ( -7.10542454e-07_wp ) * ta2 * va2 * d_tmrt * pa + & 1713 ( -1.24382300e-05_wp ) * va3 * d_tmrt * pa + & 1714 ( -7.38584400e-09_wp ) * ta * va3 * d_tmrt * pa + & 1715 ( 2.20609296e-07_wp ) * va4 * d_tmrt * pa + & 1716 ( -7.32469180e-04_wp ) * d_tmrt2 * pa + & 1717 ( -1.87381964e-05_wp ) * ta * d_tmrt2 * pa + & 1718 ( 4.80925239e-06_wp ) * ta2 * d_tmrt2 * pa + & 1719 ( -8.75492040e-08_wp ) * ta3 * d_tmrt2 * pa + & 1720 ( 2.77862930e-05_wp ) * va * d_tmrt2 * pa + & 1721 ( -5.06004592e-06_wp ) * ta * va * d_tmrt2 * pa + & 1722 ( 1.14325367e-07_wp ) * ta2 * va * d_tmrt2 * pa + & 1723 ( 2.53016723e-06_wp ) * va2 * d_tmrt2 * pa + & 1724 ( -1.72857035e-08_wp ) * ta * va2 * d_tmrt2 * pa + & 1725 ( -3.95079398e-08_wp ) * va3 * d_tmrt2 * pa + & 1726 ( -3.59413173e-07_wp ) * d_tmrt3 * pa + & 1727 ( 7.04388046e-07_wp ) * ta * d_tmrt3 * pa + & 1728 ( -1.89309167e-08_wp ) * ta2 * d_tmrt3 * pa + & 1729 ( -4.79768731e-07_wp ) * va * d_tmrt3 * pa + & 1730 ( 7.96079978e-09_wp ) * ta * va * d_tmrt3 * pa + & 1731 ( 1.62897058e-09_wp ) * va2 * d_tmrt3 * pa + & 1732 ( 3.94367674e-08_wp ) * d_tmrt4 * pa + & 1733 ( -1.18566247e-09_wp ) * ta * d_tmrt4 * pa + & 1734 ( 3.34678041e-10_wp ) * va * d_tmrt4 * pa + & 1735 ( -1.15606447e-10_wp ) * d_tmrt5 * pa 1736 1737 part_pa2 = ( -2.80626406e+00_wp ) * pa2 + & 1738 ( 5.48712484e-01_wp ) * ta * pa2 + & 1739 ( -3.99428410e-03_wp ) * ta2 * pa2 + & 1740 ( -9.54009191e-04_wp ) * ta3 * pa2 + & 1741 ( 1.93090978e-05_wp ) * ta4 * pa2 + & 1742 ( -3.08806365e-01_wp ) * va * pa2 + & 1743 ( 1.16952364e-02_wp ) * ta * va * pa2 + & 1744 ( 4.95271903e-04_wp ) * ta2 * va * pa2 + & 1745 ( -1.90710882e-05_wp ) * ta3 * va * pa2 + & 1746 ( 2.10787756e-03_wp ) * va2 * pa2 + & 1747 ( -6.98445738e-04_wp ) * ta * va2 * pa2 + & 1748 ( 2.30109073e-05_wp ) * ta2 * va2 * pa2 + & 1749 ( 4.17856590e-04_wp ) * va3 * pa2 + & 1750 ( -1.27043871e-05_wp ) * ta * va3 * pa2 + & 1751 ( -3.04620472e-06_wp ) * va4 * pa2 + & 1752 ( 5.14507424e-02_wp ) * d_tmrt * pa2 + & 1753 ( -4.32510997e-03_wp ) * ta * d_tmrt * pa2 + & 1754 ( 8.99281156e-05_wp ) * ta2 * d_tmrt * pa2 + & 1755 ( -7.14663943e-07_wp ) * ta3 * d_tmrt * pa2 + & 1756 ( -2.66016305e-04_wp ) * va * d_tmrt * pa2 + & 1757 ( 2.63789586e-04_wp ) * ta * va * d_tmrt * pa2 + & 1758 ( -7.01199003e-06_wp ) * ta2 * va * d_tmrt * pa2 + & 1759 ( -1.06823306e-04_wp ) * va2 * d_tmrt * pa2 + & 1760 ( 3.61341136e-06_wp ) * ta * va2 * d_tmrt * pa2 + & 1761 ( 2.29748967e-07_wp ) * va3 * d_tmrt * pa2 + & 1762 ( 3.04788893e-04_wp ) * d_tmrt2 * pa2 + & 1763 ( -6.42070836e-05_wp ) * ta * d_tmrt2 * pa2 + & 1764 ( 1.16257971e-06_wp ) * ta2 * d_tmrt2 * pa2 + & 1765 ( 7.68023384e-06_wp ) * va * d_tmrt2 * pa2 + & 1766 ( -5.47446896e-07_wp ) * ta * va * d_tmrt2 * pa2 + & 1767 ( -3.59937910e-08_wp ) * va2 * d_tmrt2 * pa2 + & 1768 ( -4.36497725e-06_wp ) * d_tmrt3 * pa2 + & 1769 ( 1.68737969e-07_wp ) * ta * d_tmrt3 * pa2 + & 1770 ( 2.67489271e-08_wp ) * va * d_tmrt3 * pa2 + & 1771 ( 3.23926897e-09_wp ) * d_tmrt4 * pa2 1772 1773 part_pa3 = ( -3.53874123e-02_wp ) * pa3 + & 1774 ( -2.21201190e-01_wp ) * ta * pa3 + & 1775 ( 1.55126038e-02_wp ) * ta2 * pa3 + & 1776 ( -2.63917279e-04_wp ) * ta3 * pa3 + & 1777 ( 4.53433455e-02_wp ) * va * pa3 + & 1778 ( -4.32943862e-03_wp ) * ta * va * pa3 + & 1779 ( 1.45389826e-04_wp ) * ta2 * va * pa3 + & 1780 ( 2.17508610e-04_wp ) * va2 * pa3 + & 1781 ( -6.66724702e-05_wp ) * ta * va2 * pa3 + & 1782 ( 3.33217140e-05_wp ) * va3 * pa3 + & 1783 ( -2.26921615e-03_wp ) * d_tmrt * pa3 + & 1784 ( 3.80261982e-04_wp ) * ta * d_tmrt * pa3 + & 1785 ( -5.45314314e-09_wp ) * ta2 * d_tmrt * pa3 + & 1786 ( -7.96355448e-04_wp ) * va * d_tmrt * pa3 + & 1787 ( 2.53458034e-05_wp ) * ta * va * d_tmrt * pa3 + & 1788 ( -6.31223658e-06_wp ) * va2 * d_tmrt * pa3 + & 1789 ( 3.02122035e-04_wp ) * d_tmrt2 * pa3 + & 1790 ( -4.77403547e-06_wp ) * ta * d_tmrt2 * pa3 + & 1791 ( 1.73825715e-06_wp ) * va * d_tmrt2 * pa3 + & 1792 ( -4.09087898e-07_wp ) * d_tmrt3 * pa3 1793 1794 part_pa46 = ( 6.14155345e-01_wp ) * pa4 + & 1795 ( -6.16755931e-02_wp ) * ta * pa4 + & 1796 ( 1.33374846e-03_wp ) * ta2 * pa4 + & 1797 ( 3.55375387e-03_wp ) * va * pa4 + & 1798 ( -5.13027851e-04_wp ) * ta * va * pa4 + & 1799 ( 1.02449757e-04_wp ) * va2 * pa4 + & 1800 ( -1.48526421e-03_wp ) * d_tmrt * pa4 + & 1801 ( -4.11469183e-05_wp ) * ta * d_tmrt * pa4 + & 1802 ( -6.80434415e-06_wp ) * va * d_tmrt * pa4 + & 1803 ( -9.77675906e-06_wp ) * d_tmrt2 * pa4 + & 1804 ( 8.82773108e-02_wp ) * pa5 + & 1805 ( -3.01859306e-03_wp ) * ta * pa5 + & 1806 ( 1.04452989e-03_wp ) * va * pa5 + & 1807 ( 2.47090539e-04_wp ) * d_tmrt * pa5 + & 1808 ( 1.48348065e-03_wp ) * pa6 1809 ! 1344 1810 !-- Calculate 6th order polynomial as approximation 1345 utci = ta + & 1346 ( 6.07562052e-01 ) + & 1347 ( -2.27712343e-02 ) * ta + & 1348 ( 8.06470249e-04 ) * ta * ta + & 1349 ( -1.54271372e-04 ) * ta * ta * ta + & 1350 ( -3.24651735e-06 ) * ta * ta * ta * ta + & 1351 ( 7.32602852e-08 ) * ta * ta * ta * ta * ta + & 1352 ( 1.35959073e-09 ) * ta * ta * ta * ta * ta * ta + & 1353 ( -2.25836520e+00 ) * va + & 1354 ( 8.80326035e-02 ) * ta * va + & 1355 ( 2.16844454e-03 ) * ta * ta * va + & 1356 ( -1.53347087e-05 ) * ta * ta * ta * va + & 1357 ( -5.72983704e-07 ) * ta * ta * ta * ta * va + & 1358 ( -2.55090145e-09 ) * ta * ta * ta * ta * ta * va + & 1359 ( -7.51269505e-01 ) * va * va + & 1360 ( -4.08350271e-03 ) * ta * va * va + & 1361 ( -5.21670675e-05 ) * ta * ta * va * va + & 1362 ( 1.94544667e-06 ) * ta * ta * ta * va * va + & 1363 ( 1.14099531e-08 ) * ta * ta * ta * ta * va * va + & 1364 ( 1.58137256e-01 ) * va * va * va + & 1365 ( -6.57263143e-05 ) * ta * va * va * va + & 1366 ( 2.22697524e-07 ) * ta * ta * va * va * va + & 1367 ( -4.16117031e-08 ) * ta * ta * ta * va * va * va + & 1368 ( -1.27762753e-02 ) * va * va * va * va + & 1369 ( 9.66891875e-06 ) * ta * va * va * va * va + & 1370 ( 2.52785852e-09 ) * ta * ta * va * va * va * va + & 1371 ( 4.56306672e-04 ) * va * va * va * va * va + & 1372 ( -1.74202546e-07 ) * ta * va * va * va * va * va + & 1373 ( -5.91491269e-06 ) * va * va * va * va * va * va + & 1374 ( 3.98374029e-01 ) * d_tmrt + & 1375 ( 1.83945314e-04 ) * ta * d_tmrt + & 1376 ( -1.73754510e-04 ) * ta * ta * d_tmrt + & 1377 ( -7.60781159e-07 ) * ta * ta * ta * d_tmrt + & 1378 ( 3.77830287e-08 ) * ta * ta * ta * ta * d_tmrt + & 1379 ( 5.43079673e-10 ) * ta * ta * ta * ta * ta * d_tmrt + & 1380 ( -2.00518269e-02 ) * va * d_tmrt + & 1381 ( 8.92859837e-04 ) * ta * va * d_tmrt + & 1382 ( 3.45433048e-06 ) * ta * ta * va * d_tmrt + & 1383 ( -3.77925774e-07 ) * ta * ta * ta * va * d_tmrt + & 1384 ( -1.69699377e-09 ) * ta * ta * ta * ta * va * d_tmrt + & 1385 ( 1.69992415e-04 ) * va * va * d_tmrt + & 1386 ( -4.99204314e-05 ) * ta * va * va * d_tmrt + & 1387 ( 2.47417178e-07 ) * ta * ta * va * va * d_tmrt + & 1388 ( 1.07596466e-08 ) * ta * ta * ta * va * va * d_tmrt + & 1389 ( 8.49242932e-05 ) * va * va * va * d_tmrt + & 1390 ( 1.35191328e-06 ) * ta * va * va * va * d_tmrt + & 1391 ( -6.21531254e-09 ) * ta * ta * va * va * va * d_tmrt + & 1392 ( -4.99410301e-06 ) * va * va * va * va * d_tmrt + & 1393 ( -1.89489258e-08 ) * ta * va * va * va * va * d_tmrt + & 1394 ( 8.15300114e-08 ) * va * va * va * va * va * d_tmrt + & 1395 ( 7.55043090e-04 ) * d_tmrt * d_tmrt + & 1396 ( -5.65095215e-05 ) * ta * d_tmrt * d_tmrt + & 1397 ( -4.52166564e-07 ) * ta * ta * d_tmrt * d_tmrt + & 1398 ( 2.46688878e-08 ) * ta * ta * ta * d_tmrt * d_tmrt + & 1399 ( 2.42674348e-10 ) * ta * ta * ta * ta * d_tmrt * d_tmrt + & 1400 ( 1.54547250e-04 ) * va * d_tmrt * d_tmrt + & 1401 ( 5.24110970e-06 ) * ta * va * d_tmrt * d_tmrt + & 1402 ( -8.75874982e-08 ) * ta * ta * va * d_tmrt * d_tmrt + & 1403 ( -1.50743064e-09 ) * ta * ta * ta * va * d_tmrt * d_tmrt + & 1404 ( -1.56236307e-05 ) * va * va * d_tmrt * d_tmrt + & 1405 ( -1.33895614e-07 ) * ta * va * va * d_tmrt * d_tmrt + & 1406 ( 2.49709824e-09 ) * ta * ta * va * va * d_tmrt * d_tmrt + & 1407 ( 6.51711721e-07 ) * va * va * va * d_tmrt * d_tmrt + & 1408 ( 1.94960053e-09 ) * ta * va * va * va * d_tmrt * d_tmrt + & 1409 ( -1.00361113e-08 ) * va * va * va * va * d_tmrt * d_tmrt + & 1410 ( -1.21206673e-05 ) * d_tmrt * d_tmrt * d_tmrt + & 1411 ( -2.18203660e-07 ) * ta * d_tmrt * d_tmrt * d_tmrt + & 1412 ( 7.51269482e-09 ) * ta * ta * d_tmrt * d_tmrt * d_tmrt + & 1413 ( 9.79063848e-11 ) * ta * ta * ta * d_tmrt * d_tmrt * d_tmrt + & 1414 ( 1.25006734e-06 ) * va * d_tmrt * d_tmrt * d_tmrt + & 1415 ( -1.81584736e-09 ) * ta * va * d_tmrt * d_tmrt * d_tmrt + & 1416 ( -3.52197671e-10 ) * ta * ta * va * d_tmrt * d_tmrt * d_tmrt + & 1417 ( -3.36514630e-08 ) * va * va * d_tmrt * d_tmrt * d_tmrt + & 1418 ( 1.35908359e-10 ) * ta * va * va * d_tmrt * d_tmrt * d_tmrt + & 1419 ( 4.17032620e-10 ) * va * va * va * d_tmrt * d_tmrt * d_tmrt + & 1420 ( -1.30369025e-09 ) * d_tmrt * d_tmrt * d_tmrt * d_tmrt + & 1421 ( 4.13908461e-10 ) * ta * d_tmrt * d_tmrt * d_tmrt * d_tmrt + & 1422 ( 9.22652254e-12 ) * ta * ta * d_tmrt * d_tmrt * d_tmrt * d_tmrt + & 1423 ( -5.08220384e-09 ) * va * d_tmrt * d_tmrt * d_tmrt * d_tmrt + & 1424 ( -2.24730961e-11 ) * ta * va * d_tmrt * d_tmrt * d_tmrt * d_tmrt + & 1425 ( 1.17139133e-10 ) * va * va * d_tmrt * d_tmrt * d_tmrt * d_tmrt + & 1426 ( 6.62154879e-10 ) * d_tmrt * d_tmrt * d_tmrt * d_tmrt * d_tmrt + & 1427 ( 4.03863260e-13 ) * ta * d_tmrt * d_tmrt * d_tmrt * d_tmrt * d_tmrt + & 1428 ( 1.95087203e-12 ) * va * d_tmrt * d_tmrt * d_tmrt * d_tmrt * d_tmrt + & 1429 ( -4.73602469e-12 ) * d_tmrt * d_tmrt * d_tmrt * d_tmrt * d_tmrt * & 1430 d_tmrt + & 1431 ( 5.12733497e+00 ) * pa + & 1432 ( -3.12788561e-01 ) * ta * pa + & 1433 ( -1.96701861e-02 ) * ta * ta * pa + & 1434 ( 9.99690870e-04 ) * ta * ta * ta * pa + & 1435 ( 9.51738512e-06 ) * ta * ta * ta * ta * pa + & 1436 ( -4.66426341e-07 ) * ta * ta * ta * ta * ta * pa + & 1437 ( 5.48050612e-01 ) * va * pa + & 1438 ( -3.30552823e-03 ) * ta * va * pa + & 1439 ( -1.64119440e-03 ) * ta * ta * va * pa + & 1440 ( -5.16670694e-06 ) * ta * ta * ta * va * pa + & 1441 ( 9.52692432e-07 ) * ta * ta * ta * ta * va * pa + & 1442 ( -4.29223622e-02 ) * va * va * pa + & 1443 ( 5.00845667e-03 ) * ta * va * va * pa + & 1444 ( 1.00601257e-06 ) * ta * ta * va * va * pa + & 1445 ( -1.81748644e-06 ) * ta * ta * ta * va * va * pa + & 1446 ( -1.25813502e-03 ) * va * va * va * pa + & 1447 ( -1.79330391e-04 ) * ta * va * va * va * pa + & 1448 ( 2.34994441e-06 ) * ta * ta * va * va * va * pa + & 1449 ( 1.29735808e-04 ) * va * va * va * va * pa + & 1450 ( 1.29064870e-06 ) * ta * va * va * va * va * pa + & 1451 ( -2.28558686e-06 ) * va * va * va * va * va * pa + & 1452 ( -3.69476348e-02 ) * d_tmrt * pa + & 1453 ( 1.62325322e-03 ) * ta * d_tmrt * pa + & 1454 ( -3.14279680e-05 ) * ta * ta * d_tmrt * pa + & 1455 ( 2.59835559e-06 ) * ta * ta * ta * d_tmrt * pa + & 1456 ( -4.77136523e-08 ) * ta * ta * ta * ta * d_tmrt * pa + & 1457 ( 8.64203390e-03 ) * va * d_tmrt * pa + & 1458 ( -6.87405181e-04 ) * ta * va * d_tmrt * pa + & 1459 ( -9.13863872e-06 ) * ta * ta * va * d_tmrt * pa + & 1460 ( 5.15916806e-07 ) * ta * ta * ta * va * d_tmrt * pa + & 1461 ( -3.59217476e-05 ) * va * va * d_tmrt * pa + & 1462 ( 3.28696511e-05 ) * ta * va * va * d_tmrt * pa + & 1463 ( -7.10542454e-07 ) * ta * ta * va * va * d_tmrt * pa + & 1464 ( -1.24382300e-05 ) * va * va * va * d_tmrt * pa + & 1465 ( -7.38584400e-09 ) * ta * va * va * va * d_tmrt * pa + & 1466 ( 2.20609296e-07 ) * va * va * va * va * d_tmrt * pa + & 1467 ( -7.32469180e-04 ) * d_tmrt * d_tmrt * pa + & 1468 ( -1.87381964e-05 ) * ta * d_tmrt * d_tmrt * pa + & 1469 ( 4.80925239e-06 ) * ta * ta * d_tmrt * d_tmrt * pa + & 1470 ( -8.75492040e-08 ) * ta * ta * ta * d_tmrt * d_tmrt * pa + & 1471 ( 2.77862930e-05 ) * va * d_tmrt * d_tmrt * pa + & 1472 ( -5.06004592e-06 ) * ta * va * d_tmrt * d_tmrt * pa + & 1473 ( 1.14325367e-07 ) * ta * ta * va * d_tmrt * d_tmrt * pa + & 1474 ( 2.53016723e-06 ) * va * va * d_tmrt * d_tmrt * pa + & 1475 ( -1.72857035e-08 ) * ta * va * va * d_tmrt * d_tmrt * pa + & 1476 ( -3.95079398e-08 ) * va * va * va * d_tmrt * d_tmrt * pa + & 1477 ( -3.59413173e-07 ) * d_tmrt * d_tmrt * d_tmrt * pa + & 1478 ( 7.04388046e-07 ) * ta * d_tmrt * d_tmrt * d_tmrt * pa + & 1479 ( -1.89309167e-08 ) * ta * ta * d_tmrt * d_tmrt * d_tmrt * pa + & 1480 ( -4.79768731e-07 ) * va * d_tmrt * d_tmrt * d_tmrt * pa + & 1481 ( 7.96079978e-09 ) * ta * va * d_tmrt * d_tmrt * d_tmrt * pa + & 1482 ( 1.62897058e-09 ) * va * va * d_tmrt * d_tmrt * d_tmrt * pa + & 1483 ( 3.94367674e-08 ) * d_tmrt * d_tmrt * d_tmrt * d_tmrt * pa + & 1484 ( -1.18566247e-09 ) * ta * d_tmrt * d_tmrt * d_tmrt * d_tmrt * pa + & 1485 ( 3.34678041e-10 ) * va * d_tmrt * d_tmrt * d_tmrt * d_tmrt * pa + & 1486 ( -1.15606447e-10 ) * d_tmrt * d_tmrt * d_tmrt * d_tmrt * d_tmrt * pa + & 1487 ( -2.80626406e+00 ) * pa * pa + & 1488 ( 5.48712484e-01 ) * ta * pa * pa + & 1489 ( -3.99428410e-03 ) * ta * ta * pa * pa + & 1490 ( -9.54009191e-04 ) * ta * ta * ta * pa * pa + & 1491 ( 1.93090978e-05 ) * ta * ta * ta * ta * pa * pa + & 1492 ( -3.08806365e-01 ) * va * pa * pa + & 1493 ( 1.16952364e-02 ) * ta * va * pa * pa + & 1494 ( 4.95271903e-04 ) * ta * ta * va * pa * pa + & 1495 ( -1.90710882e-05 ) * ta * ta * ta * va * pa * pa + & 1496 ( 2.10787756e-03 ) * va * va * pa * pa + & 1497 ( -6.98445738e-04 ) * ta * va * va * pa * pa + & 1498 ( 2.30109073e-05 ) * ta * ta * va * va * pa * pa + & 1499 ( 4.17856590e-04 ) * va * va * va * pa * pa + & 1500 ( -1.27043871e-05 ) * ta * va * va * va * pa * pa + & 1501 ( -3.04620472e-06 ) * va * va * va * va * pa * pa + & 1502 ( 5.14507424e-02 ) * d_tmrt * pa * pa + & 1503 ( -4.32510997e-03 ) * ta * d_tmrt * pa * pa + & 1504 ( 8.99281156e-05 ) * ta * ta * d_tmrt * pa * pa + & 1505 ( -7.14663943e-07 ) * ta * ta * ta * d_tmrt * pa * pa + & 1506 ( -2.66016305e-04 ) * va * d_tmrt * pa * pa + & 1507 ( 2.63789586e-04 ) * ta * va * d_tmrt * pa * pa + & 1508 ( -7.01199003e-06 ) * ta * ta * va * d_tmrt * pa * pa + & 1509 ( -1.06823306e-04 ) * va * va * d_tmrt * pa * pa + & 1510 ( 3.61341136e-06 ) * ta * va * va * d_tmrt * pa * pa + & 1511 ( 2.29748967e-07 ) * va * va * va * d_tmrt * pa * pa + & 1512 ( 3.04788893e-04 ) * d_tmrt * d_tmrt * pa * pa + & 1513 ( -6.42070836e-05 ) * ta * d_tmrt * d_tmrt * pa * pa + & 1514 ( 1.16257971e-06 ) * ta * ta * d_tmrt * d_tmrt * pa * pa + & 1515 ( 7.68023384e-06 ) * va * d_tmrt * d_tmrt * pa * pa + & 1516 ( -5.47446896e-07 ) * ta * va * d_tmrt * d_tmrt * pa * pa + & 1517 ( -3.59937910e-08 ) * va * va * d_tmrt * d_tmrt * pa * pa + & 1518 ( -4.36497725e-06 ) * d_tmrt * d_tmrt * d_tmrt * pa * pa + & 1519 ( 1.68737969e-07 ) * ta * d_tmrt * d_tmrt * d_tmrt * pa * pa + & 1520 ( 2.67489271e-08 ) * va * d_tmrt * d_tmrt * d_tmrt * pa * pa + & 1521 ( 3.23926897e-09 ) * d_tmrt * d_tmrt * d_tmrt * d_tmrt * pa * pa + & 1522 ( -3.53874123e-02 ) * pa * pa * pa + & 1523 ( -2.21201190e-01 ) * ta * pa * pa * pa + & 1524 ( 1.55126038e-02 ) * ta * ta * pa * pa * pa + & 1525 ( -2.63917279e-04 ) * ta * ta * ta * pa * pa * pa + & 1526 ( 4.53433455e-02 ) * va * pa * pa * pa + & 1527 ( -4.32943862e-03 ) * ta * va * pa * pa * pa + & 1528 ( 1.45389826e-04 ) * ta * ta * va * pa * pa * pa + & 1529 ( 2.17508610e-04 ) * va * va * pa * pa * pa + & 1530 ( -6.66724702e-05 ) * ta * va * va * pa * pa * pa + & 1531 ( 3.33217140e-05 ) * va * va * va * pa * pa * pa + & 1532 ( -2.26921615e-03 ) * d_tmrt * pa * pa * pa + & 1533 ( 3.80261982e-04 ) * ta * d_tmrt * pa * pa * pa + & 1534 ( -5.45314314e-09 ) * ta * ta * d_tmrt * pa * pa * pa + & 1535 ( -7.96355448e-04 ) * va * d_tmrt * pa * pa * pa + & 1536 ( 2.53458034e-05 ) * ta * va * d_tmrt * pa * pa * pa + & 1537 ( -6.31223658e-06 ) * va * va * d_tmrt * pa * pa * pa + & 1538 ( 3.02122035e-04 ) * d_tmrt * d_tmrt * pa * pa * pa + & 1539 ( -4.77403547e-06 ) * ta * d_tmrt * d_tmrt * pa * pa * pa + & 1540 ( 1.73825715e-06 ) * va * d_tmrt * d_tmrt * pa * pa * pa + & 1541 ( -4.09087898e-07 ) * d_tmrt * d_tmrt * d_tmrt * pa * pa * pa + & 1542 ( 6.14155345e-01 ) * pa * pa * pa * pa + & 1543 ( -6.16755931e-02 ) * ta * pa * pa * pa * pa + & 1544 ( 1.33374846e-03 ) * ta * ta * pa * pa * pa * pa + & 1545 ( 3.55375387e-03 ) * va * pa * pa * pa * pa + & 1546 ( -5.13027851e-04 ) * ta * va * pa * pa * pa * pa + & 1547 ( 1.02449757e-04 ) * va * va * pa * pa * pa * pa + & 1548 ( -1.48526421e-03 ) * d_tmrt * pa * pa * pa * pa + & 1549 ( -4.11469183e-05 ) * ta * d_tmrt * pa * pa * pa * pa + & 1550 ( -6.80434415e-06 ) * va * d_tmrt * pa * pa * pa * pa + & 1551 ( -9.77675906e-06 ) * d_tmrt * d_tmrt * pa * pa * pa * pa + & 1552 ( 8.82773108e-02 ) * pa * pa * pa * pa * pa + & 1553 ( -3.01859306e-03 ) * ta * pa * pa * pa * pa * pa + & 1554 ( 1.04452989e-03 ) * va * pa * pa * pa * pa * pa + & 1555 ( 2.47090539e-04 ) * d_tmrt * pa * pa * pa * pa * pa + & 1556 ( 1.48348065e-03 ) * pa * pa * pa * pa * pa * pa 1557 1811 utci_ij = ta + part_ta + part_va + part_d_tmrt + part_pa + part_pa2 + & 1812 part_pa3 + part_pa46 1813 ! 1558 1814 !-- Consider offset in result 1559 utci = utci+ offset1815 utci_ij = utci_ij + offset 1560 1816 1561 1817 END SUBROUTINE calculate_utci_static … … 1570 1826 !> Value of perct is the Perceived Temperature, degree centigrade 1571 1827 !------------------------------------------------------------------------------! 1572 SUBROUTINE calculate_perct_static( ta, vp, ws, tmrt, pair, clo, perct )1828 SUBROUTINE calculate_perct_static( ta, vp, ws, tmrt, pair, clo, perct_ij ) 1573 1829 1574 1830 IMPLICIT NONE 1575 1831 ! 1576 1832 !-- Type of input of the argument list 1577 1833 REAL(wp), INTENT ( IN ) :: ta !< Local air temperature (degC) … … 1580 1836 REAL(wp), INTENT ( IN ) :: ws !< Local wind velocitry (m/s) 1581 1837 REAL(wp), INTENT ( IN ) :: pair !< Local barometric air pressure (hPa) 1582 1838 ! 1583 1839 !-- Type of output of the argument list 1584 REAL(wp), INTENT ( OUT ) :: perct 1585 REAL(wp), INTENT ( OUT ) :: clo !< Clothing index (dimensionless)1586 1840 REAL(wp), INTENT ( OUT ) :: perct_ij !< Perceived temperature (degC) 1841 REAL(wp), INTENT ( OUT ) :: clo !< Clothing index (dimensionless) 1842 ! 1587 1843 !-- Parameters for standard "Klima-Michel" 1588 1844 REAL(wp), PARAMETER :: eta = 0._wp !< Mechanical work efficiency for walking on flat ground (compare to Fanger (1972) pp 24f) 1589 1845 REAL(wp), PARAMETER :: actlev = 134.6862_wp !< Workload by activity per standardized surface (A_Du) 1590 1846 ! 1591 1847 !-- Type of program variables 1592 1848 REAL(wp), PARAMETER :: eps = 0.0005 !< Accuracy in clothing insulation (clo) for evaluation the root of Fanger's PMV (pmva=0) … … 1614 1870 1615 1871 LOGICAL :: sultrieness 1616 1872 ! 1617 1873 !-- Initialise 1618 perct = 9999.0_wp1874 perct_ij = bio_fill_value 1619 1875 1620 1876 nerr = 0_iwp 1621 1877 ncount = 0_iwp 1622 1878 sultrieness = .FALSE. 1879 ! 1623 1880 !-- Tresholds: clothing insulation (account for model inaccuracies) 1881 ! 1624 1882 ! summer clothing 1625 1883 sclo = 0.44453_wp 1884 ! 1626 1885 ! winter clothing 1627 1886 wclo = 1.76267_wp 1628 1887 ! 1629 1888 !-- decision: firstly calculate for winter or summer clothing 1630 1889 IF ( ta <= 10._wp ) THEN 1631 1890 ! 1632 1891 !-- First guess: winter clothing insulation: cold stress 1633 1892 clo = wclo … … 1636 1895 1637 1896 IF ( pmva > 0._wp ) THEN 1638 1897 ! 1639 1898 !-- Case summer clothing insulation: heat load ? 1640 1899 clo = sclo … … 1643 1902 pmv_s = pmva 1644 1903 IF ( pmva <= 0._wp ) THEN 1904 ! 1645 1905 !-- Case: comfort achievable by varying clothing insulation 1646 ! --Between winter and summer set values1906 ! Between winter and summer set values 1647 1907 CALL iso_ridder ( ta, tmrt, vp, ws, pair, actlev, eta, sclo, & 1648 1908 pmv_s, wclo, pmv_w, eps, pmva, top, ncount, clo ) … … 1662 1922 ENDIF 1663 1923 ELSE 1664 1924 ! 1665 1925 !-- First guess: summer clothing insulation: heat load 1666 1926 clo = sclo … … 1669 1929 1670 1930 IF ( pmva < 0._wp ) THEN 1671 1931 ! 1672 1932 !-- Case winter clothing insulation: cold stress ? 1673 1933 clo = wclo … … 1677 1937 1678 1938 IF ( pmva >= 0._wp ) THEN 1679 1939 ! 1680 1940 !-- Case: comfort achievable by varying clothing insulation 1681 1941 ! between winter and summer set values … … 1698 1958 1699 1959 ENDIF 1700 1960 ! 1701 1961 !-- Determine perceived temperature by regression equation + adjustments 1702 1962 pmvs = pmva 1703 CALL perct_regression ( pmva, clo, perct )1704 ptc = perct 1963 CALL perct_regression ( pmva, clo, perct_ij ) 1964 ptc = perct_ij 1705 1965 IF ( clo >= 1.75_wp .AND. pmva <= -0.11_wp ) THEN 1966 ! 1706 1967 !-- Adjust for cold conditions according to Gagge 1986 1707 1968 CALL dpmv_cold ( pmva, ta, ws, tmrt, nerr_cold, d_pmv ) … … 1712 1973 pmvs = -0.11_wp 1713 1974 ENDIF 1714 CALL perct_regression ( pmvs, clo, perct )1975 CALL perct_regression ( pmvs, clo, perct_ij ) 1715 1976 ENDIF 1716 1977 ! clo_fanger = clo 1717 1978 clon = clo 1718 IF ( clo > 0.5_wp .AND. perct <= 8.73_wp ) THEN 1979 IF ( clo > 0.5_wp .AND. perct_ij <= 8.73_wp ) THEN 1980 ! 1719 1981 !-- Required clothing insulation (ireq) is exclusively defined for 1720 1982 ! operative temperatures (top) less 10 (C) for a 1721 1983 ! reference wind of 0.2 m/s according to 8.73 (C) for 0.1 m/s 1722 clon = ireq_neutral ( perct , ireq_minimal, nerr )1984 clon = ireq_neutral ( perct_ij, ireq_minimal, nerr ) 1723 1985 clo = clon 1724 1986 ENDIF … … 1727 1989 d_std = -99._wp 1728 1990 IF ( pmva > 0.06_wp .AND. clo <= 0.5_wp ) THEN 1991 ! 1729 1992 !-- Adjust for warm/humid conditions according to Gagge 1986 1730 1993 CALL saturation_vapor_pressure ( ta, svp_ta ) 1731 1994 d_pmv = deltapmv ( pmva, ta, vp, svp_ta, tmrt, ws, nerr ) 1732 1995 pmvs = pmva + d_pmv 1733 CALL perct_regression ( pmvs, clo, perct )1996 CALL perct_regression ( pmvs, clo, perct_ij ) 1734 1997 IF ( sult_lim < 99._wp ) THEN 1735 IF ( (perct - ptc) > sult_lim ) sultrieness = .TRUE. 1998 IF ( (perct_ij - ptc) > sult_lim ) sultrieness = .TRUE. 1999 ! 1736 2000 !-- Set factor to threshold for sultriness 1737 2001 IF ( dgtcstd /= 0_iwp ) THEN 1738 d_std = ( ( perct - ptc ) - dgtcm ) / dgtcstd2002 d_std = ( ( perct_ij - ptc ) - dgtcm ) / dgtcstd 1739 2003 ENDIF 1740 2004 ENDIF … … 1762 2026 1763 2027 IF ( ta < 0._wp ) THEN 2028 ! 1764 2029 !-- ta < 0 (degC): saturation water vapour pressure over ice 1765 2030 b = 17.84362_wp 1766 2031 c = 245.425_wp 1767 2032 ELSE 2033 ! 1768 2034 !-- ta >= 0 (degC): saturation water vapour pressure over water 1769 2035 b = 17.08085_wp 1770 2036 c = 234.175_wp 1771 2037 ENDIF 1772 2038 ! 1773 2039 !-- Saturation water vapour pressure 1774 2040 svp_ta = 6.1078_wp * EXP ( b * ta / ( c + ta ) ) … … 1788 2054 1789 2055 IMPLICIT NONE 1790 2056 ! 1791 2057 !-- Input variables of argument list: 1792 2058 REAL(wp), INTENT ( IN ) :: ta !< Ambient temperature (degC) … … 1803 2069 REAL(wp), INTENT ( IN ) :: pmv_w !< Fanger's PMV corresponding to wclo 1804 2070 REAL(wp), INTENT ( IN ) :: pmv_s !< Fanger's PMV corresponding to sclo 1805 2071 ! 1806 2072 ! Output variables of argument list: 1807 2073 REAL(wp), INTENT ( OUT ) :: pmva !< 0 (set to zero, because clo is evaluated for comfort) … … 1814 2080 ! nerr = -2: error = maximum iterations (max_iteration) exceeded 1815 2081 ! nerr = -3: error = root not bracketed between sclo and wclo 1816 2082 ! 1817 2083 !-- Type of program variables 1818 2084 INTEGER(iwp), PARAMETER :: max_iteration = 15_iwp !< max number of iterations … … 1831 2097 REAL(wp) :: sroot !< sqrt of PMV-guess 1832 2098 INTEGER(iwp) :: j !< running index 1833 2099 ! 1834 2100 !-- Initialise 1835 2101 nerr = 0_iwp 1836 2102 ! 1837 2103 !-- Set pmva = 0 (comfort): Root of PMV depending on clothing insulation 2104 x_ridder = bio_fill_value 1838 2105 pmva = 0._wp 1839 2106 clo_lower = sclo … … 1884 2151 y_lower = y_new 1885 2152 ELSE 2153 ! 1886 2154 !-- Never get here in x_ridder: singularity in y 1887 2155 nerr = -1_iwp … … 1895 2163 ENDIF 1896 2164 ENDDO 2165 ! 1897 2166 !-- x_ridder exceed maximum iterations 1898 2167 nerr = -2_iwp … … 1904 2173 x_ridder = clo_upper 1905 2174 ELSE 2175 ! 1906 2176 !-- x_ridder not bracketed by u_clo and o_clo 1907 2177 nerr = -3_iwp … … 1919 2189 !> individual and reference environment. 1920 2190 !------------------------------------------------------------------------------! 1921 SUBROUTINE perct_regression( pmv, clo, perct )2191 SUBROUTINE perct_regression( pmv, clo, perct_ij ) 1922 2192 1923 2193 IMPLICIT NONE … … 1926 2196 REAL(wp), INTENT ( IN ) :: clo !< clothing insulation index (clo) 1927 2197 1928 REAL(wp), INTENT ( OUT ) :: perct !< perct (degC) corresponding to given PMV / clo2198 REAL(wp), INTENT ( OUT ) :: perct_ij !< perct (degC) corresponding to given PMV / clo 1929 2199 1930 2200 IF ( pmv <= -0.11_wp ) THEN 1931 perct = 5.805_wp + 12.6784_wp * pmv2201 perct_ij = 5.805_wp + 12.6784_wp * pmv 1932 2202 ELSE 1933 2203 IF ( pmv >= + 0.01_wp ) THEN 1934 perct = 16.826_wp + 6.163_wp * pmv2204 perct_ij = 16.826_wp + 6.163_wp * pmv 1935 2205 ELSE 1936 perct = 21.258_wp - 9.558_wp * clo2206 perct_ij = 21.258_wp - 9.558_wp * clo 1937 2207 ENDIF 1938 2208 ENDIF … … 1953 2223 1954 2224 IMPLICIT NONE 1955 2225 ! 1956 2226 !-- Input variables of argument list: 1957 2227 REAL(wp), INTENT ( IN ) :: ta !< Ambient air temperature (degC) … … 1963 2233 REAL(wp), INTENT ( IN ) :: actlev !< Individuals activity level per unit surface area (W/m2) 1964 2234 REAL(wp), INTENT ( IN ) :: eta !< Individuals mechanical work efficiency (dimensionless) 1965 2235 ! 1966 2236 !-- Output variables of argument list: 1967 2237 REAL(wp), INTENT ( OUT ) :: pmva !< Actual Predicted Mean Vote (PMV, … … 1969 2239 ! (ta,tmrt,pa,ws,pair) and individual variables (clo, actlev, eta) 1970 2240 REAL(wp), INTENT ( OUT ) :: top !< operative temperature (degC) 1971 2241 ! 1972 2242 !-- Internal variables 1973 2243 REAL(wp) :: f_cl !< Increase in surface due to clothing (factor) … … 1992 2262 REAL(wp) :: z6 !< Heat loss through forced convection 1993 2263 INTEGER(iwp) :: i !< running index 1994 2264 ! 1995 2265 !-- Clo must be > 0. to avoid div. by 0! 1996 2266 clo = in_clo 1997 2267 IF ( clo <= 0._wp ) clo = .001_wp 1998 2268 ! 1999 2269 !-- f_cl = Increase in surface due to clothing 2000 2270 f_cl = 1._wp + .15_wp * clo 2001 2271 ! 2002 2272 !-- Case of free convection (ws < 0.1 m/s ) not considered 2003 2273 ws = in_ws … … 2005 2275 ws = .1_wp 2006 2276 ENDIF 2007 2277 ! 2008 2278 !-- Heat_convection = forced convection 2009 2279 heat_convection = 12.1_wp * SQRT ( ws * pair / 1013.25_wp ) 2010 2280 ! 2011 2281 !-- Activity = inner heat produktion per standardized surface 2012 2282 activity = actlev * ( 1._wp - eta ) 2013 2283 ! 2014 2284 !-- T_skin_aver = average skin temperature 2015 2285 t_skin_aver = 35.7_wp - .0275_wp * activity 2016 2286 ! 2017 2287 !-- Calculation of constants for evaluation below 2018 2288 bc = .155_wp * clo * 3.96_wp * 10._wp**( -8 ) * f_cl … … 2021 2291 dc = ( 1._wp + ec * cc ) / bc 2022 2292 gc = ( t_skin_aver + bc * ( tmrt + degc_to_k )**4 + ec * cc * ta ) / bc 2023 2293 ! 2024 2294 !-- Calculation of clothing surface temperature (t_clothing) based on 2025 2295 ! Newton-approximation with air temperature as initial guess … … 2029 2299 * dc - gc ) / ( 4._wp * ( t_clothing + degc_to_k )**3 + dc ) 2030 2300 ENDDO 2031 2301 ! 2032 2302 !-- Empiric factor for the adaption of the heat ballance equation 2033 2303 ! to the psycho-physical scale (Equ. 40 in FANGER) 2034 2304 z1 = ( .303_wp * EXP ( -.036_wp * actlev ) + .0275_wp ) 2035 2305 ! 2036 2306 !-- Water vapour diffution through the skin 2037 2307 z2 = .31_wp * ( 57.3_wp - .07_wp * activity-pa ) 2038 2308 ! 2039 2309 !-- Sweat evaporation from the skin surface 2040 2310 z3 = .42_wp * ( activity - 58._wp ) 2041 2311 ! 2042 2312 !-- Loss of latent heat through respiration 2043 2313 z4 = .0017_wp * actlev * ( 58.7_wp - pa ) + .0014_wp * actlev * & 2044 2314 ( 34._wp - ta ) 2045 2315 ! 2046 2316 !-- Loss of radiational heat 2047 2317 z5 = 3.96e-8_wp * f_cl * ( ( t_clothing + degc_to_k )**4 - ( tmrt + & … … 2052 2322 hr = 0._wp 2053 2323 ENDIF 2054 2324 ! 2055 2325 !-- Heat loss through forced convection cc*(t_clothing-TT) 2056 2326 z6 = cc * ( t_clothing - ta ) 2057 2327 ! 2058 2328 !-- Predicted Mean Vote 2059 2329 pmva = z1 * ( activity - z2 - z3 - z4 - z5 - z6 ) 2060 2330 ! 2061 2331 !-- Operative temperatur 2062 2332 top = ( hr * tmrt + heat_convection * ta ) / ( hr + heat_convection ) … … 2073 2343 2074 2344 IMPLICIT NONE 2075 2345 ! 2076 2346 !-- Input variables of argument list: 2077 2347 REAL(wp), INTENT ( IN ) :: pmva !< Actual Predicted Mean Vote (PMV) according to Fanger … … 2081 2351 REAL(wp), INTENT ( IN ) :: tmrt !< Mean radiant temperature (degC) at screen level 2082 2352 REAL(wp), INTENT ( IN ) :: ws !< Wind speed (m/s) 1 m above ground 2083 2353 ! 2084 2354 !-- Output variables of argument list: 2085 2355 INTEGER(iwp), INTENT ( OUT ) :: nerr !< Error status / quality flag … … 2088 2358 ! -3 = rel. humidity set to 5 % or 95 %, respectively 2089 2359 ! -4 = deltapmv set to avoid pmvs < 0 2090 2360 ! 2091 2361 !-- Internal variable types: 2092 2362 REAL(wp) :: pmv !< … … 2142 2412 +1.8686215_wp, +3.4260713_wp, +2.0116185_wp, -0.7777552_wp, -4.6715853_wp,& 2143 2413 -7.7314281_wp, -11.7602578_wp, -23.5934198_wp / 2144 2414 ! 2145 2415 !-- Test for compliance with regression range 2146 2416 IF ( pmva < -1.0_wp .OR. pmva > 7.0_wp ) THEN … … 2149 2419 nerr = 0_iwp 2150 2420 ENDIF 2151 2421 ! 2152 2422 !-- Initialise classic PMV 2153 2423 pmv = pmva 2154 2424 ! 2155 2425 !-- Water vapour pressure of air 2156 2426 p10 = 0.05_wp * svp_ta … … 2161 2431 nerr = -3_iwp 2162 2432 IF ( vp < p10 ) THEN 2433 ! 2163 2434 !-- Due to conditions of regression: r.H. >= 5 % 2164 2435 pa = p10 2165 2436 ELSE 2437 ! 2166 2438 !-- Due to conditions of regression: r.H. <= 95 % 2167 2439 pa = p95 … … 2169 2441 ENDIF 2170 2442 IF ( pa > 0._wp ) THEN 2443 ! 2171 2444 !-- Natural logarithm of pa 2172 2445 apa = LOG ( pa ) … … 2174 2447 apa = -5._wp 2175 2448 ENDIF 2176 2449 ! 2177 2450 !-- Ratio actual water vapour pressure to that of a r.H. of 50 % 2178 2451 pa_p50 = 0.5_wp * svp_ta … … 2184 2457 pa_p50 = 0._wp 2185 2458 ENDIF 2459 ! 2186 2460 !-- Square root of wind velocity 2187 2461 IF ( ws >= 0._wp ) THEN … … 2190 2464 sqvel = 0._wp 2191 2465 ENDIF 2466 ! 2192 2467 !-- Air temperature 2193 2468 ! ta = ta 2194 2469 !-- Difference mean radiation to air temperature 2195 2470 dtmrt = tmrt - ta 2196 2471 ! 2197 2472 !-- Select the valid regression coefficients 2198 2473 nreg = INT ( pmv ) 2199 2474 IF ( nreg < 0_iwp ) THEN 2475 ! 2200 2476 !-- value of the FUNCTION in the case pmv <= -1 2201 2477 deltapmv = 0._wp … … 2213 2489 ENDIF 2214 2490 ENDIF 2215 2491 ! 2216 2492 !-- Regression valid for 0. <= pmv <= 6. 2217 2493 dpmv_1 = & … … 2239 2515 + aconst ( nreg + 1_iwp ) 2240 2516 ENDIF 2241 2517 ! 2242 2518 !-- Calculate pmv modification 2243 2519 deltapmv = ( 1._wp - gew ) * dpmv_1 + gew * dpmv_2 2244 2520 pmvs = pmva + deltapmv 2245 2521 IF ( ( pmvs ) < 0._wp ) THEN 2522 ! 2246 2523 !-- Prevent negative pmv* due to problems with clothing insulation 2247 2524 nerr = -4_iwp 2248 2525 IF ( pmvs > -0.11_wp ) THEN 2526 ! 2249 2527 !-- Threshold from FUNCTION perct_regression for winter clothing insulation 2250 2528 deltapmv = deltapmv + 0.11_wp 2251 2529 ELSE 2530 ! 2252 2531 !-- Set pmvs to "0" for compliance with summer clothing insulation 2253 2532 deltapmv = -1._wp * pmva … … 2266 2545 !> perct. 2267 2546 !------------------------------------------------------------------------------! 2268 SUBROUTINE calc_sultr( perct , dperctm, dperctstd, sultr_res )2547 SUBROUTINE calc_sultr( perct_ij, dperctm, dperctstd, sultr_res ) 2269 2548 2270 2549 IMPLICIT NONE 2271 2550 ! 2272 2551 !-- Input of the argument list: 2273 REAL(wp), INTENT ( IN ) :: perct !< Classical perceived temperature: Base is Fanger's PMV2274 2552 REAL(wp), INTENT ( IN ) :: perct_ij !< Classical perceived temperature: Base is Fanger's PMV 2553 ! 2275 2554 !-- Additional output variables of argument list: 2276 2555 REAL(wp), INTENT ( OUT ) :: dperctm !< Mean deviation perct (classical gt) to gt* (rational gt … … 2279 2558 ! determining the significance to perceive sultriness 2280 2559 REAL(wp), INTENT ( OUT ) :: sultr_res 2281 2560 ! 2282 2561 !-- Types of coefficients mean deviation: third order polynomial 2283 2562 REAL(wp), PARAMETER :: dperctka = +7.5776086_wp … … 2285 2564 REAL(wp), PARAMETER :: dperctkc = +0.0213324_wp 2286 2565 REAL(wp), PARAMETER :: dperctkd = -0.00027797237_wp 2287 2566 ! 2288 2567 !-- Types of coefficients mean deviation plus standard deviation 2289 2568 ! regression coefficients: third order polynomial … … 2292 2571 REAL(wp), PARAMETER :: dperctsc = -0.00054709752_wp 2293 2572 REAL(wp), PARAMETER :: dperctsd = +0.0000063714823_wp 2294 2573 ! 2295 2574 !-- Factor to mean standard deviation defining SIGNificance for 2296 2575 ! sultriness 2297 2576 REAL(wp), PARAMETER :: faktor = 1._wp 2298 2577 ! 2299 2578 !-- Initialise 2300 2579 sultr_res = +99._wp … … 2302 2581 dperctstd = 999999._wp 2303 2582 2304 IF ( perct < 16.826_wp .OR. perct > 56._wp ) THEN 2305 !-- Unallowed classical PMV/perct 2583 IF ( perct_ij < 16.826_wp .OR. perct_ij > 56._wp ) THEN 2584 ! 2585 !-- Unallowed value of classical perct! 2306 2586 RETURN 2307 2587 ENDIF 2308 2588 ! 2309 2589 !-- Mean deviation dependent on perct 2310 dperctm = dperctka + dperctkb * perct + dperctkc * perct**2._wp + dperctkd * perct**3._wp 2311 2590 dperctm = dperctka + dperctkb * perct_ij + dperctkc * perct_ij**2._wp + & 2591 dperctkd * perct_ij**3._wp 2592 ! 2312 2593 !-- Mean deviation plus its standard deviation 2313 dperctstd = dperctsa + dperctsb * perct + dperctsc * perct**2._wp + dperctsd * perct**3._wp 2314 2594 dperctstd = dperctsa + dperctsb * perct_ij + dperctsc * perct_ij**2._wp + & 2595 dperctsd * perct_ij**3._wp 2596 ! 2315 2597 !-- Value of the FUNCTION 2316 2598 sultr_res = dperctm + faktor * dperctstd … … 2330 2612 2331 2613 IMPLICIT NONE 2332 2614 ! 2333 2615 !-- Type of input arguments 2334 2616 REAL(wp), INTENT ( IN ) :: pmva !< Fanger's classical predicted mean vote … … 2336 2618 REAL(wp), INTENT ( IN ) :: ws !< Relative wind velocity 1 m above ground (m/s) 2337 2619 REAL(wp), INTENT ( IN ) :: tmrt !< Mean radiant temperature (degC) 2338 2620 ! 2339 2621 !-- Type of output argument 2340 2622 INTEGER(iwp), INTENT ( OUT ) :: nerr !< Error indicator: 0 = o.k., +1 = denominator for intersection = 0 2341 2623 REAL(wp), INTENT ( OUT ) :: dpmv_cold_res !< Increment to adjust pmva according to the results of Gagge's 2342 2624 ! 2 node model depending on the input 2343 2625 ! 2344 2626 !-- Type of program variables 2345 2627 REAL(wp) :: delta_cold(3) … … 2352 2634 REAL(wp) :: sqrt_ws 2353 2635 INTEGER(iwp) :: i 2354 INTEGER(iwp) :: j2636 ! INTEGER(iwp) :: j 2355 2637 INTEGER(iwp) :: i_bin 2356 2638 ! 2357 2639 !-- Coefficient of the 3 regression lines 2358 2640 ! 1:const 2:*pmvc 3:*ta 4:*sqrt_ws 5:*dtmrt 2359 2641 coeff(1,1:5) = & 2360 2642 (/ +0.161_wp, +0.130_wp, -1.125E-03_wp, +1.106E-03_wp, -4.570E-04_wp /) … … 2363 2645 coeff(3,1:5) = & 2364 2646 (/ +0.05761_wp, +0.458_wp, -1.829E-02_wp, -5.577E-03_wp, -1.970E-03_wp /) 2365 2647 ! 2366 2648 !-- Initialise 2367 2649 nerr = 0_iwp … … 2384 2666 2385 2667 DO i = 1, 3 2668 ! 2386 2669 !-- Regression constant for given meteorological variables 2387 reg_a(i) = coeff(i, 1) + coeff(i, 3) * ta + coeff(i, 4) * &2670 reg_a(i) = coeff(i, 1) + coeff(i, 3) * ta + coeff(i, 4) * & 2388 2671 sqrt_ws + coeff(i,5)*dtmrt 2389 2672 delta_cold(i) = reg_a(i) + coeff(i, 2) * pmvc 2390 2673 ENDDO 2391 2674 ! 2392 2675 !-- Intersection points of regression lines in terms of Fanger's PMV 2393 2676 DO i = 1, 2 … … 2408 2691 ENDIF 2409 2692 ENDDO 2693 ! 2410 2694 !-- Adjust to operative temperature scaled according 2411 2695 ! to classical PMV (Fanger) … … 2434 2718 INTEGER(iwp) :: i, i_bin 2435 2719 2436 ! range_1 range_2range_32720 ! range_1 range_2 range_3 2437 2721 DATA (coef(i, 0), i = 1, n_bin) /0.0941540_wp, -0.1506620_wp, -0.0871439_wp/ 2438 2722 DATA (coef(i, 1), i = 1, n_bin) /0.0783162_wp, -1.0612651_wp, 0.1695040_wp/ … … 2467 2751 !> for perct < 10 (degC) 2468 2752 !------------------------------------------------------------------------------! 2469 REAL(wp) FUNCTION ireq_neutral( perct , ireq_minimal, nerr )2753 REAL(wp) FUNCTION ireq_neutral( perct_ij, ireq_minimal, nerr ) 2470 2754 2471 2755 IMPLICIT NONE 2472 2756 ! 2473 2757 !-- Type declaration of arguments 2474 REAL(wp), INTENT ( IN ) :: perct 2758 REAL(wp), INTENT ( IN ) :: perct_ij 2475 2759 REAL(wp), INTENT ( OUT ) :: ireq_minimal 2476 2760 INTEGER(iwp), INTENT ( OUT ) :: nerr 2477 2761 ! 2478 2762 !-- Type declaration for internal varables 2479 2763 REAL(wp) :: top02 2480 2764 ! 2481 2765 !-- Initialise 2482 2766 nerr = 0_iwp 2483 2767 ! 2484 2768 !-- Convert perceived temperature from basis 0.1 m/s to basis 0.2 m/s 2485 top02 = 1.8788_wp + 0.9296_wp * perct 2486 2769 top02 = 1.8788_wp + 0.9296_wp * perct_ij 2770 ! 2487 2771 !-- IREQ neutral conditions (thermal comfort) 2488 2772 ireq_neutral = 1.62_wp - 0.0564_wp * top02 2489 2773 ! 2490 2774 !-- Regression only defined for perct <= 10 (degC) 2491 2775 IF ( ireq_neutral < 0.5_wp ) THEN … … 2495 2779 ireq_neutral = 0.5_wp 2496 2780 ENDIF 2497 2781 ! 2498 2782 !-- Minimal required clothing insulation: maximal acceptable body cooling 2499 2783 ireq_minimal = 1.26_wp - 0.0588_wp * top02 … … 2523 2807 2524 2808 height = height_cm * 100._wp 2525 2809 ! 2526 2810 !-- According to Gehan-George, for children 2527 2811 IF ( age < 19_iwp ) THEN … … 2533 2817 RETURN 2534 2818 END IF 2535 2819 ! 2536 2820 !-- DuBois D, DuBois EF: A formula to estimate the approximate surface area if 2537 2821 ! height and weight be known. In: Arch. Int. Med.. 17, 1916, S. 863?871. … … 2598 2882 2599 2883 SUBROUTINE ipt_init (age, weight, height, sex, work, actlev, clo, & 2600 ta, vp, ws, tmrt, pair, dt, storage, t_clothing, &2884 ta, vp, ws, tmrt, pair, dt, storage, t_clothing, & 2601 2885 ipt ) 2602 2886 2603 2887 IMPLICIT NONE 2604 2888 ! 2605 2889 !-- Input parameters 2606 2890 REAL(wp), INTENT(in) :: age !< Persons age (years) … … 2615 2899 REAL(wp), INTENT(in) :: dt !< Timestep (s) 2616 2900 INTEGER(iwp), INTENT(in) :: sex !< Persons sex (1 = male, 2 = female) 2617 2901 ! 2618 2902 !-- Output parameters 2619 2903 REAL(wp), INTENT(out) :: actlev … … 2622 2906 REAL(wp), INTENT(out) :: t_clothing 2623 2907 REAL(wp), INTENT(out) :: ipt 2624 2908 ! 2625 2909 !-- Internal variables 2626 2910 REAL(wp), PARAMETER :: eps = 0.0005_wp … … 2644 2928 REAL(wp) :: top 2645 2929 REAL(wp) :: a_surf 2646 REAL(wp) :: acti2930 ! REAL(wp) :: acti 2647 2931 INTEGER(iwp) :: ncount 2648 2932 INTEGER(iwp) :: nerr_cold … … 2653 2937 storage = 0._wp 2654 2938 CALL persdat ( age, weight, height, sex, work, a_surf, actlev ) 2655 2939 ! 2656 2940 !-- Initialise 2657 t_clothing = -999.0_wp2658 ipt = -999.0_wp2941 t_clothing = bio_fill_value 2942 ipt = bio_fill_value 2659 2943 nerr = 0_wp 2660 2944 ncount = 0_wp 2661 2945 sultrieness = .FALSE. 2946 ! 2662 2947 !-- Tresholds: clothing insulation (account for model inaccuracies) 2663 2948 !-- Summer clothing … … 2665 2950 !-- Winter clothing 2666 2951 wclo = 1.76267_wp 2667 2952 ! 2668 2953 !-- Decision: firstly calculate for winter or summer clothing 2669 2954 IF ( ta <= 10._wp ) THEN 2670 2955 ! 2671 2956 !-- First guess: winter clothing insulation: cold stress 2672 2957 clo = wclo 2673 t_clothing = -999.0_wp! force initial run2958 t_clothing = bio_fill_value ! force initial run 2674 2959 CALL fanger_s_acti ( ta, tmrt, vp, ws, pair, clo, actlev, work, & 2675 2960 t_clothing, storage, dt, pmva ) … … 2677 2962 2678 2963 IF ( pmva > 0._wp ) THEN 2679 2964 ! 2680 2965 !-- Case summer clothing insulation: heat load ? 2681 2966 clo = sclo 2682 t_clothing = -999.0_wp! force initial run2967 t_clothing = bio_fill_value ! force initial run 2683 2968 CALL fanger_s_acti ( ta, tmrt, vp, ws, pair, clo, actlev, work, & 2684 2969 t_clothing, storage, dt, pmva ) 2685 2970 pmv_s = pmva 2686 2971 IF ( pmva <= 0._wp ) THEN 2972 ! 2687 2973 !-- Case: comfort achievable by varying clothing insulation 2688 2974 ! between winter and summer set values … … 2695 2981 ELSE IF ( pmva > 0.06_wp ) THEN 2696 2982 clo = 0.5_wp 2697 t_clothing = -999.0_wp2983 t_clothing = bio_fill_value 2698 2984 CALL fanger_s_acti ( ta, tmrt, vp, ws, pair, clo, actlev, work, & 2699 2985 t_clothing, storage, dt, pmva ) … … 2701 2987 ELSE IF ( pmva < -0.11_wp ) THEN 2702 2988 clo = 1.75_wp 2703 t_clothing = -999.0_wp2989 t_clothing = bio_fill_value 2704 2990 CALL fanger_s_acti ( ta, tmrt, vp, ws, pair, clo, actlev, work, & 2705 2991 t_clothing, storage, dt, pmva ) … … 2707 2993 2708 2994 ELSE 2709 2995 ! 2710 2996 !-- First guess: summer clothing insulation: heat load 2711 2997 clo = sclo 2712 t_clothing = -999.0_wp2998 t_clothing = bio_fill_value 2713 2999 CALL fanger_s_acti ( ta, tmrt, vp, ws, pair, clo, actlev, work, & 2714 3000 t_clothing, storage, dt, pmva ) … … 2716 3002 2717 3003 IF ( pmva < 0._wp ) THEN 2718 3004 ! 2719 3005 !-- Case winter clothing insulation: cold stress ? 2720 3006 clo = wclo 2721 t_clothing = -999.0_wp3007 t_clothing = bio_fill_value 2722 3008 CALL fanger_s_acti ( ta, tmrt, vp, ws, pair, clo, actlev, work, & 2723 3009 t_clothing, storage, dt, pmva ) … … 2725 3011 2726 3012 IF ( pmva >= 0._wp ) THEN 2727 3013 ! 2728 3014 !-- Case: comfort achievable by varying clothing insulation 2729 3015 ! between winter and summer set values … … 2736 3022 ELSE IF ( pmva < -0.11_wp ) THEN 2737 3023 clo = 1.75_wp 2738 t_clothing = -999.0_wp3024 t_clothing = bio_fill_value 2739 3025 CALL fanger_s_acti ( ta, tmrt, vp, ws, pair, clo, actlev, work, & 2740 3026 t_clothing, storage, dt, pmva ) … … 2742 3028 ELSE IF ( pmva > 0.06_wp ) THEN 2743 3029 clo = 0.5_wp 2744 t_clothing = -999.0_wp3030 t_clothing = bio_fill_value 2745 3031 CALL fanger_s_acti ( ta, tmrt, vp, ws, pair, clo, actlev, work, & 2746 3032 t_clothing, storage, dt, pmva ) … … 2748 3034 2749 3035 ENDIF 2750 3036 ! 2751 3037 !-- Determine perceived temperature by regression equation + adjustments 2752 3038 pmvs = pmva … … 2754 3040 ptc = ipt 2755 3041 IF ( clo >= 1.75_wp .AND. pmva <= -0.11_wp ) THEN 3042 ! 2756 3043 !-- Adjust for cold conditions according to Gagge 1986 2757 3044 CALL dpmv_cold ( pmva, ta, ws, tmrt, nerr_cold, d_pmv ) … … 2767 3054 clon = clo 2768 3055 IF ( clo > 0.5_wp .AND. ipt <= 8.73_wp ) THEN 3056 ! 2769 3057 !-- Required clothing insulation (ireq) is exclusively defined for 2770 3058 ! operative temperatures (top) less 10 (C) for a … … 2777 3065 d_std = -99._wp 2778 3066 IF ( pmva > 0.06_wp .AND. clo <= 0.5_wp ) THEN 3067 ! 2779 3068 !-- Adjust for warm/humid conditions according to Gagge 1986 2780 3069 CALL saturation_vapor_pressure ( ta, svp_ta ) … … 2808 3097 2809 3098 IMPLICIT NONE 2810 3099 ! 2811 3100 !-- Type of input of the argument list 2812 3101 REAL(wp), INTENT ( IN ) :: ta !< Air temperature (°C) … … 2819 3108 REAL(wp), INTENT ( IN ) :: actlev !< Internal heat production (W) 2820 3109 REAL(wp), INTENT ( IN ) :: work !< Mechanical work load (W) 2821 3110 ! 2822 3111 !-- In and output parameters 2823 3112 REAL(wp), INTENT (INOUT) :: storage !< Heat storage (W) 2824 3113 REAL(wp), INTENT (INOUT) :: t_clothing !< Clothig temperature (°C) 2825 3114 ! 2826 3115 !-- Type of output of the argument list 2827 3116 REAL(wp), INTENT ( OUT ) :: ipt !< Instationary perceived temperature (°C) 2828 3117 ! 2829 3118 !-- Type of internal variables 2830 3119 REAL(wp) :: d_pmv … … 2841 3130 2842 3131 LOGICAL :: sultrieness 2843 3132 ! 2844 3133 !-- Initialise 2845 ipt = -999.0_wp3134 ipt = bio_fill_value 2846 3135 2847 3136 nerr = 0_iwp 2848 3137 sultrieness = .FALSE. 2849 3138 ! 2850 3139 !-- Determine pmv_adjusted for current conditions 2851 3140 CALL fanger_s_acti ( ta, tmrt, vp, ws, pair, clo, actlev, work, & 2852 3141 t_clothing, storage, dt, pmva ) 2853 3142 ! 2854 3143 !-- Determine perceived temperature by regression equation + adjustments 2855 3144 CALL perct_regression ( pmva, clo, ipt ) 2856 3145 2857 3146 IF ( clo >= 1.75_wp .AND. pmva <= -0.11_wp ) THEN 3147 ! 2858 3148 !-- Adjust for cold conditions according to Gagge 1986 2859 3149 CALL dpmv_cold ( pmva, ta, ws, tmrt, nerr_cold, d_pmv ) … … 2866 3156 CALL perct_regression ( pmvs, clo, ipt ) 2867 3157 ENDIF 2868 3158 ! 2869 3159 !-- Consider sultriness if appropriate 2870 3160 ptc = ipt … … 2873 3163 d_std = -99._wp 2874 3164 IF ( pmva > 0.06_wp .AND. clo <= 0.5_wp ) THEN 3165 ! 2875 3166 !-- Adjust for warm/humid conditions according to Gagge 1986 2876 3167 CALL saturation_vapor_pressure ( ta, svp_ta ) … … 2898 3189 2899 3190 IMPLICIT NONE 2900 3191 ! 2901 3192 !-- Input argument types 2902 3193 REAL(wp), INTENT ( IN ) :: ta !< Air temperature (°C) … … 2909 3200 REAL(wp), INTENT ( IN ) :: activity !< Work load (W/m²) 2910 3201 REAL(wp), INTENT ( IN ) :: in_clo !< Clothing index (clo) (no dim) 2911 3202 ! 2912 3203 !-- Output argument types 2913 3204 REAL(wp), INTENT ( OUT ) :: pmva !< actual Predicted Mean Vote (no dim) … … 2915 3206 REAL(wp), INTENT (INOUT) :: s !< storage var. of energy balance (W/m2) 2916 3207 REAL(wp), INTENT (INOUT) :: t_cloth !< clothing temperature (°C) 2917 3208 ! 2918 3209 !-- Internal variables 2919 3210 REAL(wp), PARAMETER :: time_equil = 7200._wp … … 2928 3219 REAL(wp) :: gc !< preliminary result storage 2929 3220 REAL(wp) :: t_clothing !< clothing temperature (°C) 2930 REAL(wp) :: hr !< radiational heat resistence3221 ! REAL(wp) :: hr !< radiational heat resistence 2931 3222 REAL(wp) :: clo !< clothing insulation index (clo) 2932 3223 REAL(wp) :: ws !< wind speed (m/s) … … 2946 3237 INTEGER(iwp) :: niter !< Running index 2947 3238 2948 3239 ! 2949 3240 !-- Clo must be > 0. to avoid div. by 0! 2950 3241 clo = in_clo 2951 3242 IF ( clo < 001._wp ) clo = .001_wp 2952 3243 ! 2953 3244 !-- Increase in surface due to clothing 2954 3245 f_cl = 1._wp + .15_wp * clo 2955 3246 ! 2956 3247 !-- Case of free convection (ws < 0.1 m/s ) not considered 2957 3248 ws = in_ws … … 2959 3250 ws = .1_wp 2960 3251 ENDIF 2961 3252 ! 2962 3253 !-- Heat_convection = forced convection 2963 3254 heat_convection = 12.1_wp * SQRT ( ws * pair / 1013.25_wp ) 2964 3255 ! 2965 3256 !-- Average skin temperature 2966 3257 t_skin_aver = 35.7_wp - .0275_wp * activity 2967 3258 ! 2968 3259 !-- Calculation of constants for evaluation below 2969 3260 bc = .155_wp * clo * 3.96_wp * 10._wp**( -8._wp ) * f_cl … … 2972 3263 dc = ( 1._wp + ec * cc ) / bc 2973 3264 gc = ( t_skin_aver + bc * ( tmrt + 273.2_wp )**4._wp + ec * cc * ta ) / bc 2974 3265 ! 2975 3266 !-- Calculation of clothing surface temperature (t_clothing) based on 2976 3267 ! newton-approximation with air temperature as initial guess 2977 niter = dt 3268 niter = INT( dt * 10._wp, KIND=iwp ) 3269 IF ( niter < 1 ) niter = 1_iwp 2978 3270 adjustrate = 1._wp - EXP ( -1._wp * ( 10._wp / time_equil ) * dt ) 2979 3271 IF ( adjustrate >= 1._wp ) adjustrate = 1._wp … … 2982 3274 2983 3275 IF ( t_cloth <= -998.0_wp ) THEN ! If initial run 2984 niter = 3_ wp3276 niter = 3_iwp 2985 3277 adjustrate = 1._wp 2986 3278 adjustrate_cloth = 1._wp … … 2993 3285 dc - gc ) / ( 4._wp * ( t_clothing + 273.2_wp )**3._wp + dc ) 2994 3286 ENDDO 2995 3287 ! 2996 3288 !-- Empiric factor for the adaption of the heat ballance equation 2997 3289 ! to the psycho-physical scale (Equ. 40 in FANGER) 2998 3290 z1 = ( .303_wp * EXP ( -.036_wp * actlev ) + .0275_wp ) 2999 3291 ! 3000 3292 !-- Water vapour diffution through the skin 3001 3293 z2 = .31_wp * ( 57.3_wp - .07_wp * activity-pa ) 3002 3294 ! 3003 3295 !-- Sweat evaporation from the skin surface 3004 3296 z3 = .42_wp * ( activity - 58._wp ) 3005 3297 ! 3006 3298 !-- Loss of latent heat through respiration 3007 3299 z4 = .0017_wp * actlev * ( 58.7_wp - pa ) + .0014_wp * actlev * & 3008 3300 ( 34._wp - ta ) 3009 3301 ! 3010 3302 !-- Loss of radiational heat 3011 3303 z5 = 3.96e-8_wp * f_cl * ( ( t_clothing + 273.2_wp )**4 - ( tmrt + & 3012 3304 273.2_wp )**4 ) 3013 3305 ! 3014 3306 !-- Heat loss through forced convection 3015 3307 z6 = cc * ( t_clothing - ta ) 3016 3308 ! 3017 3309 !-- Write together as energy ballance 3018 3310 en = activity - z2 - z3 - z4 - z5 - z6 3019 3311 ! 3020 3312 !-- Manage storage 3021 3313 d_s = adjustrate * en + ( 1._wp - adjustrate ) * s 3022 3314 ! 3023 3315 !-- Predicted Mean Vote 3024 3316 pmva = z1 * d_s 3025 3317 ! 3026 3318 !-- Update storage 3027 3319 s = d_s … … 3041 3333 !------------------------------------------------------------------------------! 3042 3334 3043 SUBROUTINE calculate_pet_static( ta, vpa, v, tmrt, pair, pet )3335 SUBROUTINE calculate_pet_static( ta, vpa, v, tmrt, pair, pet_ij ) 3044 3336 3045 3337 IMPLICIT NONE 3046 3338 ! 3047 3339 !-- Input arguments: 3048 3340 REAL(wp), INTENT( IN ) :: ta !< Air temperature (°C) … … 3051 3343 REAL(wp), INTENT( IN ) :: vpa !< Vapor pressure (hPa) 3052 3344 REAL(wp), INTENT( IN ) :: pair !< Air pressure (hPa) 3053 3345 ! 3054 3346 !-- Output arguments: 3055 REAL(wp), INTENT ( OUT ) :: pet !< PET(°C)3056 3347 REAL(wp), INTENT ( OUT ) :: pet_ij !< PET (°C) 3348 ! 3057 3349 !-- Internal variables: 3058 3350 REAL(wp) :: acl !< clothing area (m²) … … 3071 3363 REAL(wp) :: tcl !< Clothing temperature (°C) 3072 3364 REAL(wp) :: wetsk !< Fraction of wet skin (factor) 3073 3365 ! 3074 3366 !-- Variables: 3075 3367 REAL(wp) :: int_heat !< Internal heat (W) 3076 3368 ! 3077 3369 !-- MEMI configuration 3078 3370 REAL(wp) :: age !< Persons age (a) … … 3085 3377 ! INTEGER(iwp) :: pos !< Posture: 1 = standing, 2 = sitting 3086 3378 ! INTEGER(iwp) :: sex !< Sex: 1 = male, 2 = female 3087 3379 ! 3088 3380 !-- Configuration, keep standard parameters! 3089 3381 age = 35._wp … … 3094 3386 clo = 0.9_wp 3095 3387 fcl = 1.15_wp 3096 3388 ! 3097 3389 !-- Call subfunctions 3098 3390 CALL in_body ( age, eta, ere, erel, ht, int_heat, mbody, pair, rtv, ta, & … … 3103 3395 vpts, wetsk ) 3104 3396 3105 CALL pet_iteration ( acl, adu, aeff, esw, facl, feff, int_heat, pair, rdcl,&3106 rd sk, rtv, ta, tcl, tsk, pet, vpts, wetsk )3397 CALL pet_iteration ( acl, adu, aeff, esw, facl, feff, int_heat, pair, & 3398 rdcl, rdsk, rtv, ta, tcl, tsk, pet_ij, vpts, wetsk ) 3107 3399 3108 3400 … … 3117 3409 SUBROUTINE in_body ( age, eta, ere, erel, ht, int_heat, mbody, pair, rtv, ta, & 3118 3410 vpa, work ) 3119 3411 ! 3120 3412 !-- Input arguments: 3121 3413 REAL(wp), INTENT( IN ) :: pair !< air pressure (hPa) … … 3127 3419 REAL(wp), INTENT( IN ) :: work !< Work load (W) 3128 3420 REAL(wp), INTENT( IN ) :: eta !< Work efficiency (dimensionless) 3129 3421 ! 3130 3422 !-- Output arguments: 3131 3423 REAL(wp), INTENT( OUT ) :: ere !< energy ballance (W) … … 3133 3425 REAL(wp), INTENT( OUT ) :: int_heat !< internal heat production (W) 3134 3426 REAL(wp), INTENT( OUT ) :: rtv !< respiratory volume 3135 3427 ! 3136 3428 !-- Constants: 3137 3429 ! REAL(wp), PARAMETER :: cair = 1010._wp !< replaced by c_p 3138 3430 ! REAL(wp), PARAMETER :: evap = 2.42_wp * 10._wp **6._wp !< replaced by l_v 3139 3431 ! 3140 3432 !-- Internal variables: 3141 3433 REAL(wp) :: eres !< Sensible respiratory heat flux (W) … … 3151 3443 3152 3444 int_heat = met * (1._wp - eta) 3153 3445 ! 3154 3446 !-- Sensible respiration energy 3155 3447 tex = 0.47_wp * ta + 21.0_wp 3156 3448 rtv = 1.44_wp * 10._wp ** (-6._wp) * met 3157 3449 eres = c_p * (ta - tex) * rtv 3158 3450 ! 3159 3451 !-- Latent respiration energy 3160 3452 vpex = 6.11_wp * 10._wp ** ( 7.45_wp * tex / ( 235._wp + tex ) ) 3161 3453 erel = 0.623_wp * l_v / pair * ( vpa - vpex ) * rtv 3162 3454 ! 3163 3455 !-- Sum of the results 3164 3456 ere = eres + erel … … 3176 3468 vpts, wetsk ) 3177 3469 3178 3470 ! 3179 3471 !-- Input arguments: 3180 3472 REAL(wp), INTENT( IN ) :: ere !< Energy ballance (W) … … 3190 3482 REAL(wp), INTENT( IN ) :: clo !< clothing insulation (clo) 3191 3483 REAL(wp), INTENT( IN ) :: fcl !< factor for surface area increase by clothing 3192 3484 ! 3193 3485 !-- Output arguments: 3194 3486 REAL(wp), INTENT( OUT ) :: acl !< Clothing surface area (m²) … … 3204 3496 REAL(wp), INTENT( OUT ) :: vpts !< Sat. vapor pressure over skin (hPa) 3205 3497 REAL(wp), INTENT( OUT ) :: wetsk !< Fraction of wet skin (dimensionless) 3206 3498 ! 3207 3499 !-- Cconstants: 3208 3500 ! REAL(wp), PARAMETER :: cair = 1010._wp !< replaced by c_p … … 3214 3506 REAL(wp), PARAMETER :: po = 1013.25_wp !< Air pressure at sea level (hPa) 3215 3507 REAL(wp), PARAMETER :: rob = 1.06_wp !< 3216 3508 ! 3217 3509 !-- Internal variables 3218 3510 REAL(wp) :: c(0:10) !< Core temperature array (°C) … … 3239 3531 REAL(wp) :: rsum !< Radiational loss or gain (W/m²) 3240 3532 REAL(wp) :: sw !< 3241 REAL(wp) :: swf !< 3533 ! REAL(wp) :: swf !< female factor, currently unused 3242 3534 REAL(wp) :: swm !< 3243 3535 REAL(wp) :: tbody !< … … 3279 3571 3280 3572 di = r2 - r1 3281 3573 ! 3282 3574 !-- Skin temperatur 3283 3575 DO j = 1, 7 … … 3295 3587 htcl = 6.28_wp * ht * y * di / ( rcl * LOG( r2 / r1 ) * acl ) 3296 3588 tsk = 1._wp / htcl * ( hc * ( tcl - ta ) + rclo2 ) + tcl 3297 3589 ! 3298 3590 !-- Radiation saldo 3299 3591 aeff = adu * feff … … 3303 3595 ( ( tmrt + degc_to_k )** 4._wp - ( tcl + degc_to_k )** 4._wp ) 3304 3596 rsum = rbare + rclo 3305 3597 ! 3306 3598 !-- Convection 3307 3599 cbare = hc * ( ta - tsk ) * adu * ( 1._wp - facl ) 3308 3600 cclo = hc * ( ta - tcl ) * acl 3309 3601 csum = cbare + cclo 3310 3602 ! 3311 3603 !-- Core temperature 3312 3604 c(0) = int_heat + ere … … 3337 3629 tcore(4) = c(0) / ( 5.28_wp * adu + c(1) * 1._wp / 40._wp ) + tsk 3338 3630 END IF 3339 3631 ! 3340 3632 !-- Transpiration 3341 3633 tbody = 0.1_wp * tsk + 0.9_wp * tcore(j) … … 3360 3652 IF ( eswdif > 0._wp ) esw = eswphy !< Limit is sweat production 3361 3653 IF ( esw > 0._wp ) esw = 0._wp !< Sweat can't be evaporated, no more cooling effect 3362 3654 ! 3363 3655 !-- Diffusion 3364 3656 rdsk = 0.79_wp * 10._wp ** 7._wp … … 3366 3658 ed = l_v / ( rdsk + rdcl ) * adu * ( 1._wp - wetsk ) * ( vpa - & 3367 3659 vpts ) 3368 3660 ! 3369 3661 !-- Max vb 3370 3662 vb1 = 34._wp - tsk … … 3374 3666 IF ( vb1 < 0._wp ) vb1 = 0._wp 3375 3667 vb = ( 6.3_wp + 75._wp * vb2 ) / ( 1._wp + 0.5_wp * vb1 ) 3376 3668 ! 3377 3669 !-- Energy ballence 3378 3670 enbal = int_heat + ed + ere + esw + csum + rsum + food 3379 3671 ! 3380 3672 !-- Clothing temperature 3381 3673 xx = 0.001_wp … … 3429 3721 IF ( ( j == 4_iwp ) .AND. ( vb < 89._wp ) ) CYCLE 3430 3722 IF ( vb > 90._wp ) vb = 90._wp 3431 3723 ! 3432 3724 !-- Loses by water 3433 3725 ws = sw * 3600._wp * 1000._wp … … 3448 3740 !------------------------------------------------------------------------------! 3449 3741 SUBROUTINE pet_iteration ( acl, adu, aeff, esw, facl, feff, int_heat, pair, & 3450 rdcl, rdsk, rtv, ta, tcl, tsk, pet , vpts, wetsk )3451 3742 rdcl, rdsk, rtv, ta, tcl, tsk, pet_ij, vpts, wetsk ) 3743 ! 3452 3744 !-- Input arguments: 3453 3745 REAL(wp), INTENT( IN ) :: acl !< clothing surface area (m²) … … 3466 3758 REAL(wp), INTENT( IN ) :: vpts !< sat. vapor pressure over skin (hPa) 3467 3759 REAL(wp), INTENT( IN ) :: wetsk !< fraction of wet skin (dimensionless) 3468 3760 ! 3469 3761 !-- Output arguments: 3470 REAL(wp), INTENT( OUT ) :: aeff !< effective surface area (m²)3471 REAL(wp), INTENT( OUT ) :: pet !< PET (°C)3472 3762 REAL(wp), INTENT( OUT ) :: aeff !< effective surface area (m²) 3763 REAL(wp), INTENT( OUT ) :: pet_ij !< PET (°C) 3764 ! 3473 3765 !-- Cconstants: 3474 ! REAL(wp), PARAMETER :: cair = 1010._wp !< replaced by c_p3475 3766 REAL(wp), PARAMETER :: emcl = 0.95_wp !< Longwave emission coef. of cloth 3476 3767 REAL(wp), PARAMETER :: emsk = 0.99_wp !< Longwave emission coef. of skin 3477 ! REAL(wp), PARAMETER :: evap = 2.42_wp * 10._wp **6._wp !< replaced by l_v3478 3768 REAL(wp), PARAMETER :: po = 1013.25_wp !< Air pressure at sea level (hPa) 3479 3769 ! 3480 3770 !-- Internal variables 3481 3771 REAL ( wp ) :: cbare !< Convection through bare skin … … 3499 3789 INTEGER ( iwp ) :: i !< running index 3500 3790 3501 pet = ta3791 pet_ij = ta 3502 3792 enbal2 = 0._wp 3503 3793 … … 3506 3796 hc = 2.67_wp + 6.5_wp * 0.1_wp ** 0.67_wp 3507 3797 hc = hc * ( pair / po ) ** 0.55_wp 3508 3798 ! 3509 3799 !-- Radiation 3510 3800 aeff = adu * feff 3511 3801 rbare = aeff * ( 1._wp - facl ) * emsk * sigma_sb * & 3512 ( ( pet + degc_to_k ) ** 4._wp - ( tsk + degc_to_k ) ** 4._wp )3802 ( ( pet_ij + degc_to_k ) ** 4._wp - ( tsk + degc_to_k ) ** 4._wp ) 3513 3803 rclo = feff * acl * emcl * sigma_sb * & 3514 ( ( pet + degc_to_k ) ** 4._wp - ( tcl + degc_to_k ) ** 4._wp )3804 ( ( pet_ij + degc_to_k ) ** 4._wp - ( tcl + degc_to_k ) ** 4._wp ) 3515 3805 rsum = rbare + rclo 3516 3806 ! 3517 3807 !-- Covection 3518 cbare = hc * ( pet - tsk ) * adu * ( 1._wp - facl )3519 cclo = hc * ( pet - tcl ) * acl3808 cbare = hc * ( pet_ij - tsk ) * adu * ( 1._wp - facl ) 3809 cclo = hc * ( pet_ij - tcl ) * acl 3520 3810 csum = cbare + cclo 3521 3811 ! 3522 3812 !-- Diffusion 3523 3813 ed = l_v / ( rdsk + rdcl ) * adu * ( 1._wp - wetsk ) * ( 12._wp - & 3524 3814 vpts ) 3525 3815 ! 3526 3816 !-- Respiration 3527 tex = 0.47_wp * pet + 21._wp3528 eres = c_p * ( pet - tex ) * rtv3817 tex = 0.47_wp * pet_ij + 21._wp 3818 eres = c_p * ( pet_ij - tex ) * rtv 3529 3819 vpex = 6.11_wp * 10._wp ** ( 7.45_wp * tex / ( 235._wp + tex ) ) 3530 3820 erel = 0.623_wp * l_v / pair * ( 12._wp - vpex ) * rtv 3531 3821 ere = eres + erel 3532 3822 ! 3533 3823 !-- Energy ballance 3534 3824 enbal = int_heat + ed + ere + esw + csum + rsum 3535 3825 ! 3536 3826 !-- Iteration concerning ta 3827 xx = 0.001_wp 3537 3828 IF ( count1 == 0_iwp ) xx = 1._wp 3538 3829 IF ( count1 == 1_iwp ) xx = 0.1_wp 3539 3830 IF ( count1 == 2_iwp ) xx = 0.01_wp 3540 IF ( count1 == 3_iwp ) xx = 0.001_wp3541 IF ( enbal > 0._wp ) pet = pet- xx3542 IF ( enbal < 0._wp ) pet = pet+ xx3831 ! IF ( count1 == 3_iwp ) xx = 0.001_wp 3832 IF ( enbal > 0._wp ) pet_ij = pet_ij - xx 3833 IF ( enbal < 0._wp ) pet_ij = pet_ij + xx 3543 3834 IF ( ( enbal <= 0._wp ) .AND. ( enbal2 > 0._wp ) ) EXIT 3544 3835 IF ( ( enbal >= 0._wp ) .AND. ( enbal2 < 0._wp ) ) EXIT … … 3549 3840 END SUBROUTINE pet_iteration 3550 3841 3842 ! 3843 !-- UVEM specific subroutines 3844 3845 !---------------------------------------------------------------------------------------------------------------------! 3846 ! Description: 3847 ! ------------ 3848 !> Module-specific routine for new module 3849 !---------------------------------------------------------------------------------------------------------------------! 3850 SUBROUTINE uvem_solar_position 3851 3852 USE date_and_time_mod, & 3853 ONLY: calc_date_and_time, day_of_year, time_utc 3854 3855 USE control_parameters, & 3856 ONLY: latitude, longitude 3857 3858 IMPLICIT NONE 3859 3860 3861 REAL(wp) :: alpha = 0.0_wp !< solar azimuth angle in radiant 3862 REAL(wp) :: doy_r = 0.0_wp !< real format of day_of_year 3863 REAL(wp) :: declination = 0.0_wp !< declination 3864 REAL(wp) :: dtor = 0.0_wp !< factor to convert degree to radiant 3865 REAL(wp) :: js = 0.0_wp !< parameter for solar position calculation 3866 REAL(wp) :: lat = 52.39_wp !< latitude 3867 REAL(wp) :: lon = 9.7_wp !< longitude 3868 REAL(wp) :: thetar = 0.0_wp !< angle for solar zenith angle calculation 3869 REAL(wp) :: thetasr = 0.0_wp !< angle for solar azimuth angle calculation 3870 REAL(wp) :: zgl = 0.0_wp !< calculated exposure by direct beam 3871 REAL(wp) :: woz = 0.0_wp !< calculated exposure by diffuse radiation 3872 REAL(wp) :: wsp = 0.0_wp !< calculated exposure by direct beam 3873 3874 3875 CALL calc_date_and_time 3876 doy_r = real(day_of_year) 3877 dtor = pi / 180.0_wp 3878 lat = latitude 3879 lon = longitude 3880 ! 3881 !-- calculation of js, necessary for calculation of equation of time (zgl) : 3882 js= 72.0_wp * ( doy_r + ( time_utc / 86400.0_wp ) ) / 73.0_wp 3883 ! 3884 !-- calculation of equation of time (zgl): 3885 zgl = 0.0066_wp + 7.3525_wp * cos( ( js + 85.9_wp ) * dtor ) + 9.9359_wp * & 3886 cos( ( 2.0_wp * js + 108.9_wp ) * dtor ) + 0.3387_wp * cos( ( 3 * js + 105.2_wp ) * dtor ) 3887 ! 3888 !-- calculation of apparent solar time woz: 3889 woz = ( ( time_utc / 3600.0_wp ) - ( 4.0_wp * ( 15.0_wp - lon ) ) / 60.0_wp ) + ( zgl / 60.0_wp ) 3890 ! 3891 !-- calculation of hour angle (wsp): 3892 wsp = ( woz - 12.0_wp ) * 15.0_wp 3893 ! 3894 !-- calculation of declination: 3895 declination = 0.3948_wp - 23.2559_wp * cos( ( js + 9.1_wp ) * dtor ) - & 3896 0.3915_wp * cos( ( 2.0_wp * js + 5.4_wp ) * dtor ) - 0.1764_wp * cos( ( 3.0_wp * js + 26.0_wp ) * dtor ) 3897 ! 3898 !-- calculation of solar zenith angle 3899 thetar = acos( sin( lat * dtor) * sin( declination * dtor ) + cos( wsp * dtor ) * & 3900 cos( lat * dtor ) * cos( declination * dtor ) ) 3901 thetasr = asin( sin( lat * dtor) * sin( declination * dtor ) + cos( wsp * dtor ) * & 3902 cos( lat * dtor ) * cos( declination * dtor ) ) 3903 sza = thetar / dtor 3904 ! 3905 !-- calculation of solar azimuth angle 3906 IF (woz .LE. 12.0_wp) alpha = pi - acos( ( sin(thetasr) * sin( lat * dtor ) - & 3907 sin( declination * dtor ) ) / ( cos(thetasr) * cos( lat * dtor ) ) ) 3908 IF (woz .GT. 12.0_wp) alpha = pi + acos( ( sin(thetasr) * sin( lat * dtor ) - & 3909 sin( declination * dtor ) ) / ( cos(thetasr) * cos( lat * dtor ) ) ) 3910 saa = alpha / dtor 3911 3912 END SUBROUTINE uvem_solar_position 3913 3914 3915 !------------------------------------------------------------------------------! 3916 ! Description: 3917 ! ------------ 3918 !> Module-specific routine for new module 3919 !---------------------------------------------------------------------------------------------------------------------! 3920 SUBROUTINE uvem_calc_exposure 3921 3922 USE indices, & 3923 ONLY: nxlg, nxrg, nyng, nysg, nys, nyn, nxl, nxr 3924 3925 3926 IMPLICIT NONE 3927 3928 INTEGER(iwp) :: i !< loop index in x direction 3929 INTEGER(iwp) :: j !< loop index in y direction 3930 INTEGER(iwp) :: szai !< loop index for different sza values 3931 3932 CALL uvem_solar_position 3933 3934 IF (sza .GE. 90) THEN 3935 vitd3_exposure(:,:) = 0.0_wp 3936 ELSE 3937 3938 DO ai = 0, 35 3939 DO zi = 0, 9 3940 projection_area_lookup_table(ai,zi) = uvem_projarea_f%var(clothing,zi,ai) 3941 ENDDO 3942 ENDDO 3943 DO ai = 0, 35 3944 DO zi = 0, 9 3945 integration_array(ai,zi) = uvem_integration_f%var(zi,ai) 3946 ENDDO 3947 ENDDO 3948 DO ai = 0, 2 3949 DO zi = 0, 90 3950 irradiance_lookup_table(ai,zi) = uvem_irradiance_f%var(zi,ai) 3951 ENDDO 3952 ENDDO 3953 DO ai = 0, 35 3954 DO zi = 0, 9 3955 DO szai = 0, 90 3956 radiance_lookup_table(ai,zi,szai) = uvem_radiance_f%var(szai,zi,ai) 3957 ENDDO 3958 ENDDO 3959 ENDDO 3960 3961 3962 3963 !-- rotate 3D-Model human to desired direction ----------------------------- 3964 projection_area_temp( 0:35,:) = projection_area_lookup_table 3965 projection_area_temp(36:71,:) = projection_area_lookup_table 3966 IF ( .NOT. turn_to_sun ) startpos_human = orientation_angle / 10.0_wp 3967 IF ( turn_to_sun ) startpos_human = saa / 10.0_wp 3968 DO ai = 0, 35 3969 xfactor = ( startpos_human ) - INT( startpos_human ) 3970 DO zi = 0, 9 3971 projection_area(ai,zi) = ( projection_area_temp( 36 - INT( startpos_human ) - 1 + ai , zi) * & 3972 ( xfactor ) ) & 3973 +( projection_area_temp( 36 - INT( startpos_human ) + ai , zi) * & 3974 ( 1.0_wp - xfactor ) ) 3975 ENDDO 3976 ENDDO 3977 ! 3978 ! 3979 !-- interpolate to accurate Solar Zenith Angle ------------------ 3980 DO ai = 0, 35 3981 xfactor = (sza)-INT(sza) 3982 DO zi = 0, 9 3983 radiance_array(ai,zi) = ( radiance_lookup_table(ai, zi, INT(sza) ) * ( 1.0_wp - xfactor) ) + & 3984 ( radiance_lookup_table(ai,zi,INT(sza) + 1) * xfactor ) 3985 ENDDO 3986 ENDDO 3987 DO iq = 0, 2 3988 irradiance(iq) = ( irradiance_lookup_table(iq, INT(sza) ) * ( 1.0_wp - xfactor)) + & 3989 (irradiance_lookup_table(iq, INT(sza) + 1) * xfactor ) 3990 ENDDO 3991 ! 3992 !-- interpolate to accurate Solar Azimuth Angle ------------------ 3993 IF ( sun_in_south ) THEN 3994 startpos_saa_float = 180.0_wp / 10.0_wp 3995 ELSE 3996 startpos_saa_float = saa / 10.0_wp 3997 ENDIF 3998 radiance_array_temp( 0:35,:) = radiance_array 3999 radiance_array_temp(36:71,:) = radiance_array 4000 xfactor = (startpos_saa_float) - INT(startpos_saa_float) 4001 DO ai = 0, 35 4002 DO zi = 0, 9 4003 radiance_array(ai,zi) = ( radiance_array_temp( 36 - INT( startpos_saa_float ) - 1 + ai , zi ) * & 4004 ( xfactor ) ) & 4005 + ( radiance_array_temp( 36 - INT( startpos_saa_float ) + ai , zi ) & 4006 * ( 1.0_wp - xfactor ) ) 4007 ENDDO 4008 ENDDO 4009 ! 4010 ! 4011 !-- calculate Projectionarea for direct beam -----------------------------' 4012 projection_area_direct_temp( 0:35,:) = projection_area 4013 projection_area_direct_temp(36:71,:) = projection_area 4014 yfactor = ( sza / 10.0_wp ) - INT( sza / 10.0_wp ) 4015 xfactor = ( startpos_saa_float ) - INT( startpos_saa_float ) 4016 projection_area_direct_beam = ( projection_area_direct_temp( INT(startpos_saa_float) ,INT(sza/10.0_wp) ) * & 4017 ( 1.0_wp - xfactor ) * ( 1.0_wp - yfactor ) ) + & 4018 ( projection_area_direct_temp( INT(startpos_saa_float) + 1,INT(sza/10.0_wp) ) * & 4019 ( xfactor ) * ( 1.0_wp - yfactor ) ) + & 4020 ( projection_area_direct_temp( INT(startpos_saa_float) ,INT(sza/10.0_wp)+1) * & 4021 ( 1.0_wp - xfactor ) * ( yfactor ) ) + & 4022 ( projection_area_direct_temp( INT(startpos_saa_float) + 1,INT(sza/10.0_wp)+1) * & 4023 ( xfactor ) * ( yfactor ) ) 4024 ! 4025 ! 4026 ! 4027 DO i = nxl, nxr !nxlg, nxrg 4028 DO j = nys, nyn !nysg, nyng 4029 ! 4030 ! !-- extract obstruction from IBSET-Integer_Array ------------------' 4031 IF (consider_obstructions ) THEN 4032 obstruction_temp1 = building_obstruction_f%var_3d(:,j,i) 4033 IF (obstruction_temp1(0) .NE. 9) THEN 4034 DO pobi = 0, 44 4035 DO bi = 0, 7 4036 IF ( btest( obstruction_temp1(pobi), bi ) .EQV. .TRUE.) THEN 4037 obstruction_temp2( ( pobi * 8 ) + bi ) = 1 4038 ELSE 4039 obstruction_temp2( ( pobi * 8 ) + bi ) = 0 4040 ENDIF 4041 ENDDO 4042 ENDDO 4043 DO zi = 0, 9 4044 obstruction(:,zi) = obstruction_temp2( zi * 36 :( zi * 36) + 35 ) 4045 ENDDO 4046 ELSE 4047 obstruction(:,:) = 0 4048 ENDIF 4049 ENDIF 4050 ! 4051 ! !-- calculated human exposure ------------------' 4052 diffuse_exposure = SUM( radiance_array * projection_area * integration_array * obstruction ) 4053 4054 obstruction_direct_beam = obstruction( nint(startpos_saa_float), nint( sza / 10.0_wp ) ) 4055 IF (sza .GE. 89.99_wp) THEN 4056 sza = 89.99999_wp 4057 ENDIF 4058 ! 4059 !-- calculate direct normal irradiance (direct beam) ------------------' 4060 direct_exposure = ( irradiance(1) / cos( pi * sza / 180.0_wp ) ) * & 4061 projection_area_direct_beam * obstruction_direct_beam 4062 4063 vitd3_exposure(j,i) = ( diffuse_exposure + direct_exposure ) / 1000.0_wp * 70.97_wp 4064 ! unit = international units vitamin D per second 4065 ENDDO 4066 ENDDO 4067 ENDIF 4068 4069 END SUBROUTINE uvem_calc_exposure 3551 4070 3552 4071 END MODULE biometeorology_mod -
palm/trunk/SOURCE/check_parameters.f90
r3545 r3569 25 25 ! ----------------- 26 26 ! $Id$ 27 ! Sort chemistry routine calls into list. 28 ! dom_dwd_user, Schrempf: 29 ! remove CALLs to uv exposure model, this is now part of biometeorology_mod 30 ! 31 ! 3545 2018-11-21 11:19:41Z gronemeier 27 32 ! Call tcm_check_parameters before other modules 28 33 ! … … 823 828 USE urban_surface_mod, & 824 829 ONLY: usm_check_data_output, usm_check_parameters 825 826 USE uv_exposure_model_mod, &827 ONLY: uvem_check_data_output828 830 829 831 USE wind_turbine_model_mod, & … … 1507 1509 ENDIF 1508 1510 1509 ! Check for chem_emission_mod parameters setting 1510 ! IF ( air_chemistry ) CALL chem_emissions_check_parameters ! forkel preliminary 1511 1512 1513 ! Check for chemitry_model_mod parameters setting 1514 IF ( air_chemistry ) CALL chem_check_parameters 1515 1516 1517 !-- Check the module settings 1511 ! 1512 !-- Check the module settings: 1513 !-- tcm_check_parameters must be called at first 1518 1514 CALL tcm_check_parameters 1515 1519 1516 IF ( biometeorology ) CALL bio_check_parameters 1520 1517 IF ( bulk_cloud_model ) CALL bcm_check_parameters 1518 IF ( air_chemistry ) CALL chem_check_parameters 1519 ! IF ( air_chemistry ) CALL chem_emissions_check_parameters ! forkel preliminary 1521 1520 IF ( gust_module_enabled ) CALL gust_check_parameters 1522 1521 IF ( large_scale_forcing .OR. nudging ) & … … 1527 1526 IF ( radiation ) CALL radiation_check_parameters 1528 1527 IF ( calculate_spectra ) CALL spectra_check_parameters 1529 CALL stg_check_parameters1528 CALL stg_check_parameters 1530 1529 IF ( urban_surface ) CALL usm_check_parameters 1531 1530 IF ( wind_turbine ) CALL wtm_check_parameters … … 3220 3219 3221 3220 IF ( unit == 'illegal' .AND. biometeorology ) THEN 3222 CALL bio_check_data_output( var, unit )3221 CALL bio_check_data_output( var, unit, i, ilen, k ) 3223 3222 ENDIF 3224 3223 … … 3244 3243 CALL usm_check_data_output( var, unit ) 3245 3244 ENDIF 3246 3247 IF ( unit == 'illegal' .AND. uv_exposure &3248 .AND. var(1:5) == 'uvem_' ) THEN3249 CALL uvem_check_data_output( var, unit, i, ilen, k )3250 ENDIF3251 3252 3245 ! 3253 3246 !-- Finally, check for user-defined quantities -
palm/trunk/SOURCE/data_output_2d.f90
r3554 r3569 19 19 ! 20 20 ! Current revisions: 21 ! ----------------- 21 ! ------------------ 22 22 ! 23 23 ! … … 25 25 ! ----------------- 26 26 ! $Id$ 27 ! Remove fill_value from bio_data_output_2d call 28 ! dom_dwd_user, Schrempf: 29 ! Clean up of biometeorology calls, 30 ! remove uv exposure model code, this is now part of biometeorology_mod. 31 ! 32 ! 3554 2018-11-22 11:24:52Z gronemeier 27 33 ! - add variable description 28 34 ! - rename variable 'if' into 'ivar' … … 289 295 ntdim_2d_xy, ntdim_2d_xz, ntdim_2d_yz, & 290 296 ocean_mode, psolver, section, simulated_time, & 291 time_since_reference_point , uv_exposure297 time_since_reference_point 292 298 293 299 USE cpulog, & … … 336 342 USE turbulence_closure_mod, & 337 343 ONLY: tcm_data_output_2d 338 339 USE uv_exposure_model_mod, &340 ONLY: uvem_data_output_2d341 344 342 345 … … 1328 1331 ENDIF 1329 1332 1330 IF ( .NOT. found .AND. biometeorology & 1331 .AND. mode == 'xy' ) THEN 1333 IF ( .NOT. found .AND. biometeorology ) THEN 1332 1334 CALL bio_data_output_2d( av, do2d(av,ivar), found, grid, & 1333 local_pf, two_d, nzb_do, nzt_do, & 1334 fill_value ) 1335 local_pf, two_d, nzb_do, nzt_do ) 1335 1336 ENDIF 1336 1337 … … 1355 1356 mode, local_pf, two_d ) 1356 1357 ENDIF 1357 1358 IF ( .NOT. found .AND. uv_exposure ) THEN1359 CALL uvem_data_output_2d( av, do2d(av,ivar), found, grid, &1360 local_pf, two_d, nzb_do, nzt_do )1361 ENDIF1362 1358 1363 1359 IF ( .NOT. found .AND. air_chemistry ) THEN -
palm/trunk/SOURCE/init_3d_model.f90
r3547 r3569 25 25 ! ----------------- 26 26 ! $Id$ 27 ! dom_dwd_user, Schrempf: 28 ! Remove uv exposure model code, this is now part of biometeorology_mod, 29 ! remove bio_init_arrays. 30 ! 31 ! 3547 2018-11-21 13:21:24Z suehring 27 32 ! variables documented 28 33 ! … … 560 565 561 566 USE biometeorology_mod, & 562 ONLY: bio_init , bio_init_arrays567 ONLY: bio_init 563 568 564 569 USE bulk_cloud_model_mod, & … … 677 682 ONLY: usm_init_urban_surface, usm_allocate_surface 678 683 679 USE uv_exposure_model_mod, &680 ONLY: uvem_init, uvem_init_arrays681 682 684 USE virtual_measurement_mod, & 683 685 ONLY: vm_init … … 708 710 INTEGER(iwp) :: l !< loop variable 709 711 INTEGER(iwp) :: nzt_l !< index of top PE boundary for multigrid level 710 711 712 REAL(wp) :: dx_l !< grid spacing along x on different multigrid level 712 713 REAL(wp) :: dy_l !< grid spacing along y on different multigrid level … … 1091 1092 ! 1092 1093 !-- Allocate arrays for other modules 1093 IF ( biometeorology ) CALL bio_init_arrays1094 1094 IF ( bulk_cloud_model ) CALL bcm_init_arrays 1095 1095 IF ( gust_module_enabled ) CALL gust_init_arrays … … 1098 1098 IF ( salsa ) CALL salsa_init_arrays 1099 1099 IF ( wind_turbine ) CALL wtm_init_arrays 1100 IF ( uv_exposure ) CALL uvem_init_arrays1101 1100 1102 1101 ! … … 1106 1105 ENDIF 1107 1106 1108 !1109 !-- Read uv exposure input data1110 IF ( uv_exposure ) THEN1111 CALL uvem_init1112 ENDIF1113 1107 1114 1108 ! -
palm/trunk/SOURCE/modules.f90
r3543 r3569 25 25 ! ----------------- 26 26 ! $Id$ 27 ! dom_dwd_user, Schrempf: 28 ! -uv_exposure flag, UV model is now part of biometeorology_mod 29 ! 30 ! 3543 2018-11-20 17:06:15Z suehring 27 31 ! +type_x_byte, type_y_byte 28 32 ! … … 1407 1411 LOGICAL :: use_ug_for_galilei_tr = .TRUE. !< namelist parameter 1408 1412 LOGICAL :: use_upstream_for_tke = .FALSE. !< namelist parameter 1409 LOGICAL :: uv_exposure = .FALSE. !< switch for uv exposure model1410 1413 LOGICAL :: virtual_flight = .FALSE. !< use virtual flight model 1411 1414 LOGICAL :: virtual_measurement = .FALSE. !< control parameter to switch-on virtual measurements -
palm/trunk/SOURCE/netcdf_interface_mod.f90
r3529 r3569 25 25 ! ----------------- 26 26 ! $Id$ 27 ! dom_dwd_user, Schrempf: 28 ! Remove uv exposure model code, this is now part of biometeorology_mod 29 ! 30 ! 3529 2018-11-15 21:03:15Z gronemeier 27 31 ! - set time units 28 32 ! - add additional global attributes, … … 646 650 skip_time_do2d_xy, skip_time_do2d_xz, skip_time_do2d_yz, & 647 651 skip_time_do3d, topography, num_leg, num_var_fl, & 648 urban_surface , uv_exposure652 urban_surface 649 653 650 654 USE grid_variables, & … … 694 698 ONLY: usm_define_netcdf_grid 695 699 696 USE uv_exposure_model_mod, &697 ONLY: uvem_define_netcdf_grid698 700 699 701 … … 2871 2873 ENDIF 2872 2874 ! 2873 !-- Check for human thermal comfortquantities2875 !-- Check for biometeorology quantities 2874 2876 IF ( .NOT. found .AND. biometeorology ) THEN 2875 2877 CALL bio_define_netcdf_grid( do2d( av, i), found, & … … 2885 2887 ENDIF 2886 2888 2887 !2888 !-- Check for UV exposure quantities2889 IF ( .NOT. found .AND. uv_exposure ) THEN2890 CALL uvem_define_netcdf_grid( do2d(av,i), found, &2891 grid_x, grid_y, grid_z )2892 ENDIF2893 2889 2894 2890 ! -
palm/trunk/SOURCE/parin.f90
r3545 r3569 25 25 ! ----------------- 26 26 ! $Id$ 27 ! dom_dwd_user, Schrempf: 28 ! Remove uv exposure model code, this is now part of biometeorology_mod. 29 ! 30 ! 3545 2018-11-21 11:19:41Z gronemeier 27 31 ! remove rans_mode from initialization_parameters 28 32 ! … … 561 565 USE urban_surface_mod, & 562 566 ONLY: usm_parin 563 564 USE uv_exposure_model_mod, &565 ONLY: uvem_parin566 567 567 568 USE vertical_nesting_mod, & … … 933 934 CALL stg_parin 934 935 CALL chem_parin 935 CALL uvem_parin936 936 CALL im_parin 937 937 CALL salsa_parin -
palm/trunk/SOURCE/sum_up_3d_data.f90
r3553 r3569 19 19 ! 20 20 ! Current revisions: 21 ! ----------------- 21 ! ------------------ 22 22 ! 23 23 ! … … 25 25 ! ----------------- 26 26 ! $Id$ 27 ! dom_dwd_user, Schrempf: 28 ! Remove CALLs to uv exposure model, this is now part of biometeorology_mod 29 ! 30 ! 3553 2018-11-22 10:30:48Z suehring 27 31 ! variables documented 28 32 ! … … 260 264 ONLY: air_chemistry, average_count_3d, biometeorology, doav, doav_n, & 261 265 land_surface, ocean_mode, rho_surface, urban_surface, & 262 uv_exposure,varnamelength266 varnamelength 263 267 264 268 USE cpulog, & … … 297 301 USE urban_surface_mod, & 298 302 ONLY: usm_average_3d_data 299 300 USE uv_exposure_model_mod, &301 ONLY: uvem_3d_data_averaging302 303 303 304 … … 536 537 ENDIF 537 538 538 IF ( biometeorology .AND. trimvar(1:4) == 'bio_') THEN539 IF ( biometeorology ) THEN 539 540 CALL bio_3d_data_averaging( 'allocate', doav(ii) ) 540 541 ENDIF … … 562 563 IF ( urban_surface .AND. trimvar(1:4) == 'usm_' ) THEN 563 564 CALL usm_average_3d_data( 'allocate', doav(ii) ) 564 ENDIF565 566 IF ( uv_exposure .AND. trimvar(1:5) == 'uvem_') THEN567 CALL uvem_3d_data_averaging( 'allocate', doav(ii) )568 565 ENDIF 569 566 … … 1190 1187 ENDIF 1191 1188 1192 IF ( biometeorology .AND. trimvar(1:4) == 'bio_') THEN1189 IF ( biometeorology ) THEN 1193 1190 CALL bio_3d_data_averaging( 'sum', doav(ii) ) 1194 1191 ENDIF … … 1222 1219 ENDIF 1223 1220 1224 IF ( uv_exposure ) THEN1225 CALL uvem_3d_data_averaging( 'sum', doav(ii) )1226 ENDIF1227 1228 1221 ! 1229 1222 !-- User-defined quantities -
palm/trunk/SOURCE/time_integration.f90
r3525 r3569 25 25 ! ----------------- 26 26 ! $Id$ 27 ! dom_dwd_user, Schrempf: 28 ! Changes due to merge of uv exposure model into biometeorology_mod. 29 ! 30 ! 3525 2018-11-14 16:06:14Z kanani 27 31 ! Changes related to clean-up of biometeorology (dom_dwd_user) 28 32 ! … … 415 419 416 420 USE biometeorology_mod, & 417 ONLY: bio_calculate_thermal_index_maps, time_bio_results 421 ONLY: bio_calculate_thermal_index_maps, time_bio_results, & 422 thermal_comfort, uvem_calc_exposure, uv_exposure 418 423 419 424 USE bulk_cloud_model_mod, & … … 468 473 turbulent_inflow, turbulent_outflow, urban_surface, & 469 474 use_initial_profile_as_reference, & 470 use_single_reference_value, u v_exposure, u_gtrans, v_gtrans,&475 use_single_reference_value, u_gtrans, v_gtrans, & 471 476 virtual_flight, virtual_measurement, wind_turbine, & 472 477 ws_scheme_mom, ws_scheme_sca … … 579 584 ONLY: user_actions 580 585 581 USE uv_exposure_model_mod, &582 ONLY: uvem_calc_exposure583 586 584 587 USE wind_turbine_model_mod, & … … 1254 1257 ENDIF 1255 1258 ENDIF 1256 ! 1257 !-- If required, do UV exposure calculations 1258 IF ( uv_exposure ) THEN 1259 CALL uvem_calc_exposure 1260 ENDIF 1259 1261 1260 ! 1262 1261 !-- If required, calculate indoor temperature, waste heat, heat flux … … 1366 1365 !-- Biometeorology calculation of stationary thermal indices 1367 1366 IF ( biometeorology .AND. time_do3d >= dt_do3d ) THEN 1368 CALL bio_calculate_thermal_index_maps ( .FALSE. ) 1369 time_bio_results = time_since_reference_point 1367 ! 1368 !-- If required, do thermal comfort calculations 1369 IF ( thermal_comfort ) THEN 1370 CALL bio_calculate_thermal_index_maps ( .FALSE. ) 1371 time_bio_results = time_since_reference_point 1372 ENDIF 1373 ! 1374 !-- If required, do UV exposure calculations 1375 IF ( uv_exposure ) THEN 1376 CALL uvem_calc_exposure 1377 ENDIF 1370 1378 ENDIF 1371 1379
Note: See TracChangeset
for help on using the changeset viewer.