Home

Context Navigation

source: palm/trunk/SOURCE/biometeorology_mod.f90 @ 3461

Last change on this file since 3461 was 3448, checked in by kanani, 6 years ago
Implementation of human thermal indices (from branch biomet_p2 at r3444)
Property svn:keywords set to `Id`
File size: 49.6 KB

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

Note: See TracBrowser for help on using the repository browser.

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |