Home

Context Navigation

biometeorology_mod.f90 @ 3475

Last change on this file since 3475 was 3475, checked in by kanani, 6 years ago

Add option for using MRT from RTM instead of simplified global value

Property svn:keywords set to Id

Property svn:mergeinfo set to (toggle deleted branches)

/palm/branches/chemistry/SOURCE/biometeorology_mod.f90	2047-3190,3218-3297
/palm/branches/resler/SOURCE/biometeorology_mod.f90	2023-3474
/palm/branches/salsa/SOURCE/biometeorology_mod.f90	2503-3460
/palm/branches/forwind/SOURCE/biometeorology_mod.f90	1564-1913
/palm/branches/fricke/SOURCE/biometeorology_mod.f90	942-977
/palm/branches/hoffmann/SOURCE/biometeorology_mod.f90	989-1052
/palm/branches/letzel/masked_output/SOURCE/biometeorology_mod.f90	296-409
/palm/branches/palm4u/SOURCE/biometeorology_mod.f90	2540-2692
/palm/branches/rans/SOURCE/biometeorology_mod.f90	2078-3128
/palm/branches/suehring/biometeorology_mod.f90	423-666

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

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

Context Navigation

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

Download in other formats: