Home

Context Navigation

source: palm/trunk/SOURCE/diagnostic_output_quantities_mod.f90 @ 4868

Last change on this file since 4868 was 4866, checked in by suehring, 4 years ago
Implemented vertical passive scalar flux
Property svn:keywords set to `Id`
File size: 97.2 KB

Rev	Line
[3994]	1	!> @file diagnostic_output_quantities_mod.f90
[4583]	2	!--------------------------------------------------------------------------------------------------!
[3994]	3	! This file is part of the PALM model system.
	4	!
[4583]	5	! PALM is free software: you can redistribute it and/or modify it under the terms of the GNU General
	6	! Public License as published by the Free Software Foundation, either version 3 of the License, or
	7	! (at your option) any later version.
[3994]	8	!
[4583]	9	! PALM is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the
	10	! implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General
	11	! Public License for more details.
[3994]	12	!
[4583]	13	! You should have received a copy of the GNU General Public License along with PALM. If not, see
	14	! <http://www.gnu.org/licenses/>.
[3994]	15	!
[4828]	16	! Copyright 1997-2021 Leibniz Universitaet Hannover
[4583]	17	!--------------------------------------------------------------------------------------------------!
[3994]	18	!
	19	! Current revisions:
	20	! ------------------
[4671]	21	!
[4866]	22	!
[3994]	23	! Former revisions:
	24	! -----------------
	25	! $Id: diagnostic_output_quantities_mod.f90 4866 2021-02-05 17:50:14Z raasch $
[4866]	26	! Implemented vertical passive scalar flux
	27	! Implemented air temperature as diagnostic output quantity (previous commit)
	28	!
	29	! 4861 2021-02-01 10:51:45Z raasch
[4768]	30	! Enable 3D data output also with 64-bit precision
	31	!
	32	! 4757 2020-10-26 10:23:38Z schwenkel
[4757]	33	! Implement relative humidity as diagnostic output quantity
	34	!
	35	! 4671 2020-09-09 20:27:58Z pavelkrc
[4671]	36	! Implementation of downward facing USM and LSM surfaces
[4861]	37	!
[4671]	38	! 4583 2020-06-29 12:36:47Z raasch
[4583]	39	! file re-formatted to follow the PALM coding standard
	40	!
	41	! 4560 2020-06-11 12:19:47Z suehring
[4560]	42	! - Bugfix in calculation of vertical momentum and scalar fluxes
	43	! - remove averaged output variables from PUBLIC list
[4583]	44	!
[4560]	45	! 4535 2020-05-15 12:07:23Z raasch
[4535]	46	! bugfix for restart data format query
[4583]	47	!
[4535]	48	! 4518 2020-05-04 15:44:28Z suehring
[4583]	49	! * Define arrays over ghost points in order to allow for standard mpi-io treatment. By this
	50	! modularization of restart-data input is possible with the module interface.
[4518]	51	! * Bugfix: add missing restart input of wtheta_av, wq_av, wu_av, and wv_av.
[4583]	52	!
[4518]	53	! 4517 2020-05-03 14:29:30Z raasch
[4457]	54	! use statement for exchange horiz added,
	55	! bugfix for call of exchange horiz 2d
[4583]	56	!
[4457]	57	! 4431 2020-02-27 23:23:01Z gronemeier
[4431]	58	! added wspeed and wdir output; bugfix: set fill_value in case of masked output
	59	!
	60	! 4360 2020-01-07 11:25:50Z suehring
[4583]	61	! added output of wu, wv, wtheta and wq to enable covariance calculation according to temporal EC
	62	! method
[4431]	63	!
[4351]	64	! 4346 2019-12-18 11:55:56Z motisi
[4583]	65	! Introduction of wall_flags_total_0, which currently sets bits based on static topography
	66	! information used in wall_flags_static_0
[4431]	67	!
[4346]	68	! 4331 2019-12-10 18:25:02Z suehring
[4331]	69	! - Modularize 2-m potential temperature output
	70	! - New output for 10-m wind speed
[4431]	71	!
[4331]	72	! 4329 2019-12-10 15:46:36Z motisi
[4329]	73	! Renamed wall_flags_0 to wall_flags_static_0
[4431]	74	!
[4329]	75	! 4182 2019-08-22 15:20:23Z scharf
[4182]	76	! Corrected "Former revisions" section
[4431]	77	!
[4182]	78	! 4167 2019-08-16 11:01:48Z suehring
[4583]	79	! Changed behaviour of masked output over surface to follow terrain and ignore buildings
	80	! (J.Resler, T.Gronemeier)
[4431]	81	!
[4167]	82	! 4157 2019-08-14 09:19:12Z suehring
[4583]	83	! Initialization restructured, in order to work also when data output during spin-up is enabled.
[4431]	84	!
[4157]	85	! 4132 2019-08-02 12:34:17Z suehring
[4132]	86	! Bugfix in masked data output
[4431]	87	!
[4132]	88	! 4069 2019-07-01 14:05:51Z Giersch
[4583]	89	! Masked output running index mid has been introduced as a local variable to avoid runtime error
	90	! (Loop variable has been modified) in time_integration
[4431]	91	!
[4069]	92	! 4039 2019-06-18 10:32:41Z suehring
[4583]	93	! - Add output of uu, vv, ww to enable variance calculation according temporal EC method
[4039]	94	! - Allocate arrays only when they are required
	95	! - Formatting adjustment
	96	! - Rename subroutines
	97	! - Further modularization
[4431]	98	!
[4039]	99	! 3998 2019-05-23 13:38:11Z suehring
[4431]	100	! Bugfix in gathering all output strings
	101	!
[3998]	102	! 3995 2019-05-22 18:59:54Z suehring
[4583]	103	! Avoid compiler warnings about unused variable and fix string operation which is not allowed with
	104	! PGI compiler
[4431]	105	!
[3995]	106	! 3994 2019-05-22 18:08:09Z suehring
[3994]	107	! Initial revision
	108	!
[4182]	109	! Authors:
	110	! --------
[3994]	111	! @author Farah Kanani-Suehring
[4182]	112	!
[4431]	113	!
[3994]	114	! Description:
	115	! ------------
	116	!> ...
[4583]	117	!--------------------------------------------------------------------------------------------------!
[3994]	118	MODULE diagnostic_output_quantities_mod
[4431]	119
[4583]	120	USE arrays_3d, &
	121	ONLY: ddzu, &
[4757]	122	exner, &
	123	hyp, &
[4583]	124	pt, &
	125	q, &
[4866]	126	s, &
[4583]	127	u, &
	128	v, &
	129	w, &
	130	zu, &
[4331]	131	zw
	132
[4583]	133	USE basic_constants_and_equations_mod, &
[4861]	134	ONLY: degc_to_k, &
	135	kappa, &
[4757]	136	pi, &
	137	magnus, &
	138	rd_d_rv
[4331]	139
[4583]	140	USE control_parameters, &
	141	ONLY: current_timestep_number, &
	142	data_output, &
[4757]	143	data_output_pr, &
	144	humidity, &
[4583]	145	message_string, &
[4866]	146	passive_scalar, &
[4583]	147	restart_data_format_output, &
[4331]	148	varnamelength
[4431]	149	!
[4583]	150	! USE cpulog, &
[3994]	151	! ONLY: cpu_log, log_point
[4039]	152
[4583]	153	USE exchange_horiz_mod, &
[4457]	154	ONLY: exchange_horiz_2d
	155
[4583]	156	USE grid_variables, &
[4039]	157	ONLY: ddx, ddy
[4331]	158
[4583]	159	USE indices, &
	160	ONLY: nbgp, &
	161	nxl, &
	162	nxlg, &
	163	nxr, &
	164	nxrg, &
	165	nyn, &
	166	nyng, &
	167	nys, &
	168	nysg, &
	169	nzb, &
	170	nzt, &
[4346]	171	wall_flags_total_0
[4331]	172
[3994]	173	USE kinds
	174
[4757]	175	USE profil_parameter, &
	176	ONLY: dopr_index
	177
[4518]	178	USE restart_data_mpi_io_mod, &
	179	ONLY: rd_mpi_io_check_array, &
	180	rrd_mpi_io, &
	181	wrd_mpi_io
	182
[4757]	183	USE statistics, &
	184	ONLY: hom, &
	185	pr_palm, &
	186	rmask, &
	187	statistic_regions, &
	188	sums_l
	189
[4583]	190	USE surface_mod, &
	191	ONLY: surf_def_h, &
	192	surf_lsm_h, &
	193	surf_type, &
[4331]	194	surf_usm_h
[3994]	195
	196	IMPLICIT NONE
	197
	198	CHARACTER(LEN=varnamelength), DIMENSION(500) :: do_all = ' '
	199
	200	INTEGER(iwp) :: timestep_number_at_prev_calc = 0 !< ...at previous diagnostic output calculation
[4431]	201
[4039]	202	LOGICAL :: initialized_diagnostic_output_quantities = .FALSE. !< flag indicating whether output is initialized
	203	LOGICAL :: prepared_diagnostic_output_quantities = .FALSE. !< flag indicating whether output is p
[3994]	204
[4331]	205	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: pt_2m !< 2-m air potential temperature
	206	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: pt_2m_av !< averaged 2-m air potential temperature
	207	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: uv_10m !< horizontal wind speed at 10m
	208	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: uv_10m_av !< averaged horizontal wind speed at 10m
	209
[4757]	210	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET :: rh !< relative humidity
	211	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET :: rh_av !< avg. relative humidity
[4861]	212	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET :: ta !< air temperature
	213	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET :: ta_av !< avg. air temperature
[3994]	214	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET :: ti !< rotation(u,v,w) aka turbulence intensity
[4431]	215	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET :: ti_av !< avg. rotation(u,v,w) aka turbulence intensity
	216	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET :: u_center !< u at center of grid box
	217	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET :: u_center_av !< mean of u_center
	218	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET :: uu !< uu
	219	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET :: uu_av !< mean of uu
	220	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET :: wspeed !< horizontal wind speed
	221	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET :: wspeed_av !< mean of horizotal wind speed
	222	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET :: v_center !< v at center of grid box
	223	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET :: v_center_av !< mean of v_center
	224	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET :: vv !< vv
	225	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET :: vv_av !< mean of vv
	226	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET :: wdir !< wind direction
	227	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET :: wdir_av !< mean wind direction
	228	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET :: ww !< ww
	229	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET :: ww_av !< mean of ww
	230	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET :: wu !< wu
	231	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET :: wu_av !< mean of wu
	232	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET :: wv !< wv
	233	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET :: wv_av !< mean of wv
	234	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET :: wtheta !< wtheta
	235	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET :: wtheta_av !< mean of wtheta
	236	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET :: wq !< wq
	237	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET :: wq_av !< mean of wq
[4866]	238	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET :: ws !< ws
	239	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET :: ws_av !< mean of ws
[3994]	240
	241
	242	SAVE
	243
	244	PRIVATE
	245
	246	!
	247	!-- Public variables
[4583]	248	PUBLIC do_all, &
	249	initialized_diagnostic_output_quantities, &
	250	prepared_diagnostic_output_quantities, &
[4560]	251	timestep_number_at_prev_calc
[4431]	252	!
	253	!-- Public routines
[4583]	254	PUBLIC doq_3d_data_averaging, &
	255	doq_calculate, &
	256	doq_check_data_output, &
[4757]	257	doq_check_data_output_pr, &
[4583]	258	doq_define_netcdf_grid, &
	259	doq_init, &
	260	doq_output_2d, &
	261	doq_output_3d, &
	262	doq_output_mask, &
[4757]	263	doq_statistics, &
[4583]	264	doq_rrd_local, &
[4039]	265	doq_wrd_local
[3994]	266
	267
[4039]	268	INTERFACE doq_3d_data_averaging
	269	MODULE PROCEDURE doq_3d_data_averaging
[4431]	270	END INTERFACE doq_3d_data_averaging
[3994]	271
[4039]	272	INTERFACE doq_calculate
	273	MODULE PROCEDURE doq_calculate
	274	END INTERFACE doq_calculate
[3994]	275
[4039]	276	INTERFACE doq_check_data_output
	277	MODULE PROCEDURE doq_check_data_output
	278	END INTERFACE doq_check_data_output
[4431]	279
[4757]	280	INTERFACE doq_check_data_output_pr
	281	MODULE PROCEDURE doq_check_data_output_pr
	282	END INTERFACE doq_check_data_output_pr
	283
[4039]	284	INTERFACE doq_define_netcdf_grid
	285	MODULE PROCEDURE doq_define_netcdf_grid
	286	END INTERFACE doq_define_netcdf_grid
[4431]	287
[4039]	288	INTERFACE doq_output_2d
	289	MODULE PROCEDURE doq_output_2d
	290	END INTERFACE doq_output_2d
[4431]	291
[4039]	292	INTERFACE doq_output_3d
	293	MODULE PROCEDURE doq_output_3d
	294	END INTERFACE doq_output_3d
[4431]	295
[4039]	296	INTERFACE doq_output_mask
	297	MODULE PROCEDURE doq_output_mask
	298	END INTERFACE doq_output_mask
[4431]	299
[4039]	300	INTERFACE doq_init
	301	MODULE PROCEDURE doq_init
	302	END INTERFACE doq_init
[3994]	303
[4757]	304	INTERFACE doq_statistics
	305	MODULE PROCEDURE doq_statistics
	306	END INTERFACE doq_statistics
	307
[4039]	308	INTERFACE doq_prepare
	309	MODULE PROCEDURE doq_prepare
	310	END INTERFACE doq_prepare
[4431]	311
[4518]	312	INTERFACE doq_rrd_local
	313	MODULE PROCEDURE doq_rrd_local_ftn
	314	MODULE PROCEDURE doq_rrd_local_mpi
	315	END INTERFACE doq_rrd_local
[4431]	316
[4039]	317	INTERFACE doq_wrd_local
	318	MODULE PROCEDURE doq_wrd_local
	319	END INTERFACE doq_wrd_local
[3994]	320
	321	CONTAINS
[4431]	322
[4583]	323	!--------------------------------------------------------------------------------------------------!
[4039]	324	! Description:
	325	! ------------
[4583]	326	!> Sum up and time-average diagnostic output quantities as well as allocate the array necessary for
	327	!> storing the average.
	328	!--------------------------------------------------------------------------------------------------!
[4039]	329	SUBROUTINE doq_3d_data_averaging( mode, variable )
[3994]	330
[4583]	331	USE control_parameters, &
[4039]	332	ONLY: average_count_3d
	333
[4431]	334	CHARACTER (LEN=*) :: mode !<
	335	CHARACTER (LEN=*) :: variable !<
[4039]	336
	337	INTEGER(iwp) :: i !<
	338	INTEGER(iwp) :: j !<
	339	INTEGER(iwp) :: k !<
	340
[4583]	341
[4039]	342	IF ( mode == 'allocate' ) THEN
	343
	344	SELECT CASE ( TRIM( variable ) )
	345
[4757]	346	CASE ( 'rh' )
	347	IF ( .NOT. ALLOCATED( rh_av ) ) THEN
	348	ALLOCATE( rh_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
	349	ENDIF
	350	rh_av = 0.0_wp
	351
[4861]	352	CASE ( 'ta' )
	353	IF ( .NOT. ALLOCATED( ta_av ) ) THEN
	354	ALLOCATE( ta_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
	355	ENDIF
	356	ta_av = 0.0_wp
	357
[4039]	358	CASE ( 'ti' )
	359	IF ( .NOT. ALLOCATED( ti_av ) ) THEN
[4518]	360	ALLOCATE( ti_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4039]	361	ENDIF
	362	ti_av = 0.0_wp
[4431]	363
[4039]	364	CASE ( 'uu' )
	365	IF ( .NOT. ALLOCATED( uu_av ) ) THEN
[4518]	366	ALLOCATE( uu_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4039]	367	ENDIF
	368	uu_av = 0.0_wp
[4431]	369
[4039]	370	CASE ( 'vv' )
	371	IF ( .NOT. ALLOCATED( vv_av ) ) THEN
[4518]	372	ALLOCATE( vv_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4039]	373	ENDIF
	374	vv_av = 0.0_wp
[4431]	375
[4039]	376	CASE ( 'ww' )
	377	IF ( .NOT. ALLOCATED( ww_av ) ) THEN
[4518]	378	ALLOCATE( ww_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4039]	379	ENDIF
	380	ww_av = 0.0_wp
[4431]	381
[4351]	382	CASE ( 'wu' )
	383	IF ( .NOT. ALLOCATED( wu_av ) ) THEN
[4518]	384	ALLOCATE( wu_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4351]	385	ENDIF
	386	wu_av = 0.0_wp
[4431]	387
[4351]	388	CASE ( 'wv' )
	389	IF ( .NOT. ALLOCATED( wv_av ) ) THEN
[4518]	390	ALLOCATE( wv_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4351]	391	ENDIF
	392	wv_av = 0.0_wp
[4431]	393
[4351]	394	CASE ( 'wtheta' )
	395	IF ( .NOT. ALLOCATED( wtheta_av ) ) THEN
[4518]	396	ALLOCATE( wtheta_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4351]	397	ENDIF
	398	wtheta_av = 0.0_wp
[4431]	399
[4351]	400	CASE ( 'wq' )
	401	IF ( .NOT. ALLOCATED( wq_av ) ) THEN
[4518]	402	ALLOCATE( wq_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4351]	403	ENDIF
	404	wq_av = 0.0_wp
[4431]	405
[4866]	406	CASE ( 'ws' )
	407	IF ( .NOT. ALLOCATED( ws_av ) ) THEN
	408	ALLOCATE( ws_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
	409	ENDIF
	410	ws_av = 0.0_wp
	411
[4331]	412	CASE ( 'theta_2m*' )
	413	IF ( .NOT. ALLOCATED( pt_2m_av ) ) THEN
	414	ALLOCATE( pt_2m_av(nysg:nyng,nxlg:nxrg) )
	415	ENDIF
	416	pt_2m_av = 0.0_wp
[4431]	417
[4331]	418	CASE ( 'wspeed_10m*' )
	419	IF ( .NOT. ALLOCATED( uv_10m_av ) ) THEN
	420	ALLOCATE( uv_10m_av(nysg:nyng,nxlg:nxrg) )
	421	ENDIF
	422	uv_10m_av = 0.0_wp
[4039]	423
[4431]	424	CASE ( 'wspeed' )
	425	IF ( .NOT. ALLOCATED( wspeed_av ) ) THEN
[4518]	426	ALLOCATE( wspeed_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4431]	427	ENDIF
	428	wspeed_av = 0.0_wp
	429
	430	CASE ( 'wdir' )
	431	IF ( .NOT. ALLOCATED( u_center_av ) ) THEN
[4518]	432	ALLOCATE( u_center_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4431]	433	ENDIF
	434	IF ( .NOT. ALLOCATED( v_center_av ) ) THEN
[4518]	435	ALLOCATE( v_center_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4431]	436	ENDIF
	437	u_center_av = 0.0_wp
	438	v_center_av = 0.0_wp
	439
[4039]	440	CASE DEFAULT
	441	CONTINUE
	442
	443	END SELECT
	444
	445	ELSEIF ( mode == 'sum' ) THEN
	446
	447	SELECT CASE ( TRIM( variable ) )
[4431]	448
[4757]	449	CASE ( 'rh' )
	450	IF ( ALLOCATED( rh_av ) ) THEN
	451	DO i = nxl, nxr
	452	DO j = nys, nyn
	453	DO k = nzb, nzt+1
	454	rh_av(k,j,i) = rh_av(k,j,i) + rh(k,j,i)
	455	ENDDO
	456	ENDDO
	457	ENDDO
	458	ENDIF
	459
[4861]	460	CASE ( 'ta' )
	461	IF ( ALLOCATED( ta_av ) ) THEN
	462	DO i = nxl, nxr
	463	DO j = nys, nyn
	464	DO k = nzb, nzt+1
	465	ta_av(k,j,i) = ta_av(k,j,i) + ta(k,j,i)
	466	ENDDO
	467	ENDDO
	468	ENDDO
	469	ENDIF
	470
[4039]	471	CASE ( 'ti' )
[4583]	472	IF ( ALLOCATED( ti_av ) ) THEN
[4039]	473	DO i = nxl, nxr
	474	DO j = nys, nyn
	475	DO k = nzb, nzt+1
	476	ti_av(k,j,i) = ti_av(k,j,i) + ti(k,j,i)
	477	ENDDO
	478	ENDDO
	479	ENDDO
	480	ENDIF
[4431]	481
[4039]	482	CASE ( 'uu' )
[4583]	483	IF ( ALLOCATED( uu_av ) ) THEN
[4039]	484	DO i = nxl, nxr
	485	DO j = nys, nyn
	486	DO k = nzb, nzt+1
	487	uu_av(k,j,i) = uu_av(k,j,i) + uu(k,j,i)
	488	ENDDO
	489	ENDDO
	490	ENDDO
	491	ENDIF
[4431]	492
[4039]	493	CASE ( 'vv' )
[4583]	494	IF ( ALLOCATED( vv_av ) ) THEN
[4039]	495	DO i = nxl, nxr
	496	DO j = nys, nyn
	497	DO k = nzb, nzt+1
	498	vv_av(k,j,i) = vv_av(k,j,i) + vv(k,j,i)
	499	ENDDO
	500	ENDDO
	501	ENDDO
	502	ENDIF
[4431]	503
[4039]	504	CASE ( 'ww' )
[4583]	505	IF ( ALLOCATED( ww_av ) ) THEN
[4039]	506	DO i = nxl, nxr
	507	DO j = nys, nyn
	508	DO k = nzb, nzt+1
	509	ww_av(k,j,i) = ww_av(k,j,i) + ww(k,j,i)
	510	ENDDO
	511	ENDDO
	512	ENDDO
	513	ENDIF
[4431]	514
[4351]	515	CASE ( 'wu' )
[4583]	516	IF ( ALLOCATED( wu_av ) ) THEN
[4351]	517	DO i = nxl, nxr
	518	DO j = nys, nyn
	519	DO k = nzb, nzt+1
	520	wu_av(k,j,i) = wu_av(k,j,i) + wu(k,j,i)
	521	ENDDO
	522	ENDDO
	523	ENDDO
	524	ENDIF
[4431]	525
[4351]	526	CASE ( 'wv' )
[4583]	527	IF ( ALLOCATED( wv_av ) ) THEN
[4351]	528	DO i = nxl, nxr
	529	DO j = nys, nyn
	530	DO k = nzb, nzt+1
	531	wv_av(k,j,i) = wv_av(k,j,i) + wv(k,j,i)
	532	ENDDO
	533	ENDDO
	534	ENDDO
	535	ENDIF
[4431]	536
[4351]	537	CASE ( 'wtheta' )
[4583]	538	IF ( ALLOCATED( wtheta_av ) ) THEN
[4351]	539	DO i = nxl, nxr
	540	DO j = nys, nyn
	541	DO k = nzb, nzt+1
	542	wtheta_av(k,j,i) = wtheta_av(k,j,i) + wtheta(k,j,i)
	543	ENDDO
	544	ENDDO
	545	ENDDO
	546	ENDIF
[4431]	547
[4351]	548	CASE ( 'wq' )
[4583]	549	IF ( ALLOCATED( wq_av ) ) THEN
[4351]	550	DO i = nxl, nxr
	551	DO j = nys, nyn
	552	DO k = nzb, nzt+1
	553	wq_av(k,j,i) = wq_av(k,j,i) + wq(k,j,i)
	554	ENDDO
	555	ENDDO
	556	ENDDO
	557	ENDIF
[4431]	558
[4866]	559	CASE ( 'ws' )
	560	IF ( ALLOCATED( ws_av ) ) THEN
	561	DO i = nxl, nxr
	562	DO j = nys, nyn
	563	DO k = nzb, nzt+1
	564	ws_av(k,j,i) = ws_av(k,j,i) + ws(k,j,i)
	565	ENDDO
	566	ENDDO
	567	ENDDO
	568	ENDIF
	569
[4331]	570	CASE ( 'theta_2m*' )
[4583]	571	IF ( ALLOCATED( pt_2m_av ) ) THEN
[4331]	572	DO i = nxl, nxr
	573	DO j = nys, nyn
	574	pt_2m_av(j,i) = pt_2m_av(j,i) + pt_2m(j,i)
	575	ENDDO
	576	ENDDO
	577	ENDIF
[4039]	578
[4331]	579	CASE ( 'wspeed_10m*' )
[4583]	580	IF ( ALLOCATED( uv_10m_av ) ) THEN
[4331]	581	DO i = nxl, nxr
	582	DO j = nys, nyn
	583	uv_10m_av(j,i) = uv_10m_av(j,i) + uv_10m(j,i)
	584	ENDDO
	585	ENDDO
	586	ENDIF
	587
[4431]	588	CASE ( 'wspeed' )
[4583]	589	IF ( ALLOCATED( wspeed_av ) ) THEN
[4431]	590	DO i = nxl, nxr
	591	DO j = nys, nyn
	592	DO k = nzb, nzt+1
	593	wspeed_av(k,j,i) = wspeed_av(k,j,i) + wspeed(k,j,i)
	594	ENDDO
	595	ENDDO
	596	ENDDO
	597	ENDIF
	598
	599	CASE ( 'wdir' )
[4583]	600	IF ( ALLOCATED( u_center_av ) .AND. ALLOCATED( v_center_av ) ) THEN
[4431]	601	DO i = nxl, nxr
	602	DO j = nys, nyn
	603	DO k = nzb, nzt+1
	604	u_center_av(k,j,i) = u_center_av(k,j,i) + u_center(k,j,i)
	605	v_center_av(k,j,i) = v_center_av(k,j,i) + v_center(k,j,i)
	606	ENDDO
	607	ENDDO
	608	ENDDO
	609	ENDIF
	610
[4039]	611	CASE DEFAULT
	612	CONTINUE
	613
	614	END SELECT
	615
	616	ELSEIF ( mode == 'average' ) THEN
	617
	618	SELECT CASE ( TRIM( variable ) )
	619
[4757]	620	CASE ( 'rh' )
	621	IF ( ALLOCATED( rh_av ) ) THEN
	622	DO i = nxl, nxr
	623	DO j = nys, nyn
	624	DO k = nzb, nzt+1
	625	rh_av(k,j,i) = rh_av(k,j,i) / REAL( average_count_3d, KIND=wp )
	626	ENDDO
	627	ENDDO
	628	ENDDO
	629	ENDIF
	630
[4861]	631	CASE ( 'ta' )
	632	IF ( ALLOCATED( ta_av ) ) THEN
	633	DO i = nxl, nxr
	634	DO j = nys, nyn
	635	DO k = nzb, nzt+1
	636	ta_av(k,j,i) = ta_av(k,j,i) / REAL( average_count_3d, KIND=wp )
	637	ENDDO
	638	ENDDO
	639	ENDDO
	640	ENDIF
	641
[4039]	642	CASE ( 'ti' )
[4583]	643	IF ( ALLOCATED( ti_av ) ) THEN
[4039]	644	DO i = nxl, nxr
	645	DO j = nys, nyn
	646	DO k = nzb, nzt+1
	647	ti_av(k,j,i) = ti_av(k,j,i) / REAL( average_count_3d, KIND=wp )
	648	ENDDO
	649	ENDDO
	650	ENDDO
	651	ENDIF
[4431]	652
[4039]	653	CASE ( 'uu' )
[4583]	654	IF ( ALLOCATED( uu_av ) ) THEN
[4039]	655	DO i = nxl, nxr
	656	DO j = nys, nyn
	657	DO k = nzb, nzt+1
	658	uu_av(k,j,i) = uu_av(k,j,i) / REAL( average_count_3d, KIND=wp )
	659	ENDDO
	660	ENDDO
	661	ENDDO
	662	ENDIF
[4431]	663
[4039]	664	CASE ( 'vv' )
[4583]	665	IF ( ALLOCATED( vv_av ) ) THEN
[4039]	666	DO i = nxl, nxr
	667	DO j = nys, nyn
	668	DO k = nzb, nzt+1
	669	vv_av(k,j,i) = vv_av(k,j,i) / REAL( average_count_3d, KIND=wp )
	670	ENDDO
	671	ENDDO
	672	ENDDO
	673	ENDIF
[4431]	674
[4039]	675	CASE ( 'ww' )
[4583]	676	IF ( ALLOCATED( ww_av ) ) THEN
[4039]	677	DO i = nxl, nxr
	678	DO j = nys, nyn
	679	DO k = nzb, nzt+1
	680	ww_av(k,j,i) = ww_av(k,j,i) / REAL( average_count_3d, KIND=wp )
	681	ENDDO
	682	ENDDO
	683	ENDDO
	684	ENDIF
	685
[4351]	686	CASE ( 'wu' )
[4583]	687	IF ( ALLOCATED( wu_av ) ) THEN
[4351]	688	DO i = nxl, nxr
	689	DO j = nys, nyn
	690	DO k = nzb, nzt+1
	691	wu_av(k,j,i) = wu_av(k,j,i) / REAL( average_count_3d, KIND=wp )
	692	ENDDO
	693	ENDDO
	694	ENDDO
	695	ENDIF
[4431]	696
[4351]	697	CASE ( 'wv' )
[4583]	698	IF ( ALLOCATED( wv_av ) ) THEN
[4351]	699	DO i = nxl, nxr
	700	DO j = nys, nyn
	701	DO k = nzb, nzt+1
	702	wv_av(k,j,i) = wv_av(k,j,i) / REAL( average_count_3d, KIND=wp )
	703	ENDDO
	704	ENDDO
	705	ENDDO
	706	ENDIF
[4431]	707
[4351]	708	CASE ( 'wtheta' )
[4583]	709	IF ( ALLOCATED( wtheta_av ) ) THEN
[4351]	710	DO i = nxl, nxr
	711	DO j = nys, nyn
	712	DO k = nzb, nzt+1
	713	wtheta_av(k,j,i) = wtheta_av(k,j,i) / REAL( average_count_3d, KIND=wp )
	714	ENDDO
	715	ENDDO
	716	ENDDO
	717	ENDIF
[4431]	718
[4351]	719	CASE ( 'wq' )
[4583]	720	IF ( ALLOCATED( wq_av ) ) THEN
[4351]	721	DO i = nxl, nxr
	722	DO j = nys, nyn
	723	DO k = nzb, nzt+1
	724	wq_av(k,j,i) = wq_av(k,j,i) / REAL( average_count_3d, KIND=wp )
	725	ENDDO
	726	ENDDO
	727	ENDDO
	728	ENDIF
	729
[4866]	730	CASE ( 'ws' )
	731	IF ( ALLOCATED( ws_av ) ) THEN
	732	DO i = nxl, nxr
	733	DO j = nys, nyn
	734	DO k = nzb, nzt+1
	735	ws_av(k,j,i) = ws_av(k,j,i) / REAL( average_count_3d, KIND=wp )
	736	ENDDO
	737	ENDDO
	738	ENDDO
	739	ENDIF
	740
[4331]	741	CASE ( 'theta_2m*' )
[4583]	742	IF ( ALLOCATED( pt_2m_av ) ) THEN
[4331]	743	DO i = nxlg, nxrg
	744	DO j = nysg, nyng
	745	pt_2m_av(j,i) = pt_2m_av(j,i) / REAL( average_count_3d, KIND=wp )
	746	ENDDO
	747	ENDDO
[4457]	748	CALL exchange_horiz_2d( pt_2m_av )
[4331]	749	ENDIF
	750
	751	CASE ( 'wspeed_10m*' )
[4583]	752	IF ( ALLOCATED( uv_10m_av ) ) THEN
[4331]	753	DO i = nxlg, nxrg
	754	DO j = nysg, nyng
	755	uv_10m_av(j,i) = uv_10m_av(j,i) / REAL( average_count_3d, KIND=wp )
	756	ENDDO
	757	ENDDO
[4457]	758	CALL exchange_horiz_2d( uv_10m_av )
[4331]	759	ENDIF
	760
[4431]	761	CASE ( 'wspeed' )
[4583]	762	IF ( ALLOCATED( wspeed_av ) ) THEN
[4431]	763	DO i = nxl, nxr
	764	DO j = nys, nyn
	765	DO k = nzb, nzt+1
	766	wspeed_av(k,j,i) = wspeed_av(k,j,i) / REAL( average_count_3d, KIND=wp )
	767	ENDDO
	768	ENDDO
	769	ENDDO
	770	ENDIF
	771
	772	CASE ( 'wdir' )
[4583]	773	IF ( ALLOCATED( u_center_av ) .AND. ALLOCATED( v_center_av ) ) THEN
[4431]	774
	775	IF ( .NOT. ALLOCATED( wdir_av ) ) THEN
[4518]	776	ALLOCATE( wdir_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4431]	777	ENDIF
	778	wdir_av = 0.0_wp
	779
	780	DO i = nxl, nxr
	781	DO j = nys, nyn
	782	DO k = nzb, nzt+1
	783	u_center_av(k,j,i) = u_center_av(k,j,i) / REAL( average_count_3d, KIND=wp )
	784	v_center_av(k,j,i) = v_center_av(k,j,i) / REAL( average_count_3d, KIND=wp )
[4583]	785	wdir_av(k,j,i) = ATAN2( u_center_av(k,j,i), v_center_av(k,j,i) ) &
	786	/ pi * 180.0_wp + 180.0_wp
[4431]	787	ENDDO
	788	ENDDO
	789	ENDDO
	790	ENDIF
	791
[4039]	792	END SELECT
	793
	794	ENDIF
	795
	796
[4431]	797	END SUBROUTINE doq_3d_data_averaging
	798
[4583]	799	!--------------------------------------------------------------------------------------------------!
[3994]	800	! Description:
	801	! ------------
[4039]	802	!> Check data output for diagnostic output
[4583]	803	!--------------------------------------------------------------------------------------------------!
[4331]	804	SUBROUTINE doq_check_data_output( var, unit, i, ilen, k )
[3994]	805
[4039]	806	IMPLICIT NONE
[3994]	807
[4431]	808	CHARACTER (LEN=*) :: unit !<
[4039]	809	CHARACTER (LEN=*) :: var !<
[3994]	810
[4331]	811	INTEGER(iwp), OPTIONAL, INTENT(IN) :: i !< Current element of data_output
	812	INTEGER(iwp), OPTIONAL, INTENT(IN) :: ilen !< Length of current entry in data_output
	813	INTEGER(iwp), OPTIONAL, INTENT(IN) :: k !< Output is xy mode? 0 = no, 1 = yes
	814
[4583]	815
[4039]	816	SELECT CASE ( TRIM( var ) )
	817
[4757]	818	CASE ( 'rh' )
	819	unit = '%'
	820
[4861]	821	CASE ( 'ta' )
	822	unit = 'degree_C'
	823
[4039]	824	CASE ( 'ti' )
	825	unit = '1/s'
[4431]	826
[4039]	827	CASE ( 'uu' )
	828	unit = 'm2/s2'
[4431]	829
[4039]	830	CASE ( 'vv' )
	831	unit = 'm2/s2'
[4431]	832
[4039]	833	CASE ( 'ww' )
	834	unit = 'm2/s2'
[4431]	835
[4351]	836	CASE ( 'wu' )
	837	unit = 'm2/s2'
[4431]	838
[4351]	839	CASE ( 'wv' )
	840	unit = 'm2/s2'
[4431]	841
[4351]	842	CASE ( 'wtheta' )
	843	unit = 'Km/s'
[4431]	844
[4351]	845	CASE ( 'wq' )
	846	unit = 'm/s'
[4431]	847
[4866]	848	CASE ( 'ws' )
	849	IF ( .NOT. passive_scalar ) THEN
	850	message_string = 'data_output = ' // TRIM( var ) // &
	851	' is not implemented for passive_scalar = .FALSE.'
	852	CALL message( 'diagnostic_output', 'PA0093', 1, 2, 0, 6, 0 )
	853	ENDIF
	854	unit = 'm/s'
	855
[4431]	856	CASE ( 'wspeed' )
	857	unit = 'm/s'
	858
	859	CASE ( 'wdir' )
	860	unit = 'degree'
[4331]	861	!
	862	!-- Treat horizotal cross-section output quanatities
	863	CASE ( 'theta_2m', 'wspeed_10m' )
	864	!
	865	!-- Check if output quantity is _xy only.
	866	IF ( k == 0 .OR. data_output(i)(ilen-2:ilen) /= '_xy' ) THEN
[4583]	867	message_string = 'illegal value for data_output: "' // &
	868	TRIM( var ) // '" & only 2d-horizontal ' // &
[4331]	869	'cross sections are allowed for this value'
	870	CALL message( 'diagnostic_output', 'PA0111', 1, 2, 0, 6, 0 )
	871	ENDIF
	872
	873	IF ( TRIM( var ) == 'theta_2m*' ) unit = 'K'
	874	IF ( TRIM( var ) == 'wspeed_10m*' ) unit = 'm/s'
	875
[4039]	876	CASE DEFAULT
	877	unit = 'illegal'
	878
	879	END SELECT
	880
	881
	882	END SUBROUTINE doq_check_data_output
[4431]	883
[4757]	884
[4583]	885	!--------------------------------------------------------------------------------------------------!
[4757]	886	! Description:
	887	! ------------
	888	!> Set the unit of user defined profile output quantities. For those variables not recognized by the
	889	!> user, the parameter unit is set to "illegal", which tells the calling routine that the
	890	!> output variable is not defined and leads to a program abort.
	891	!--------------------------------------------------------------------------------------------------!
	892	SUBROUTINE doq_check_data_output_pr( variable, var_count, unit, dopr_unit )
	893
	894	CHARACTER (LEN=*) :: unit !<
	895	CHARACTER (LEN=*) :: variable !<
	896	CHARACTER (LEN=*) :: dopr_unit !< local value of dopr_unit
	897
	898	INTEGER(iwp) :: pr_index !<
	899	INTEGER(iwp) :: var_count !<
	900
	901	SELECT CASE ( TRIM( variable ) )
	902
	903	CASE ( 'rh' )
	904	IF ( .NOT. humidity ) THEN
	905	message_string = 'data_output_pr = ' // TRIM( data_output_pr(var_count) ) // &
	906	' requires humidity'
	907	CALL message( 'check_parameters', 'PA0358', 1, 2, 0, 6, 0 )
	908	ENDIF
	909	pr_index = 130
	910	dopr_index(var_count) = pr_index
	911	dopr_unit = '%'
	912	unit = dopr_unit
	913	hom(:,2,pr_index,:) = SPREAD( zu, 2, statistic_regions+1 )
	914
	915	CASE DEFAULT
	916	unit = 'illegal'
	917
	918	END SELECT
	919
	920
	921	END SUBROUTINE doq_check_data_output_pr
	922
	923
	924	!--------------------------------------------------------------------------------------------------!
[4039]	925	!
	926	! Description:
	927	! ------------
	928	!> Subroutine defining appropriate grid for netcdf variables.
[4583]	929	!--------------------------------------------------------------------------------------------------!
[4039]	930	SUBROUTINE doq_define_netcdf_grid( variable, found, grid_x, grid_y, grid_z )
[4431]	931
[3994]	932	IMPLICIT NONE
[4039]	933
	934	CHARACTER (LEN=*), INTENT(IN) :: variable !<
[4431]	935	CHARACTER (LEN=*), INTENT(OUT) :: grid_x !<
	936	CHARACTER (LEN=*), INTENT(OUT) :: grid_y !<
	937	CHARACTER (LEN=*), INTENT(OUT) :: grid_z !<
[4039]	938
[4583]	939	LOGICAL, INTENT(OUT) :: found !<
	940
	941
[4039]	942	found = .TRUE.
[4431]	943
[4039]	944	SELECT CASE ( TRIM( variable ) )
[3995]	945	!
[4039]	946	!-- s grid
[4757]	947	CASE ( 'rh', 'rh_xy', 'rh_xz', 'rh_yz', &
[4861]	948	'ta', 'ta_xy', 'ta_xz', 'ta_yz', &
[4757]	949	'ti', 'ti_xy', 'ti_xz', 'ti_yz', &
[4583]	950	'wspeed', 'wspeed_xy', 'wspeed_xz', 'wspeed_yz', &
[4560]	951	'wdir', 'wdir_xy', 'wdir_xz', 'wdir_yz' )
[3994]	952
[4039]	953	grid_x = 'x'
	954	grid_y = 'y'
[4331]	955	grid_z = 'zu'
[4039]	956	!
[4331]	957	!-- s grid surface variables
	958	CASE ( 'theta_2m_xy', 'wspeed_10m' )
	959
	960	grid_x = 'x'
	961	grid_y = 'y'
	962	grid_z = 'zu'
	963	!
[4039]	964	!-- u grid
	965	CASE ( 'uu', 'uu_xy', 'uu_xz', 'uu_yz' )
[3994]	966
[4039]	967	grid_x = 'xu'
	968	grid_y = 'y'
	969	grid_z = 'zu'
	970	!
	971	!-- v grid
	972	CASE ( 'vv', 'vv_xy', 'vv_xz', 'vv_yz' )
	973
	974	grid_x = 'x'
	975	grid_y = 'yv'
	976	grid_z = 'zu'
	977	!
	978	!-- w grid
[4583]	979	CASE ( 'ww', 'ww_xy', 'ww_xz', 'ww_yz', &
	980	'wu', 'wu_xy', 'wu_xz', 'wu_yz', &
	981	'wv', 'wv_xy', 'wv_xz', 'wv_yz', &
	982	'wtheta', 'wtheta_xy', 'wtheta_xz', 'wtheta_yz', &
[4866]	983	'wq', 'wq_xy', 'wq_xz', 'wq_yz', &
	984	'ws', 'ws_xy', 'ws_xz', 'ws_yz' )
[4039]	985
	986	grid_x = 'x'
	987	grid_y = 'y'
	988	grid_z = 'zw'
	989
	990	CASE DEFAULT
	991	found = .FALSE.
	992	grid_x = 'none'
	993	grid_y = 'none'
	994	grid_z = 'none'
	995
	996	END SELECT
	997
	998
	999	END SUBROUTINE doq_define_netcdf_grid
[4431]	1000
[4583]	1001	!--------------------------------------------------------------------------------------------------!
[4039]	1002	!
	1003	! Description:
	1004	! ------------
	1005	!> Subroutine defining 2D output variables
[4583]	1006	!--------------------------------------------------------------------------------------------------!
	1007	SUBROUTINE doq_output_2d( av, variable, found, grid, mode, local_pf, two_d, nzb_do, nzt_do, &
	1008	fill_value )
[4039]	1009
	1010
	1011	IMPLICIT NONE
	1012
[4431]	1013	CHARACTER (LEN=*) :: grid !<
	1014	CHARACTER (LEN=*) :: mode !<
	1015	CHARACTER (LEN=*) :: variable !<
[4039]	1016
	1017	INTEGER(iwp) :: av !< value indicating averaged or non-averaged output
	1018	INTEGER(iwp) :: flag_nr !< number of the topography flag (0: scalar, 1: u, 2: v, 3: w)
	1019	INTEGER(iwp) :: i !< grid index x-direction
	1020	INTEGER(iwp) :: j !< grid index y-direction
	1021	INTEGER(iwp) :: k !< grid index z-direction
[4431]	1022	INTEGER(iwp) :: nzb_do !<
	1023	INTEGER(iwp) :: nzt_do !<
[4039]	1024
	1025	LOGICAL :: found !< true if variable is in list
	1026	LOGICAL :: resorted !< true if array is resorted
	1027	LOGICAL :: two_d !< flag parameter that indicates 2D variables (horizontal cross sections)
	1028
	1029	REAL(wp) :: fill_value !< value for the _FillValue attribute
[4431]	1030
	1031	REAL(wp), DIMENSION(nxl:nxr,nys:nyn,nzb_do:nzt_do) :: local_pf !<
[4039]	1032	REAL(wp), DIMENSION(:,:,:), POINTER :: to_be_resorted !< points to array which needs to be resorted for output
[4431]	1033
[4583]	1034
[4039]	1035	flag_nr = 0
	1036	found = .TRUE.
	1037	resorted = .FALSE.
	1038	two_d = .FALSE.
	1039
	1040	SELECT CASE ( TRIM( variable ) )
	1041
[4757]	1042	CASE ( 'rh_xy', 'rh_xz', 'rh_yz' )
	1043	IF ( av == 0 ) THEN
	1044	to_be_resorted => rh
	1045	ELSE
	1046	IF ( .NOT. ALLOCATED( rh_av ) ) THEN
	1047	ALLOCATE( rh_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
	1048	rh_av = REAL( fill_value, KIND = wp )
	1049	ENDIF
	1050	to_be_resorted => rh_av
	1051	ENDIF
	1052	flag_nr = 0
	1053
	1054	IF ( mode == 'xy' ) grid = 'zu'
	1055
[4861]	1056	CASE ( 'ta_xy', 'ta_xz', 'ta_yz' )
	1057	IF ( av == 0 ) THEN
	1058	to_be_resorted => ta
	1059	ELSE
	1060	IF ( .NOT. ALLOCATED( ta_av ) ) THEN
	1061	ALLOCATE( ta_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
	1062	ta_av = REAL( fill_value, KIND = wp )
	1063	ENDIF
	1064	to_be_resorted => ta_av
	1065	ENDIF
	1066	flag_nr = 0
	1067
	1068	IF ( mode == 'xy' ) grid = 'zu'
	1069
[4039]	1070	CASE ( 'ti_xy', 'ti_xz', 'ti_yz' )
	1071	IF ( av == 0 ) THEN
	1072	to_be_resorted => ti
	1073	ELSE
[4583]	1074	IF ( .NOT. ALLOCATED( ti_av ) ) THEN
[4518]	1075	ALLOCATE( ti_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4039]	1076	ti_av = REAL( fill_value, KIND = wp )
	1077	ENDIF
	1078	to_be_resorted => ti_av
	1079	ENDIF
	1080	flag_nr = 0
[4431]	1081
[4039]	1082	IF ( mode == 'xy' ) grid = 'zu'
[4431]	1083
[4039]	1084	CASE ( 'uu_xy', 'uu_xz', 'uu_yz' )
	1085	IF ( av == 0 ) THEN
	1086	to_be_resorted => uu
	1087	ELSE
[4583]	1088	IF ( .NOT. ALLOCATED( uu_av ) ) THEN
[4518]	1089	ALLOCATE( uu_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4039]	1090	uu_av = REAL( fill_value, KIND = wp )
	1091	ENDIF
	1092	to_be_resorted => uu_av
	1093	ENDIF
	1094	flag_nr = 1
[4431]	1095
[4039]	1096	IF ( mode == 'xy' ) grid = 'zu'
[4431]	1097
[4039]	1098	CASE ( 'vv_xy', 'vv_xz', 'vv_yz' )
	1099	IF ( av == 0 ) THEN
	1100	to_be_resorted => vv
	1101	ELSE
[4583]	1102	IF ( .NOT. ALLOCATED( vv_av ) ) THEN
[4518]	1103	ALLOCATE( vv_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4039]	1104	vv_av = REAL( fill_value, KIND = wp )
	1105	ENDIF
	1106	to_be_resorted => vv_av
	1107	ENDIF
	1108	flag_nr = 2
[4431]	1109
[4039]	1110	IF ( mode == 'xy' ) grid = 'zu'
[4431]	1111
[4039]	1112	CASE ( 'ww_xy', 'ww_xz', 'ww_yz' )
	1113	IF ( av == 0 ) THEN
	1114	to_be_resorted => ww
	1115	ELSE
[4583]	1116	IF ( .NOT. ALLOCATED( ww_av ) ) THEN
[4518]	1117	ALLOCATE( ww_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4039]	1118	ww_av = REAL( fill_value, KIND = wp )
	1119	ENDIF
	1120	to_be_resorted => ww_av
	1121	ENDIF
	1122	flag_nr = 3
[4431]	1123
[4039]	1124	IF ( mode == 'xy' ) grid = 'zw'
[4431]	1125
[4351]	1126	CASE ( 'wu_xy', 'wu_xz', 'wu_yz' )
	1127	IF ( av == 0 ) THEN
	1128	to_be_resorted => wu
	1129	ELSE
[4583]	1130	IF ( .NOT. ALLOCATED( wu_av ) ) THEN
[4518]	1131	ALLOCATE( wu_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4351]	1132	wu_av = REAL( fill_value, KIND = wp )
	1133	ENDIF
	1134	to_be_resorted => wu_av
	1135	ENDIF
	1136	flag_nr = 0
[4431]	1137
[4351]	1138	IF ( mode == 'xy' ) grid = 'zw'
[4431]	1139
[4351]	1140	CASE ( 'wv_xy', 'wv_xz', 'wv_yz' )
	1141	IF ( av == 0 ) THEN
	1142	to_be_resorted => wv
	1143	ELSE
[4583]	1144	IF ( .NOT. ALLOCATED( wv_av ) ) THEN
[4518]	1145	ALLOCATE( wv_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4351]	1146	wv_av = REAL( fill_value, KIND = wp )
	1147	ENDIF
	1148	to_be_resorted => wv_av
	1149	ENDIF
	1150	flag_nr = 0
[4431]	1151
[4351]	1152	IF ( mode == 'xy' ) grid = 'zw'
[4431]	1153
[4351]	1154	CASE ( 'wtheta_xy', 'wtheta_xz', 'wtheta_yz' )
	1155	IF ( av == 0 ) THEN
	1156	to_be_resorted => wtheta
	1157	ELSE
[4583]	1158	IF ( .NOT. ALLOCATED( wtheta_av ) ) THEN
[4518]	1159	ALLOCATE( wtheta_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4351]	1160	wtheta_av = REAL( fill_value, KIND = wp )
	1161	ENDIF
	1162	to_be_resorted => wtheta_av
	1163	ENDIF
	1164	flag_nr = 0
[4431]	1165
[4351]	1166	IF ( mode == 'xy' ) grid = 'zw'
[4431]	1167
[4351]	1168	CASE ( 'wq_xy', 'wq_xz', 'wq_yz' )
	1169	IF ( av == 0 ) THEN
	1170	to_be_resorted => wq
	1171	ELSE
[4583]	1172	IF ( .NOT. ALLOCATED( wq_av ) ) THEN
[4518]	1173	ALLOCATE( wq_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4351]	1174	wq_av = REAL( fill_value, KIND = wp )
	1175	ENDIF
	1176	to_be_resorted => wq_av
	1177	ENDIF
	1178	flag_nr = 0
[4431]	1179
[4351]	1180	IF ( mode == 'xy' ) grid = 'zw'
[4431]	1181
[4866]	1182	CASE ( 'ws_xy', 'ws_xz', 'ws_yz' )
	1183	IF ( av == 0 ) THEN
	1184	to_be_resorted => ws
	1185	ELSE
	1186	IF ( .NOT. ALLOCATED( ws_av ) ) THEN
	1187	ALLOCATE( ws_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
	1188	ws_av = REAL( fill_value, KIND = wp )
	1189	ENDIF
	1190	to_be_resorted => ws_av
	1191	ENDIF
	1192	flag_nr = 0
	1193
	1194	IF ( mode == 'xy' ) grid = 'zw'
	1195
[4331]	1196	CASE ( 'theta_2m*_xy' ) ! 2d-array
	1197	IF ( av == 0 ) THEN
	1198	DO i = nxl, nxr
	1199	DO j = nys, nyn
	1200	local_pf(i,j,nzb+1) = pt_2m(j,i)
	1201	ENDDO
	1202	ENDDO
	1203	ELSE
[4583]	1204	IF ( .NOT. ALLOCATED( pt_2m_av ) ) THEN
[4331]	1205	ALLOCATE( pt_2m_av(nysg:nyng,nxlg:nxrg) )
	1206	pt_2m_av = REAL( fill_value, KIND = wp )
	1207	ENDIF
	1208	DO i = nxl, nxr
	1209	DO j = nys, nyn
	1210	local_pf(i,j,nzb+1) = pt_2m_av(j,i)
	1211	ENDDO
	1212	ENDDO
	1213	ENDIF
	1214	resorted = .TRUE.
	1215	two_d = .TRUE.
	1216	grid = 'zu1'
[4431]	1217
[4331]	1218	CASE ( 'wspeed_10m*_xy' ) ! 2d-array
	1219	IF ( av == 0 ) THEN
	1220	DO i = nxl, nxr
	1221	DO j = nys, nyn
	1222	local_pf(i,j,nzb+1) = uv_10m(j,i)
	1223	ENDDO
	1224	ENDDO
	1225	ELSE
[4583]	1226	IF ( .NOT. ALLOCATED( uv_10m_av ) ) THEN
[4331]	1227	ALLOCATE( uv_10m_av(nysg:nyng,nxlg:nxrg) )
	1228	uv_10m_av = REAL( fill_value, KIND = wp )
	1229	ENDIF
	1230	DO i = nxl, nxr
	1231	DO j = nys, nyn
	1232	local_pf(i,j,nzb+1) = uv_10m_av(j,i)
	1233	ENDDO
	1234	ENDDO
	1235	ENDIF
	1236	resorted = .TRUE.
	1237	two_d = .TRUE.
	1238	grid = 'zu1'
[4039]	1239
[4431]	1240	CASE ( 'wspeed_xy', 'wspeed_xz', 'wspeed_yz' )
	1241	IF ( av == 0 ) THEN
	1242	to_be_resorted => wspeed
	1243	ELSE
[4583]	1244	IF ( .NOT. ALLOCATED( wspeed_av ) ) THEN
[4518]	1245	ALLOCATE( wspeed_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4431]	1246	wspeed_av = REAL( fill_value, KIND = wp )
	1247	ENDIF
	1248	to_be_resorted => wspeed_av
	1249	ENDIF
	1250	flag_nr = 0
	1251
	1252	IF ( mode == 'xy' ) grid = 'zu'
	1253
	1254	CASE ( 'wdir_xy', 'wdir_xz', 'wdir_yz' )
	1255	IF ( av == 0 ) THEN
	1256	to_be_resorted => wdir
	1257	ELSE
[4583]	1258	IF ( .NOT. ALLOCATED( wdir_av ) ) THEN
[4518]	1259	ALLOCATE( wdir_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4431]	1260	wdir_av = REAL( fill_value, KIND = wp )
	1261	ENDIF
	1262	to_be_resorted => wdir_av
	1263	ENDIF
	1264	flag_nr = 0
	1265
	1266	IF ( mode == 'xy' ) grid = 'zu'
	1267
[4039]	1268	CASE DEFAULT
	1269	found = .FALSE.
	1270	grid = 'none'
	1271
	1272	END SELECT
[4431]	1273
	1274	IF ( found .AND. .NOT. resorted ) THEN
[4039]	1275	DO i = nxl, nxr
	1276	DO j = nys, nyn
	1277	DO k = nzb_do, nzt_do
[4583]	1278	local_pf(i,j,k) = MERGE( to_be_resorted(k,j,i), &
	1279	REAL( fill_value, KIND = wp ), &
	1280	BTEST( wall_flags_total_0(k,j,i), flag_nr ) )
[4039]	1281	ENDDO
	1282	ENDDO
	1283	ENDDO
[3994]	1284	ENDIF
[4431]	1285
[4039]	1286	END SUBROUTINE doq_output_2d
[4431]	1287
	1288
[4583]	1289	!--------------------------------------------------------------------------------------------------!
[4039]	1290	!
	1291	! Description:
	1292	! ------------
	1293	!> Subroutine defining 3D output variables
[4583]	1294	!--------------------------------------------------------------------------------------------------!
	1295	SUBROUTINE doq_output_3d( av, variable, found, local_pf, fill_value, nzb_do, nzt_do )
[4431]	1296
[4039]	1297	IMPLICIT NONE
[3994]	1298
[4431]	1299	CHARACTER (LEN=*) :: variable !<
[3994]	1300
[4039]	1301	INTEGER(iwp) :: av !< index indicating averaged or instantaneous output
	1302	INTEGER(iwp) :: flag_nr !< number of the topography flag (0: scalar, 1: u, 2: v, 3: w)
	1303	INTEGER(iwp) :: i !< index variable along x-direction
	1304	INTEGER(iwp) :: j !< index variable along y-direction
	1305	INTEGER(iwp) :: k !< index variable along z-direction
	1306	INTEGER(iwp) :: nzb_do !< lower limit of the data output (usually 0)
	1307	INTEGER(iwp) :: nzt_do !< vertical upper limit of the data output (usually nz_do3d)
[3994]	1308
[4039]	1309	LOGICAL :: found !< true if variable is in list
	1310	LOGICAL :: resorted !< true if array is resorted
[3994]	1311
[4039]	1312	REAL(wp) :: fill_value !< value for the _FillValue attribute
	1313
[4768]	1314	REAL(wp), DIMENSION(nxl:nxr,nys:nyn,nzb_do:nzt_do) :: local_pf !<
[4039]	1315	REAL(wp), DIMENSION(:,:,:), POINTER :: to_be_resorted !< points to array which needs to be resorted for output
	1316
[4583]	1317
[4039]	1318	flag_nr = 0
	1319	found = .TRUE.
	1320	resorted = .FALSE.
[4431]	1321
[4039]	1322	SELECT CASE ( TRIM( variable ) )
	1323
[4757]	1324	CASE ( 'rh' )
	1325	IF ( av == 0 ) THEN
	1326	to_be_resorted => rh
	1327	ELSE
	1328	IF ( .NOT. ALLOCATED( rh ) ) THEN
	1329	ALLOCATE( rh(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
	1330	rh = REAL( fill_value, KIND = wp )
	1331	ENDIF
	1332	to_be_resorted => rh
	1333	ENDIF
	1334	flag_nr = 0
	1335
[4861]	1336	CASE ( 'ta' )
	1337	IF ( av == 0 ) THEN
	1338	to_be_resorted => ta
	1339	ELSE
	1340	IF ( .NOT. ALLOCATED( ta ) ) THEN
	1341	ALLOCATE( ta(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
	1342	ta = REAL( fill_value, KIND = wp )
	1343	ENDIF
	1344	to_be_resorted => ta
	1345	ENDIF
	1346	flag_nr = 0
	1347
[4039]	1348	CASE ( 'ti' )
	1349	IF ( av == 0 ) THEN
	1350	to_be_resorted => ti
	1351	ELSE
[4583]	1352	IF ( .NOT. ALLOCATED( ti_av ) ) THEN
[4518]	1353	ALLOCATE( ti_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4039]	1354	ti_av = REAL( fill_value, KIND = wp )
	1355	ENDIF
	1356	to_be_resorted => ti_av
	1357	ENDIF
	1358	flag_nr = 0
[4431]	1359
[4039]	1360	CASE ( 'uu' )
	1361	IF ( av == 0 ) THEN
	1362	to_be_resorted => uu
	1363	ELSE
[4583]	1364	IF ( .NOT. ALLOCATED( uu_av ) ) THEN
[4518]	1365	ALLOCATE( uu_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4039]	1366	uu_av = REAL( fill_value, KIND = wp )
	1367	ENDIF
	1368	to_be_resorted => uu_av
	1369	ENDIF
	1370	flag_nr = 1
[4431]	1371
[4039]	1372	CASE ( 'vv' )
	1373	IF ( av == 0 ) THEN
	1374	to_be_resorted => vv
	1375	ELSE
[4583]	1376	IF ( .NOT. ALLOCATED( vv_av ) ) THEN
[4518]	1377	ALLOCATE( vv_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4039]	1378	vv_av = REAL( fill_value, KIND = wp )
	1379	ENDIF
	1380	to_be_resorted => vv_av
	1381	ENDIF
	1382	flag_nr = 2
[4431]	1383
[4039]	1384	CASE ( 'ww' )
	1385	IF ( av == 0 ) THEN
	1386	to_be_resorted => ww
	1387	ELSE
[4583]	1388	IF ( .NOT. ALLOCATED( ww_av ) ) THEN
[4518]	1389	ALLOCATE( ww_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4039]	1390	ww_av = REAL( fill_value, KIND = wp )
	1391	ENDIF
	1392	to_be_resorted => ww_av
	1393	ENDIF
	1394	flag_nr = 3
	1395
[4351]	1396	CASE ( 'wu' )
	1397	IF ( av == 0 ) THEN
	1398	to_be_resorted => wu
	1399	ELSE
[4583]	1400	IF ( .NOT. ALLOCATED( wu_av ) ) THEN
[4518]	1401	ALLOCATE( wu_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4351]	1402	wu_av = REAL( fill_value, KIND = wp )
	1403	ENDIF
	1404	to_be_resorted => wu_av
	1405	ENDIF
	1406	flag_nr = 0
	1407
	1408	CASE ( 'wv' )
	1409	IF ( av == 0 ) THEN
	1410	to_be_resorted => wv
	1411	ELSE
[4583]	1412	IF ( .NOT. ALLOCATED( wv_av ) ) THEN
[4518]	1413	ALLOCATE( wv_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4351]	1414	wv_av = REAL( fill_value, KIND = wp )
	1415	ENDIF
	1416	to_be_resorted => wv_av
	1417	ENDIF
	1418	flag_nr = 0
[4431]	1419
[4351]	1420	CASE ( 'wtheta' )
	1421	IF ( av == 0 ) THEN
	1422	to_be_resorted => wtheta
	1423	ELSE
[4583]	1424	IF ( .NOT. ALLOCATED( wtheta_av ) ) THEN
[4518]	1425	ALLOCATE( wtheta_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4351]	1426	wtheta_av = REAL( fill_value, KIND = wp )
	1427	ENDIF
	1428	to_be_resorted => wtheta_av
	1429	ENDIF
	1430	flag_nr = 0
[4431]	1431
[4351]	1432	CASE ( 'wq' )
	1433	IF ( av == 0 ) THEN
	1434	to_be_resorted => wq
	1435	ELSE
[4583]	1436	IF ( .NOT. ALLOCATED( wq_av ) ) THEN
[4518]	1437	ALLOCATE( wq_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4351]	1438	wq_av = REAL( fill_value, KIND = wp )
	1439	ENDIF
	1440	to_be_resorted => wq_av
	1441	ENDIF
	1442	flag_nr = 0
	1443
[4866]	1444	CASE ( 'ws' )
	1445	IF ( av == 0 ) THEN
	1446	to_be_resorted => ws
	1447	ELSE
	1448	IF ( .NOT. ALLOCATED( ws_av ) ) THEN
	1449	ALLOCATE( ws_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
	1450	ws_av = REAL( fill_value, KIND = wp )
	1451	ENDIF
	1452	to_be_resorted => ws_av
	1453	ENDIF
	1454	flag_nr = 0
	1455
[4431]	1456	CASE ( 'wspeed' )
	1457	IF ( av == 0 ) THEN
	1458	to_be_resorted => wspeed
	1459	ELSE
[4583]	1460	IF ( .NOT. ALLOCATED( wspeed_av ) ) THEN
[4518]	1461	ALLOCATE( wspeed_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4431]	1462	wspeed_av = REAL( fill_value, KIND = wp )
	1463	ENDIF
	1464	to_be_resorted => wspeed_av
	1465	ENDIF
	1466	flag_nr = 0
	1467
	1468	CASE ( 'wdir' )
	1469	IF ( av == 0 ) THEN
	1470	to_be_resorted => wdir
	1471	ELSE
[4583]	1472	IF ( .NOT. ALLOCATED( wdir_av ) ) THEN
[4518]	1473	ALLOCATE( wdir_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4431]	1474	wdir_av = REAL( fill_value, KIND = wp )
	1475	ENDIF
	1476	to_be_resorted => wdir_av
	1477	ENDIF
	1478	flag_nr = 0
	1479
[4039]	1480	CASE DEFAULT
	1481	found = .FALSE.
	1482
	1483	END SELECT
[4431]	1484
	1485	IF ( found .AND. .NOT. resorted ) THEN
[4039]	1486	DO i = nxl, nxr
	1487	DO j = nys, nyn
	1488	DO k = nzb_do, nzt_do
[4583]	1489	local_pf(i,j,k) = MERGE( to_be_resorted(k,j,i), &
	1490	REAL( fill_value, KIND = wp ), &
	1491	BTEST( wall_flags_total_0(k,j,i), flag_nr ) )
[4039]	1492	ENDDO
	1493	ENDDO
	1494	ENDDO
	1495	ENDIF
	1496
	1497	END SUBROUTINE doq_output_3d
[4431]	1498
[4583]	1499
	1500	!--------------------------------------------------------------------------------------------------!
	1501	!
[3994]	1502	! Description:
	1503	! ------------
[4583]	1504	!> Resorts the user-defined output quantity with indices (k,j,i) to a temporary array with indices
	1505	!> (i,j,k) for masked data output.
	1506	!--------------------------------------------------------------------------------------------------!
[4069]	1507	SUBROUTINE doq_output_mask( av, variable, found, local_pf, mid )
[4431]	1508
[4039]	1509	USE control_parameters
[4431]	1510
[4039]	1511	USE indices
[4167]	1512
[4039]	1513	IMPLICIT NONE
[3994]	1514
[4583]	1515	REAL(wp), PARAMETER :: fill_value = -9999.0_wp !< value for the _FillValue attribute
	1516
[4167]	1517	CHARACTER (LEN=*) :: variable !<
	1518	CHARACTER (LEN=5) :: grid !< flag to distinquish between staggered grids
[3994]	1519
[4167]	1520	INTEGER(iwp) :: av !< index indicating averaged or instantaneous output
	1521	INTEGER(iwp) :: flag_nr !< number of the topography flag (0: scalar, 1: u, 2: v, 3: w)
	1522	INTEGER(iwp) :: i !< index variable along x-direction
	1523	INTEGER(iwp) :: j !< index variable along y-direction
	1524	INTEGER(iwp) :: k !< index variable along z-direction
	1525	INTEGER(iwp) :: im !< loop index for masked variables
	1526	INTEGER(iwp) :: jm !< loop index for masked variables
	1527	INTEGER(iwp) :: kk !< masked output index variable along z-direction
	1528	INTEGER(iwp) :: mid !< masked output running index
	1529	INTEGER(iwp) :: ktt !< k index of highest horizontal surface
[3994]	1530
[4167]	1531	LOGICAL :: found !< true if variable is in list
	1532	LOGICAL :: resorted !< true if array is resorted
[3994]	1533
[4583]	1534	REAL(wp), DIMENSION(mask_size_l(mid,1),mask_size_l(mid,2),mask_size_l(mid,3)) :: local_pf !<
[3994]	1535
[4583]	1536	REAL(wp), DIMENSION(:,:,:), POINTER :: to_be_resorted !< points to array which needs to be resorted for output
[4167]	1537
[4583]	1538
[4039]	1539	flag_nr = 0
	1540	found = .TRUE.
	1541	resorted = .FALSE.
	1542	grid = 's'
	1543
	1544	SELECT CASE ( TRIM( variable ) )
	1545
[4757]	1546	CASE ( 'rh' )
	1547	IF ( av == 0 ) THEN
	1548	to_be_resorted => rh
	1549	ELSE
	1550	to_be_resorted => rh_av
	1551	ENDIF
	1552	grid = 's'
	1553	flag_nr = 0
	1554
[4861]	1555	CASE ( 'ta' )
	1556	IF ( av == 0 ) THEN
	1557	to_be_resorted => ta
	1558	ELSE
	1559	to_be_resorted => ta_av
	1560	ENDIF
	1561	grid = 's'
	1562	flag_nr = 0
	1563
[4039]	1564	CASE ( 'ti' )
	1565	IF ( av == 0 ) THEN
	1566	to_be_resorted => ti
	1567	ELSE
	1568	to_be_resorted => ti_av
	1569	ENDIF
	1570	grid = 's'
	1571	flag_nr = 0
[4431]	1572
[4039]	1573	CASE ( 'uu' )
	1574	IF ( av == 0 ) THEN
	1575	to_be_resorted => uu
	1576	ELSE
	1577	to_be_resorted => uu_av
	1578	ENDIF
	1579	grid = 'u'
	1580	flag_nr = 1
[4431]	1581
[4039]	1582	CASE ( 'vv' )
	1583	IF ( av == 0 ) THEN
	1584	to_be_resorted => vv
	1585	ELSE
	1586	to_be_resorted => vv_av
	1587	ENDIF
	1588	grid = 'v'
	1589	flag_nr = 2
[4431]	1590
[4039]	1591	CASE ( 'ww' )
	1592	IF ( av == 0 ) THEN
	1593	to_be_resorted => ww
	1594	ELSE
	1595	to_be_resorted => ww_av
	1596	ENDIF
	1597	grid = 'w'
	1598	flag_nr = 3
[4431]	1599
[4351]	1600	CASE ( 'wu' )
	1601	IF ( av == 0 ) THEN
	1602	to_be_resorted => wu
	1603	ELSE
	1604	to_be_resorted => wu_av
	1605	ENDIF
	1606	grid = 's'
	1607	flag_nr = 0
[4431]	1608
[4351]	1609	CASE ( 'wv' )
	1610	IF ( av == 0 ) THEN
	1611	to_be_resorted => wv
	1612	ELSE
	1613	to_be_resorted => wv_av
	1614	ENDIF
	1615	grid = 's'
	1616	flag_nr = 0
[4431]	1617
[4351]	1618	CASE ( 'wtheta' )
	1619	IF ( av == 0 ) THEN
	1620	to_be_resorted => wtheta
	1621	ELSE
	1622	to_be_resorted => wtheta_av
	1623	ENDIF
	1624	grid = 's'
	1625	flag_nr = 0
[4431]	1626
[4351]	1627	CASE ( 'wq' )
	1628	IF ( av == 0 ) THEN
	1629	to_be_resorted => wq
	1630	ELSE
	1631	to_be_resorted => wq_av
	1632	ENDIF
	1633	grid = 's'
	1634	flag_nr = 0
[4039]	1635
[4866]	1636	CASE ( 'ws' )
	1637	IF ( av == 0 ) THEN
	1638	to_be_resorted => ws
	1639	ELSE
	1640	to_be_resorted => ws_av
	1641	ENDIF
	1642	grid = 's'
	1643	flag_nr = 0
	1644
[4431]	1645	CASE ( 'wspeed' )
	1646	IF ( av == 0 ) THEN
	1647	to_be_resorted => wspeed
	1648	ELSE
	1649	to_be_resorted => wspeed_av
	1650	ENDIF
	1651	grid = 's'
	1652	flag_nr = 0
	1653
	1654	CASE ( 'wdir' )
	1655	IF ( av == 0 ) THEN
	1656	to_be_resorted => wdir
	1657	ELSE
	1658	to_be_resorted => wdir_av
	1659	ENDIF
	1660	grid = 's'
	1661	flag_nr = 0
	1662
[4039]	1663	CASE DEFAULT
	1664	found = .FALSE.
	1665
	1666	END SELECT
[4431]	1667
[4039]	1668	IF ( found .AND. .NOT. resorted ) THEN
	1669	IF ( .NOT. mask_surface(mid) ) THEN
	1670	!
	1671	!-- Default masked output
	1672	DO i = 1, mask_size_l(mid,1)
	1673	DO j = 1, mask_size_l(mid,2)
	1674	DO k = 1, mask_size_l(mid,3)
[4583]	1675	local_pf(i,j,k) = MERGE( to_be_resorted(mask_k(mid,k), &
	1676	mask_j(mid,j), &
	1677	mask_i(mid,i)), &
	1678	REAL( fill_value, KIND = wp ), &
	1679	BTEST( wall_flags_total_0(mask_k(mid,k), &
	1680	mask_j(mid,j), &
	1681	mask_i(mid,i)), &
[4431]	1682	flag_nr ) )
[4039]	1683	ENDDO
	1684	ENDDO
	1685	ENDDO
	1686
	1687	ELSE
	1688	!
	1689	!-- Terrain-following masked output
	1690	DO i = 1, mask_size_l(mid,1)
	1691	DO j = 1, mask_size_l(mid,2)
	1692	!
[4167]	1693	!-- Get k index of the highest terraing surface
	1694	im = mask_i(mid,i)
	1695	jm = mask_j(mid,j)
[4583]	1696	ktt = MINLOC( MERGE( 1, 0, BTEST( wall_flags_total_0(:,jm,im), 5 ) ), DIM=1 ) - 1
[4167]	1697	DO k = 1, mask_size_l(mid,3)
	1698	kk = MIN( ktt+mask_k(mid,k), nzt+1 )
[4039]	1699	!
[4167]	1700	!-- Set value if not in building
[4346]	1701	IF ( BTEST( wall_flags_total_0(kk,jm,im), 6 ) ) THEN
[4167]	1702	local_pf(i,j,k) = fill_value
	1703	ELSE
	1704	local_pf(i,j,k) = to_be_resorted(kk,jm,im)
	1705	ENDIF
[4039]	1706	ENDDO
	1707	ENDDO
	1708	ENDDO
	1709
	1710	ENDIF
	1711	ENDIF
[4431]	1712
[4039]	1713	END SUBROUTINE doq_output_mask
	1714
[4583]	1715	!--------------------------------------------------------------------------------------------------!
[4039]	1716	! Description:
	1717	! ------------
	1718	!> Allocate required arrays
[4583]	1719	!--------------------------------------------------------------------------------------------------!
[4039]	1720	SUBROUTINE doq_init
	1721
[3994]	1722	IMPLICIT NONE
[4431]	1723
[4039]	1724	INTEGER(iwp) :: ivar !< loop index over all 2d/3d/mask output quantities
[4157]	1725
[4583]	1726
[4039]	1727	!
	1728	!-- Next line is to avoid compiler warnings about unused variables
	1729	IF ( timestep_number_at_prev_calc == 0 ) CONTINUE
[4157]	1730	!
	1731	!-- Preparatory steps and initialization of output arrays
	1732	IF ( .NOT. prepared_diagnostic_output_quantities ) CALL doq_prepare
[3994]	1733
[4039]	1734	initialized_diagnostic_output_quantities = .FALSE.
[4431]	1735
[4039]	1736	ivar = 1
	1737
[4431]	1738	DO WHILE ( ivar <= SIZE( do_all ) )
	1739
[4039]	1740	SELECT CASE ( TRIM( do_all(ivar) ) )
	1741	!
[4757]	1742	!-- Allocate array for 'relative humidity'
	1743	CASE ( 'rh' )
	1744	IF ( .NOT. ALLOCATED( rh ) ) THEN
	1745	ALLOCATE( rh(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
	1746	rh = 0.0_wp
	1747	ENDIF
	1748	!
[4861]	1749	!-- Allocate array for 'air temperature'
	1750	CASE ( 'ta' )
	1751	IF ( .NOT. ALLOCATED( ta ) ) THEN
	1752	ALLOCATE( ta(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
	1753	ta = 0.0_wp
	1754	ENDIF
	1755	!
[4039]	1756	!-- Allocate array for 'turbulence intensity'
	1757	CASE ( 'ti' )
	1758	IF ( .NOT. ALLOCATED( ti ) ) THEN
[4518]	1759	ALLOCATE( ti(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4039]	1760	ti = 0.0_wp
	1761	ENDIF
	1762	!
	1763	!-- Allocate array for uu
	1764	CASE ( 'uu' )
	1765	IF ( .NOT. ALLOCATED( uu ) ) THEN
[4518]	1766	ALLOCATE( uu(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4039]	1767	uu = 0.0_wp
	1768	ENDIF
	1769	!
	1770	!-- Allocate array for vv
	1771	CASE ( 'vv' )
	1772	IF ( .NOT. ALLOCATED( vv ) ) THEN
[4518]	1773	ALLOCATE( vv(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4039]	1774	vv = 0.0_wp
	1775	ENDIF
	1776	!
	1777	!-- Allocate array for ww
	1778	CASE ( 'ww' )
	1779	IF ( .NOT. ALLOCATED( ww ) ) THEN
[4518]	1780	ALLOCATE( ww(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4039]	1781	ww = 0.0_wp
	1782	ENDIF
[4331]	1783	!
[4351]	1784	!-- Allocate array for wu
	1785	CASE ( 'wu' )
	1786	IF ( .NOT. ALLOCATED( wu ) ) THEN
[4518]	1787	ALLOCATE( wu(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4351]	1788	wu = 0.0_wp
	1789	ENDIF
	1790	!
	1791	!-- Allocate array for wv
	1792	CASE ( 'wv' )
	1793	IF ( .NOT. ALLOCATED( wv ) ) THEN
[4518]	1794	ALLOCATE( wv(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4351]	1795	wv = 0.0_wp
	1796	ENDIF
	1797	!
	1798	!-- Allocate array for wtheta
	1799	CASE ( 'wtheta' )
	1800	IF ( .NOT. ALLOCATED( wtheta ) ) THEN
[4518]	1801	ALLOCATE( wtheta(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4351]	1802	wtheta = 0.0_wp
	1803	ENDIF
	1804	!
	1805	!-- Allocate array for wq
	1806	CASE ( 'wq' )
	1807	IF ( .NOT. ALLOCATED( wq ) ) THEN
[4518]	1808	ALLOCATE( wq(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4351]	1809	wq = 0.0_wp
	1810	ENDIF
	1811	!
[4866]	1812	!-- Allocate array for ws
	1813	CASE ( 'ws' )
	1814	IF ( .NOT. ALLOCATED( ws ) ) THEN
	1815	ALLOCATE( ws(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
	1816	ws = 0.0_wp
	1817	ENDIF
	1818	!
[4331]	1819	!-- Allocate array for 2-m potential temperature
	1820	CASE ( 'theta_2m*' )
	1821	IF ( .NOT. ALLOCATED( pt_2m ) ) THEN
[4518]	1822	ALLOCATE( pt_2m(nysg:nyng,nxlg:nxrg) )
[4331]	1823	pt_2m = 0.0_wp
	1824	ENDIF
	1825	!
	1826	!-- Allocate array for 10-m wind speed
	1827	CASE ( 'wspeed_10m*' )
	1828	IF ( .NOT. ALLOCATED( uv_10m ) ) THEN
[4518]	1829	ALLOCATE( uv_10m(nysg:nyng,nxlg:nxrg) )
[4331]	1830	uv_10m = 0.0_wp
	1831	ENDIF
[4431]	1832	!
	1833	!-- Allocate array for wspeed
	1834	CASE ( 'wspeed' )
	1835	IF ( .NOT. ALLOCATED( wspeed ) ) THEN
[4518]	1836	ALLOCATE( wspeed(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4431]	1837	wspeed = 0.0_wp
	1838	ENDIF
[4039]	1839
[4431]	1840	!
	1841	!-- Allocate array for wdir
	1842	CASE ( 'wdir' )
	1843	IF ( .NOT. ALLOCATED( u_center ) ) THEN
[4518]	1844	ALLOCATE( u_center(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4431]	1845	u_center = 0.0_wp
	1846	ENDIF
	1847	IF ( .NOT. ALLOCATED( v_center ) ) THEN
[4518]	1848	ALLOCATE( v_center(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4431]	1849	v_center = 0.0_wp
	1850	ENDIF
	1851	IF ( .NOT. ALLOCATED( wdir ) ) THEN
[4518]	1852	ALLOCATE( wdir(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4431]	1853	wdir = 0.0_wp
	1854	ENDIF
	1855
[4039]	1856	END SELECT
	1857
	1858	ivar = ivar + 1
	1859	ENDDO
	1860
	1861	initialized_diagnostic_output_quantities = .TRUE.
[4431]	1862
[4039]	1863	END SUBROUTINE doq_init
	1864
	1865
	1866	!--------------------------------------------------------------------------------------------------!
	1867	! Description:
	1868	! ------------
[4757]	1869	!> Calculation of module-specific statistics, i.e. horizontally averaged profiles and time series.
	1870	!> This is called for every statistic region sr, but at least for the region "total domain" (sr=0).
	1871	!--------------------------------------------------------------------------------------------------!
	1872	SUBROUTINE doq_statistics( mode, sr, tn )
	1873
	1874
	1875	CHARACTER (LEN=*) :: mode !<
	1876
	1877	! INTEGER(iwp) :: i !<
	1878	! INTEGER(iwp) :: j !<
	1879	! INTEGER(iwp) :: k !<
	1880	INTEGER(iwp) :: sr !<
	1881	INTEGER(iwp) :: tn !<
	1882	!
	1883	!-- Next line is to avoid compiler warning about unused variables. Please remove.
	1884	IF ( sr == 0 .OR. tn == 0 ) CONTINUE
	1885
	1886	IF ( mode == 'profiles' ) THEN
	1887
	1888	ELSEIF ( mode == 'time_series' ) THEN
	1889
	1890	ENDIF
	1891
	1892	END SUBROUTINE doq_statistics
	1893
	1894
	1895	!--------------------------------------------------------------------------------------------------!
	1896	! Description:
	1897	! ------------
[4039]	1898	!> Calculate diagnostic quantities
	1899	!--------------------------------------------------------------------------------------------------!
	1900	SUBROUTINE doq_calculate
	1901
	1902	IMPLICIT NONE
	1903
[4331]	1904	INTEGER(iwp) :: i !< grid index x-dimension
[4431]	1905	INTEGER(iwp) :: j !< grid index y-dimension
[4331]	1906	INTEGER(iwp) :: k !< grid index z-dimension
[4039]	1907	INTEGER(iwp) :: ivar !< loop index over all 2d/3d/mask output quantities
[4431]	1908
[4757]	1909	REAL(wp) :: temp !< temperature
	1910	REAL(wp) :: e_s !< saturation vapor pressure
	1911	REAL(wp) :: q_s !< saturation mixing ratio
	1912
[4331]	1913	TYPE(surf_type), POINTER :: surf !< surf-type array, used to generalize subroutines
[3994]	1914
	1915
	1916	! CALL cpu_log( log_point(41), 'calculate_quantities', 'start' )
[4157]	1917
[3994]	1918	!
[4583]	1919	!-- Save timestep number to check in time_integration if doq_calculate has been called already,
	1920	!-- since the CALL occurs at two locations, but the calculations need to be done only once per
	1921	!-- timestep.
[3994]	1922	timestep_number_at_prev_calc = current_timestep_number
	1923
[4039]	1924	ivar = 1
[3994]	1925
[4431]	1926	DO WHILE ( ivar <= SIZE( do_all ) )
[3994]	1927
[4039]	1928	SELECT CASE ( TRIM( do_all(ivar) ) )
[3994]	1929	!
[4757]	1930	!-- rh (relative humidity)
	1931	CASE ( 'rh' )
	1932	DO i = nxl, nxr
	1933	DO j = nys, nyn
	1934	DO k = nzb+1, nzt
	1935	IF ( humidity ) THEN
	1936	temp = exner(k) * pt(k,j,i)
	1937	e_s = magnus( temp )
	1938	q_s = rd_d_rv * e_s / ( hyp(k) - e_s )
	1939	rh(k,j,i) = q(k,j,i) / q_s * 100.0_wp &
	1940	* MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i), 1 ) )
	1941	ENDIF
	1942	ENDDO
	1943	ENDDO
	1944	ENDDO
	1945	!
[4861]	1946	!-- ta (air temperature)
	1947	CASE ( 'ta' )
	1948	DO i = nxl, nxr
	1949	DO j = nys, nyn
	1950	DO k = nzb+1, nzt
	1951	ta(k,j,i) = exner(k) * pt(k,j,i) - degc_to_k
	1952	ENDDO
	1953	ENDDO
	1954	ENDDO
	1955	!
[3994]	1956	!-- Calculate 'turbulence intensity' from rot[(u,v,w)] at scalar grid point
	1957	CASE ( 'ti' )
	1958	DO i = nxl, nxr
	1959	DO j = nys, nyn
	1960	DO k = nzb+1, nzt
[4583]	1961	ti(k,j,i) = 0.25_wp * SQRT( &
	1962	( ( w(k,j+1,i) + w(k-1,j+1,i) &
	1963	- w(k,j-1,i) - w(k-1,j-1,i) ) * ddy &
	1964	- ( v(k+1,j,i) + v(k+1,j+1,i) &
	1965	- v(k-1,j,i) - v(k-1,j+1,i) ) * ddzu(k) )**2 &
	1966	+ ( ( u(k+1,j,i) + u(k+1,j,i+1) &
	1967	- u(k-1,j,i) - u(k-1,j,i+1) ) * ddzu(k) &
	1968	- ( w(k,j,i+1) + w(k-1,j,i+1) &
	1969	- w(k,j,i-1) - w(k-1,j,i-1) ) * ddx )**2 &
	1970	+ ( ( v(k,j,i+1) + v(k,j+1,i+1) &
	1971	- v(k,j,i-1) - v(k,j+1,i-1) ) * ddx &
	1972	- ( u(k,j+1,i) + u(k,j+1,i+1) &
	1973	- u(k,j-1,i) - u(k,j-1,i+1) ) * ddy )**2 ) &
	1974	* MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i), 0 ) )
[3994]	1975	ENDDO
	1976	ENDDO
[4431]	1977	ENDDO
[4039]	1978	!
	1979	!-- uu
	1980	CASE ( 'uu' )
	1981	DO i = nxl, nxr
	1982	DO j = nys, nyn
	1983	DO k = nzb+1, nzt
[4583]	1984	uu(k,j,i) = u(k,j,i) * u(k,j,i) &
	1985	* MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i), 1 ) )
[4039]	1986	ENDDO
	1987	ENDDO
[3994]	1988	ENDDO
[4039]	1989	!
	1990	!-- vv
	1991	CASE ( 'vv' )
	1992	DO i = nxl, nxr
	1993	DO j = nys, nyn
	1994	DO k = nzb+1, nzt
[4583]	1995	vv(k,j,i) = v(k,j,i) * v(k,j,i) &
	1996	* MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i), 2 ) )
[4039]	1997	ENDDO
	1998	ENDDO
	1999	ENDDO
	2000	!
	2001	!-- ww
	2002	CASE ( 'ww' )
	2003	DO i = nxl, nxr
	2004	DO j = nys, nyn
	2005	DO k = nzb+1, nzt-1
[4583]	2006	ww(k,j,i) = w(k,j,i) * w(k,j,i) &
	2007	* MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i), 3 ) )
[4039]	2008	ENDDO
	2009	ENDDO
	2010	ENDDO
[4331]	2011	!
[4351]	2012	!-- wu
	2013	CASE ( 'wu' )
	2014	DO i = nxl, nxr
	2015	DO j = nys, nyn
	2016	DO k = nzb+1, nzt-1
[4583]	2017	wu(k,j,i) = w(k,j,i) &
	2018	* 0.25_wp * ( u(k,j,i) + u(k,j,i+1) + u(k+1,j,i) + u(k+1,j,i+1) )&
	2019	* MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i), 0 ) )
[4351]	2020	ENDDO
	2021	ENDDO
	2022	ENDDO
	2023	!
	2024	!-- wv
	2025	CASE ( 'wv' )
	2026	DO i = nxl, nxr
	2027	DO j = nys, nyn
	2028	DO k = nzb+1, nzt-1
[4583]	2029	wv(k,j,i) = w(k,j,i) &
	2030	* 0.25_wp * ( v(k,j,i) + v(k,j+1,i) + v(k+1,j,i) + v(k+1,j+1,i) )&
	2031	* MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i), 0 ) )
[4351]	2032	ENDDO
	2033	ENDDO
	2034	ENDDO
	2035	!
	2036	!-- wtheta
	2037	CASE ( 'wtheta' )
	2038	DO i = nxl, nxr
	2039	DO j = nys, nyn
	2040	DO k = nzb+1, nzt-1
[4583]	2041	wtheta(k,j,i) = w(k,j,i) &
	2042	* 0.5_wp * ( pt(k,j,i) + pt(k+1,j,i) ) &
	2043	* MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i), 3 ))
[4351]	2044	ENDDO
	2045	ENDDO
	2046	ENDDO
	2047	!
	2048	!-- wq
	2049	CASE ( 'wq' )
	2050	DO i = nxl, nxr
	2051	DO j = nys, nyn
	2052	DO k = nzb+1, nzt-1
[4583]	2053	wq(k,j,i) = w(k,j,i) * 0.5_wp * ( q(k,j,i) + q(k+1,j,i) ) &
	2054	* MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i), 3 ) )
[4351]	2055	ENDDO
	2056	ENDDO
	2057	ENDDO
	2058	!
[4866]	2059	!-- ws
	2060	CASE ( 'ws' )
	2061	DO i = nxl, nxr
	2062	DO j = nys, nyn
	2063	DO k = nzb+1, nzt-1
	2064	ws(k,j,i) = w(k,j,i) * 0.5_wp * ( s(k,j,i) + s(k+1,j,i) ) &
	2065	* MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i), 3 ) )
	2066	ENDDO
	2067	ENDDO
	2068	ENDDO
	2069	!
[4331]	2070	!-- 2-m potential temperature
	2071	CASE ( 'theta_2m*' )
	2072	!
[4583]	2073	!-- 2-m potential temperature is caluclated from surface arrays. In case the 2m level is
	2074	!-- below the first grid point, MOST is applied, else, linear interpolation between two
	2075	!-- vertical grid levels is applied. To access all surfaces, iterate over all horizontally-
[4331]	2076	!-- upward facing surface types.
	2077	surf => surf_def_h(0)
	2078	CALL calc_pt_2m
[4671]	2079	surf => surf_lsm_h(0)
[4331]	2080	CALL calc_pt_2m
[4671]	2081	surf => surf_usm_h(0)
[4331]	2082	CALL calc_pt_2m
	2083	!
	2084	!-- 10-m wind speed
	2085	CASE ( 'wspeed_10m*' )
	2086	!
[4583]	2087	!-- 10-m wind speed is caluclated from surface arrays. In case the 10m level is below the
	2088	!-- first grid point, MOST is applied, else, linear interpolation between two vertical grid
	2089	!-- levels is applied. To access all surfaces, iterate over all horizontally-upward facing
	2090	!-- surface types.
[4331]	2091	surf => surf_def_h(0)
	2092	CALL calc_wind_10m
[4671]	2093	surf => surf_lsm_h(0)
[4331]	2094	CALL calc_wind_10m
[4671]	2095	surf => surf_usm_h(0)
[4331]	2096	CALL calc_wind_10m
[4431]	2097	!
[4583]	2098	!-- Horizontal wind speed
[4431]	2099	CASE ( 'wspeed' )
	2100	DO i = nxl, nxr
	2101	DO j = nys, nyn
	2102	DO k = nzb, nzt+1
[4583]	2103	wspeed(k,j,i) = SQRT( ( 0.5_wp * ( u(k,j,i) + u(k,j,i+1) ) )**2 &
	2104	+ ( 0.5_wp * ( v(k,j,i) + v(k,j+1,i) ) )**2 ) &
	2105	* MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i), 0 ))
[4431]	2106	ENDDO
	2107	ENDDO
	2108	ENDDO
[3994]	2109
[4431]	2110	!
[4583]	2111	!-- Horizontal wind direction
[4431]	2112	CASE ( 'wdir' )
	2113	DO i = nxl, nxr
	2114	DO j = nys, nyn
	2115	DO k = nzb, nzt+1
	2116	u_center(k,j,i) = 0.5_wp * ( u(k,j,i) + u(k,j,i+1) )
	2117	v_center(k,j,i) = 0.5_wp * ( v(k,j,i) + v(k,j+1,i) )
	2118
[4583]	2119	wdir(k,j,i) = ATAN2( u_center(k,j,i), v_center(k,j,i) ) &
	2120	/ pi * 180.0_wp + 180.0_wp
[4431]	2121	ENDDO
	2122	ENDDO
	2123	ENDDO
	2124
[4039]	2125	END SELECT
[3994]	2126
[4039]	2127	ivar = ivar + 1
[3994]	2128	ENDDO
	2129
	2130	! CALL cpu_log( log_point(41), 'calculate_quantities', 'stop' )
	2131
[4331]	2132	!
[4583]	2133	!-- The following block contains subroutines to calculate diagnostic surface quantities.
[4331]	2134	CONTAINS
[4583]	2135	!--------------------------------------------------------------------------------------------------!
[4331]	2136	! Description:
	2137	! ------------
	2138	!> Calculation of 2-m potential temperature.
[4583]	2139	!--------------------------------------------------------------------------------------------------!
[4331]	2140	SUBROUTINE calc_pt_2m
	2141
[4583]	2142	USE surface_layer_fluxes_mod, &
[4331]	2143	ONLY: psi_h
	2144
	2145	IMPLICIT NONE
	2146
	2147	INTEGER(iwp) :: kk !< running index along the z-dimension
	2148	INTEGER(iwp) :: m !< running index for surface elements
[4431]	2149
[4331]	2150	DO m = 1, surf%ns
	2151
	2152	i = surf%i(m)
	2153	j = surf%j(m)
	2154	k = surf%k(m)
	2155	!
[4583]	2156	!-- If 2-m level is below the first grid level, MOST is used for calculation of
	2157	!-- 2-m temperature.
[4331]	2158	IF ( surf%z_mo(m) > 2.0_wp ) THEN
[4583]	2159	pt_2m(j,i) = surf%pt_surface(m) + surf%ts(m) / kappa &
	2160	* ( LOG( 2.0_wp / surf%z0h(m) ) &
	2161	- psi_h( 2.0_wp / surf%ol(m) ) &
	2162	+ psi_h( surf%z0h(m) / surf%ol(m) ) )
[4331]	2163	!
[4583]	2164	!-- If 2-m level is above the first grid level, 2-m temperature is linearly interpolated
	2165	!-- between the two nearest vertical grid levels. Note, since 2-m temperature is only
	2166	!-- computed for horizontal upward-facing surfaces, only a vertical interpolation is
	2167	!-- necessary.
[4331]	2168	ELSE
	2169	!
[4431]	2170	!-- zw(k-1) defines the height of the surface.
[4331]	2171	kk = k
	2172	DO WHILE ( zu(kk) - zw(k-1) < 2.0_wp .AND. kk <= nzt )
[4431]	2173	kk = kk + 1
[4331]	2174	ENDDO
	2175	!
[4431]	2176	!-- kk defines the index of the first grid level >= 2m.
[4583]	2177	pt_2m(j,i) = pt(kk-1,j,i) + &
	2178	( zw(k-1) + 2.0_wp - zu(kk-1) ) * &
	2179	( pt(kk,j,i) - pt(kk-1,j,i) ) / &
[4331]	2180	( zu(kk) - zu(kk-1) )
	2181	ENDIF
	2182
	2183	ENDDO
	2184
	2185	END SUBROUTINE calc_pt_2m
[4431]	2186
[4583]	2187	!--------------------------------------------------------------------------------------------------!
[4331]	2188	! Description:
	2189	! ------------
	2190	!> Calculation of 10-m wind speed.
[4583]	2191	!--------------------------------------------------------------------------------------------------!
[4331]	2192	SUBROUTINE calc_wind_10m
	2193
[4583]	2194	USE surface_layer_fluxes_mod, &
[4331]	2195	ONLY: psi_m
	2196
	2197	IMPLICIT NONE
	2198
	2199	INTEGER(iwp) :: kk !< running index along the z-dimension
	2200	INTEGER(iwp) :: m !< running index for surface elements
	2201
	2202	REAL(wp) :: uv_l !< wind speed at lower grid point
	2203	REAL(wp) :: uv_u !< wind speed at upper grid point
[4431]	2204
[4583]	2205
[4331]	2206	DO m = 1, surf%ns
	2207
	2208	i = surf%i(m)
	2209	j = surf%j(m)
	2210	k = surf%k(m)
	2211	!
[4583]	2212	!-- If 10-m level is below the first grid level, MOST is used for calculation of 10-m
	2213	!-- temperature.
[4331]	2214	IF ( surf%z_mo(m) > 10.0_wp ) THEN
[4583]	2215	uv_10m(j,i) = surf%us(m) / kappa &
	2216	* ( LOG( 10.0_wp / surf%z0(m) ) &
	2217	- psi_m( 10.0_wp / surf%ol(m) ) &
	2218	+ psi_m( surf%z0(m) / surf%ol(m) ) )
[4331]	2219	!
[4583]	2220	!-- If 10-m level is above the first grid level, 10-m wind speed is linearly interpolated
	2221	!-- between the two nearest vertical grid levels. Note, since 10-m temperature is only
	2222	!-- computed for horizontal upward-facing surfaces, only a vertical interpolation is
	2223	!-- necessary.
[4331]	2224	ELSE
	2225	!
[4431]	2226	!-- zw(k-1) defines the height of the surface.
[4331]	2227	kk = k
	2228	DO WHILE ( zu(kk) - zw(k-1) < 10.0_wp .AND. kk <= nzt )
[4431]	2229	kk = kk + 1
[4331]	2230	ENDDO
	2231	!
	2232	!-- kk defines the index of the first grid level >= 10m.
[4583]	2233	uv_l = SQRT( ( 0.5_wp * ( u(kk-1,j,i) + u(kk-1,j,i+1) ) )**2 &
[4331]	2234	+ ( 0.5_wp * ( v(kk-1,j,i) + v(kk-1,j+1,i) ) )**2 )
	2235
[4583]	2236	uv_u = SQRT( ( 0.5_wp * ( u(kk,j,i) + u(kk,j,i+1) ) )**2 &
[4331]	2237	+ ( 0.5_wp * ( v(kk,j,i) + v(kk,j+1,i) ) )**2 )
	2238
[4583]	2239	uv_10m(j,i) = uv_l + ( zw(k-1) + 10.0_wp - zu(kk-1) ) * &
	2240	( uv_u - uv_l ) / &
[4331]	2241	( zu(kk) - zu(kk-1) )
	2242
	2243	ENDIF
	2244
	2245	ENDDO
	2246
	2247	END SUBROUTINE calc_wind_10m
	2248
[4039]	2249	END SUBROUTINE doq_calculate
[3994]	2250
	2251
[4583]	2252	!--------------------------------------------------------------------------------------------------!
[3994]	2253	! Description:
	2254	! ------------
[4583]	2255	!> Preparation of the diagnostic output, counting of the module-specific output quantities and
	2256	!> gathering of the output names.
	2257	!--------------------------------------------------------------------------------------------------!
[4039]	2258	SUBROUTINE doq_prepare
[3994]	2259
[4583]	2260	USE control_parameters, &
[4069]	2261	ONLY: do2d, do3d, domask, masks
[3994]	2262
	2263	IMPLICIT NONE
	2264
[4583]	2265	CHARACTER (LEN=varnamelength), DIMENSION(0:1,500) :: do2d_var = ' ' !< label array for 2d output quantities
[3994]	2266
	2267	INTEGER(iwp) :: av !< index defining type of output, av=0 instantaneous, av=1 averaged
	2268	INTEGER(iwp) :: ivar !< loop index
	2269	INTEGER(iwp) :: ivar_all !< loop index
	2270	INTEGER(iwp) :: l !< index for cutting string
[4431]	2271	INTEGER(iwp) :: mid !< masked output running index
[3994]	2272
[4583]	2273
[3994]	2274	prepared_diagnostic_output_quantities = .FALSE.
	2275
	2276	ivar = 1
	2277	ivar_all = 1
	2278
	2279	DO av = 0, 1
	2280	!
	2281	!-- Remove _xy, _xz, or _yz from string
	2282	l = MAX( 3, LEN_TRIM( do2d(av,ivar) ) )
[3998]	2283	do2d_var(av,ivar)(1:l-3) = do2d(av,ivar)(1:l-3)
[3994]	2284	!
[4039]	2285	!-- Gather 2d output quantity names.
[4583]	2286	!-- Check for double occurrence of output quantity, e.g. by _xy, _yz, _xz.
[3994]	2287	DO WHILE ( do2d_var(av,ivar)(1:1) /= ' ' )
	2288	IF ( .NOT. ANY( do_all == do2d_var(av,ivar) ) ) THEN
	2289	do_all(ivar_all) = do2d_var(av,ivar)
	2290	ENDIF
	2291	ivar = ivar + 1
	2292	ivar_all = ivar_all + 1
	2293	l = MAX( 3, LEN_TRIM( do2d(av,ivar) ) )
[3998]	2294	do2d_var(av,ivar)(1:l-3) = do2d(av,ivar)(1:l-3)
[3994]	2295	ENDDO
	2296
	2297	ivar = 1
	2298	!
[4039]	2299	!-- Gather 3d output quantity names
[3994]	2300	DO WHILE ( do3d(av,ivar)(1:1) /= ' ' )
	2301	do_all(ivar_all) = do3d(av,ivar)
	2302	ivar = ivar + 1
	2303	ivar_all = ivar_all + 1
	2304	ENDDO
	2305
	2306	ivar = 1
	2307	!
[4583]	2308	!-- Gather masked output quantity names. Also check for double output e.g. by different masks.
[3994]	2309	DO mid = 1, masks
	2310	DO WHILE ( domask(mid,av,ivar)(1:1) /= ' ' )
[4039]	2311	IF ( .NOT. ANY( do_all == domask(mid,av,ivar) ) ) THEN
[4132]	2312	do_all(ivar_all) = domask(mid,av,ivar)
[4039]	2313	ENDIF
[3994]	2314
	2315	ivar = ivar + 1
	2316	ivar_all = ivar_all + 1
	2317	ENDDO
	2318	ivar = 1
	2319	ENDDO
	2320
	2321	ENDDO
	2322
	2323	prepared_diagnostic_output_quantities = .TRUE.
	2324
[4039]	2325	END SUBROUTINE doq_prepare
[4431]	2326
[4583]	2327	!--------------------------------------------------------------------------------------------------!
[4039]	2328	! Description:
	2329	! ------------
	2330	!> Subroutine reads local (subdomain) restart data
[4583]	2331	!--------------------------------------------------------------------------------------------------!
	2332	SUBROUTINE doq_rrd_local_ftn( k, nxlf, nxlc, nxl_on_file, nxrf, nxrc, nxr_on_file, nynf, nync, &
	2333	nyn_on_file, nysf, nysc, nys_on_file, tmp_2d, tmp_3d, found )
[4431]	2334
[4518]	2335	USE control_parameters
	2336
	2337	USE indices
	2338
	2339	USE kinds
	2340
	2341	USE pegrid
	2342
	2343
	2344	IMPLICIT NONE
	2345
	2346	INTEGER(iwp) :: k !<
	2347	INTEGER(iwp) :: nxlc !<
	2348	INTEGER(iwp) :: nxlf !<
	2349	INTEGER(iwp) :: nxl_on_file !<
	2350	INTEGER(iwp) :: nxrc !<
	2351	INTEGER(iwp) :: nxrf !<
	2352	INTEGER(iwp) :: nxr_on_file !<
	2353	INTEGER(iwp) :: nync !<
	2354	INTEGER(iwp) :: nynf !<
	2355	INTEGER(iwp) :: nyn_on_file !<
	2356	INTEGER(iwp) :: nysc !<
	2357	INTEGER(iwp) :: nysf !<
	2358	INTEGER(iwp) :: nys_on_file !<
	2359
	2360	LOGICAL, INTENT(OUT) :: found
	2361
	2362	REAL(wp), DIMENSION(nys_on_file-nbgp:nyn_on_file+nbgp,nxl_on_file-nbgp:nxr_on_file+nbgp) :: tmp_2d !<
	2363
	2364	REAL(wp), DIMENSION(nzb:nzt+1,nys_on_file-nbgp:nyn_on_file+nbgp,nxl_on_file-nbgp:nxr_on_file+nbgp) :: tmp_3d !<
	2365
	2366
	2367	found = .TRUE.
	2368
	2369	SELECT CASE ( restart_string(1:length) )
	2370
	2371	CASE ( 'pt_2m_av' )
	2372	IF ( .NOT. ALLOCATED( pt_2m_av ) ) THEN
	2373	ALLOCATE( pt_2m_av(nysg:nyng,nxlg:nxrg) )
	2374	ENDIF
	2375	IF ( k == 1 ) READ ( 13 ) tmp_2d
	2376	pt_2m_av(nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = &
	2377	tmp_2d(nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp)
[4583]	2378
[4757]	2379	CASE ( 'rh_av' )
	2380	IF ( .NOT. ALLOCATED( rh_av ) ) THEN
	2381	ALLOCATE( rh_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
	2382	ENDIF
	2383	IF ( k == 1 ) READ ( 13 ) tmp_3d
	2384	rh_av(:,nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = &
	2385	tmp_3d(:,nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp)
	2386
[4861]	2387	CASE ( 'ta_av' )
	2388	IF ( .NOT. ALLOCATED( ta_av ) ) THEN
	2389	ALLOCATE( ta_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
	2390	ENDIF
	2391	IF ( k == 1 ) READ ( 13 ) tmp_3d
	2392	ta_av(:,nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = &
	2393	tmp_3d(:,nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp)
	2394
[4518]	2395	CASE ( 'ti_av' )
	2396	IF ( .NOT. ALLOCATED( ti_av ) ) THEN
	2397	ALLOCATE( ti_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
	2398	ENDIF
	2399	IF ( k == 1 ) READ ( 13 ) tmp_3d
	2400	ti_av(:,nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = &
	2401	tmp_3d(:,nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp)
	2402
	2403	CASE ( 'u_center_av' )
	2404	IF ( .NOT. ALLOCATED( u_center_av ) ) THEN
	2405	ALLOCATE( u_center_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
	2406	ENDIF
	2407	IF ( k == 1 ) READ ( 13 ) tmp_3d
	2408	u_center_av(:,nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = &
	2409	tmp_3d(:,nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp)
	2410
	2411	CASE ( 'uu_av' )
	2412	IF ( .NOT. ALLOCATED( uu_av ) ) THEN
	2413	ALLOCATE( uu_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
	2414	ENDIF
	2415	IF ( k == 1 ) READ ( 13 ) tmp_3d
	2416	uu_av(:,nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = &
	2417	tmp_3d(:,nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp)
	2418
	2419	CASE ( 'uv_10m_av' )
	2420	IF ( .NOT. ALLOCATED( uv_10m_av ) ) THEN
	2421	ALLOCATE( uv_10m_av(nysg:nyng,nxlg:nxrg) )
	2422	ENDIF
	2423	IF ( k == 1 ) READ ( 13 ) tmp_2d
	2424	uv_10m_av(nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = &
	2425	tmp_2d(nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp)
	2426
	2427	CASE ( 'v_center_av' )
	2428	IF ( .NOT. ALLOCATED( v_center_av ) ) THEN
	2429	ALLOCATE( v_center_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
	2430	ENDIF
	2431	IF ( k == 1 ) READ ( 13 ) tmp_3d
	2432	v_center_av(:,nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = &
	2433	tmp_3d(:,nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp)
	2434
	2435	CASE ( 'vv_av' )
	2436	IF ( .NOT. ALLOCATED( vv_av ) ) THEN
	2437	ALLOCATE( vv_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
	2438	ENDIF
	2439	IF ( k == 1 ) READ ( 13 ) tmp_3d
	2440	vv_av(:,nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = &
	2441	tmp_3d(:,nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp)
	2442
	2443	CASE ( 'wtheta_av' )
	2444	IF ( .NOT. ALLOCATED( wtheta_av ) ) THEN
	2445	ALLOCATE( wtheta_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
	2446	ENDIF
	2447	IF ( k == 1 ) READ ( 13 ) tmp_3d
	2448	wtheta_av(:,nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = &
	2449	tmp_3d(:,nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp)
	2450
	2451	CASE ( 'wq_av' )
	2452	IF ( .NOT. ALLOCATED( wq_av ) ) THEN
	2453	ALLOCATE( wq_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
	2454	ENDIF
	2455	IF ( k == 1 ) READ ( 13 ) tmp_3d
	2456	wq_av(:,nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = &
	2457	tmp_3d(:,nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp)
	2458
[4866]	2459	CASE ( 'ws_av' )
	2460	IF ( .NOT. ALLOCATED( ws_av ) ) THEN
	2461	ALLOCATE( ws_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
	2462	ENDIF
	2463	IF ( k == 1 ) READ ( 13 ) tmp_3d
	2464	ws_av(:,nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = &
	2465	tmp_3d(:,nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp)
	2466
[4518]	2467	CASE ( 'wspeed_av' )
	2468	IF ( .NOT. ALLOCATED( wspeed_av ) ) THEN
	2469	ALLOCATE( wspeed_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
	2470	ENDIF
	2471	IF ( k == 1 ) READ ( 13 ) tmp_3d
	2472	wspeed_av(:,nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = &
	2473	tmp_3d(:,nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp)
	2474
	2475	CASE ( 'wu_av' )
	2476	IF ( .NOT. ALLOCATED( wu_av ) ) THEN
	2477	ALLOCATE( wu_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
	2478	ENDIF
	2479	IF ( k == 1 ) READ ( 13 ) tmp_3d
	2480	wu_av(:,nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = &
	2481	tmp_3d(:,nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp)
	2482
	2483	CASE ( 'wv_av' )
	2484	IF ( .NOT. ALLOCATED( wv_av ) ) THEN
	2485	ALLOCATE( wv_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
	2486	ENDIF
	2487	IF ( k == 1 ) READ ( 13 ) tmp_3d
	2488	wv_av(:,nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = &
	2489	tmp_3d(:,nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp)
	2490
	2491	CASE ( 'ww_av' )
	2492	IF ( .NOT. ALLOCATED( ww_av ) ) THEN
	2493	ALLOCATE( ww_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
	2494	ENDIF
	2495	IF ( k == 1 ) READ ( 13 ) tmp_3d
	2496	ww_av(:,nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = &
	2497	tmp_3d(:,nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp)
	2498
	2499
	2500	CASE DEFAULT
	2501
	2502	found = .FALSE.
	2503
	2504	END SELECT
	2505
	2506	END SUBROUTINE doq_rrd_local_ftn
	2507
[4583]	2508	!--------------------------------------------------------------------------------------------------!
[4039]	2509	! Description:
	2510	! ------------
[4518]	2511	!> Read module-specific local restart data arrays (MPI-IO).
[4583]	2512	!--------------------------------------------------------------------------------------------------!
[4518]	2513	SUBROUTINE doq_rrd_local_mpi
	2514
	2515	LOGICAL :: array_found !< control flad indicating whether the array is found in the file or not
	2516
	2517
	2518	CALL rd_mpi_io_check_array( 'pt_2m_av' , found = array_found )
	2519	IF ( array_found ) THEN
	2520	IF ( .NOT. ALLOCATED( pt_2m_av ) ) ALLOCATE( pt_2m_av(nysg:nyng,nxlg:nxrg) )
	2521	CALL rrd_mpi_io( 'pt_2m_av', pt_2m_av )
	2522	ENDIF
	2523
[4861]	2524	CALL rd_mpi_io_check_array( 'ta_av' , found = array_found )
	2525	IF ( array_found ) THEN
	2526	IF ( .NOT. ALLOCATED( ta_av ) ) ALLOCATE( ta_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
	2527	CALL rrd_mpi_io( 'ta_av', ta_av )
	2528	ENDIF
	2529
[4518]	2530	CALL rd_mpi_io_check_array( 'ti_av' , found = array_found )
	2531	IF ( array_found ) THEN
	2532	IF ( .NOT. ALLOCATED( ti_av ) ) ALLOCATE( ti_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
	2533	CALL rrd_mpi_io( 'ti_av', ti_av )
	2534	ENDIF
	2535
[4757]	2536	CALL rd_mpi_io_check_array( 'rh_av' , found = array_found )
	2537	IF ( array_found ) THEN
	2538	IF ( .NOT. ALLOCATED( rh_av ) ) ALLOCATE( ti_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
	2539	CALL rrd_mpi_io( 'rh_av', rh_av )
	2540	ENDIF
	2541
[4518]	2542	CALL rd_mpi_io_check_array( 'u_center_av' , found = array_found )
	2543	IF ( array_found ) THEN
	2544	IF ( .NOT. ALLOCATED( u_center_av ) ) ALLOCATE( u_center_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
	2545	CALL rrd_mpi_io( 'u_center_av', u_center_av )
	2546	ENDIF
	2547
	2548	CALL rd_mpi_io_check_array( 'uu_av' , found = array_found )
	2549	IF ( array_found ) THEN
	2550	IF ( .NOT. ALLOCATED( uu_av ) ) ALLOCATE( uu_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
	2551	CALL rrd_mpi_io( 'uu_av', uu_av )
	2552	ENDIF
	2553
	2554	CALL rd_mpi_io_check_array( 'uv_10m_av' , found = array_found )
	2555	IF ( array_found ) THEN
	2556	IF ( .NOT. ALLOCATED( uv_10m_av ) ) ALLOCATE( uv_10m_av(nysg:nyng,nxlg:nxrg) )
	2557	CALL rrd_mpi_io( 'uv_10m_av', uv_10m_av )
	2558	ENDIF
	2559
	2560	CALL rd_mpi_io_check_array( 'v_center_av' , found = array_found )
	2561	IF ( array_found ) THEN
	2562	IF ( .NOT. ALLOCATED( v_center_av ) ) ALLOCATE( v_center_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
	2563	CALL rrd_mpi_io( 'v_center_av', v_center_av )
	2564	ENDIF
	2565
	2566	CALL rd_mpi_io_check_array( 'vv_av' , found = array_found )
	2567	IF ( array_found ) THEN
	2568	IF ( .NOT. ALLOCATED( vv_av ) ) ALLOCATE( vv_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
	2569	CALL rrd_mpi_io( 'vv_av', vv_av )
	2570	ENDIF
	2571
	2572	CALL rd_mpi_io_check_array( 'ww_av' , found = array_found )
	2573	IF ( array_found ) THEN
	2574	IF ( .NOT. ALLOCATED( ww_av ) ) ALLOCATE( ww_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
	2575	CALL rrd_mpi_io( 'ww_av', ww_av )
	2576	ENDIF
	2577
	2578	CALL rd_mpi_io_check_array( 'wu_av' , found = array_found )
	2579	IF ( array_found ) THEN
	2580	IF ( .NOT. ALLOCATED( wu_av ) ) ALLOCATE( wu_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
	2581	CALL rrd_mpi_io( 'wu_av', wu_av )
	2582	ENDIF
	2583
	2584	CALL rd_mpi_io_check_array( 'wv_av' , found = array_found )
	2585	IF ( array_found ) THEN
	2586	IF ( .NOT. ALLOCATED( wv_av ) ) ALLOCATE( wv_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
	2587	CALL rrd_mpi_io( 'wv_av', wv_av )
	2588	ENDIF
	2589
	2590	CALL rd_mpi_io_check_array( 'wtheta_av' , found = array_found )
	2591	IF ( array_found ) THEN
	2592	IF ( .NOT. ALLOCATED( wtheta_av ) ) ALLOCATE( wtheta_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
	2593	CALL rrd_mpi_io( 'wtheta_av', wtheta_av )
	2594	ENDIF
	2595
	2596	CALL rd_mpi_io_check_array( 'wq_av' , found = array_found )
	2597	IF ( array_found ) THEN
	2598	IF ( .NOT. ALLOCATED( wq_av ) ) ALLOCATE( wq_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
	2599	CALL rrd_mpi_io( 'wq_av', wq_av )
	2600	ENDIF
	2601
[4866]	2602	CALL rd_mpi_io_check_array( 'ws_av' , found = array_found )
	2603	IF ( array_found ) THEN
	2604	IF ( .NOT. ALLOCATED( ws_av ) ) ALLOCATE( ws_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
	2605	CALL rrd_mpi_io( 'ws_av', ws_av )
	2606	ENDIF
	2607
[4518]	2608	CALL rd_mpi_io_check_array( 'wspeed_av' , found = array_found )
	2609	IF ( array_found ) THEN
	2610	IF ( .NOT. ALLOCATED( wspeed_av ) ) ALLOCATE( wspeed_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
	2611	CALL rrd_mpi_io( 'wspeed_av', wspeed_av )
	2612	ENDIF
	2613
	2614	END SUBROUTINE doq_rrd_local_mpi
	2615
[4583]	2616	!--------------------------------------------------------------------------------------------------!
[4518]	2617	! Description:
	2618	! ------------
[4039]	2619	!> Subroutine writes local (subdomain) restart data
[4583]	2620	!--------------------------------------------------------------------------------------------------!
[4039]	2621	SUBROUTINE doq_wrd_local
[3994]	2622
[4039]	2623	IMPLICIT NONE
	2624
[4583]	2625
[4517]	2626	IF ( TRIM( restart_data_format_output ) == 'fortran_binary' ) THEN
[4331]	2627
[4517]	2628	IF ( ALLOCATED( pt_2m_av ) ) THEN
	2629	CALL wrd_write_string( 'pt_2m_av' )
	2630	WRITE ( 14 ) pt_2m_av
	2631	ENDIF
[4331]	2632
[4757]	2633	IF ( ALLOCATED( rh_av ) ) THEN
	2634	CALL wrd_write_string( 'rh_av' )
	2635	WRITE ( 14 ) rh_av
	2636	ENDIF
	2637
[4861]	2638	IF ( ALLOCATED( ta_av ) ) THEN
	2639	CALL wrd_write_string( 'ta_av' )
	2640	WRITE ( 14 ) ta_av
	2641	ENDIF
	2642
[4517]	2643	IF ( ALLOCATED( ti_av ) ) THEN
	2644	CALL wrd_write_string( 'ti_av' )
	2645	WRITE ( 14 ) ti_av
	2646	ENDIF
[4331]	2647
[4518]	2648	IF ( ALLOCATED( u_center_av ) ) THEN
	2649	CALL wrd_write_string( 'u_center_av' )
	2650	WRITE ( 14 ) u_center_av
	2651	ENDIF
	2652
[4517]	2653	IF ( ALLOCATED( uu_av ) ) THEN
	2654	CALL wrd_write_string( 'uu_av' )
	2655	WRITE ( 14 ) uu_av
	2656	ENDIF
[4331]	2657
[4517]	2658	IF ( ALLOCATED( uv_10m_av ) ) THEN
	2659	CALL wrd_write_string( 'uv_10m_av' )
	2660	WRITE ( 14 ) uv_10m_av
	2661	ENDIF
[4331]	2662
[4518]	2663	IF ( ALLOCATED( v_center_av ) ) THEN
	2664	CALL wrd_write_string( 'v_center_av' )
	2665	WRITE ( 14 ) v_center_av
	2666	ENDIF
	2667
[4517]	2668	IF ( ALLOCATED( vv_av ) ) THEN
	2669	CALL wrd_write_string( 'vv_av' )
	2670	WRITE ( 14 ) vv_av
	2671	ENDIF
[4039]	2672
[4517]	2673	IF ( ALLOCATED( ww_av ) ) THEN
	2674	CALL wrd_write_string( 'ww_av' )
	2675	WRITE ( 14 ) ww_av
	2676	ENDIF
[4039]	2677
[4517]	2678	IF ( ALLOCATED( wu_av ) ) THEN
	2679	CALL wrd_write_string( 'wu_av' )
	2680	WRITE ( 14 ) wu_av
	2681	ENDIF
[4431]	2682
[4517]	2683	IF ( ALLOCATED( wv_av ) ) THEN
	2684	CALL wrd_write_string( 'wv_av' )
	2685	WRITE ( 14 ) wv_av
	2686	ENDIF
[4431]	2687
[4517]	2688	IF ( ALLOCATED( wtheta_av ) ) THEN
	2689	CALL wrd_write_string( 'wtheta_av' )
	2690	WRITE ( 14 ) wtheta_av
	2691	ENDIF
[4351]	2692
[4517]	2693	IF ( ALLOCATED( wq_av ) ) THEN
	2694	CALL wrd_write_string( 'wq_av' )
	2695	WRITE ( 14 ) wq_av
	2696	ENDIF
[4431]	2697
[4866]	2698	IF ( ALLOCATED( ws_av ) ) THEN
	2699	CALL wrd_write_string( 'ws_av' )
	2700	WRITE ( 14 ) ws_av
	2701	ENDIF
	2702
[4517]	2703	IF ( ALLOCATED( wspeed_av ) ) THEN
	2704	CALL wrd_write_string( 'wspeed_av' )
	2705	WRITE ( 14 ) wspeed_av
	2706	ENDIF
[4431]	2707
[4535]	2708	ELSEIF ( restart_data_format_output(1:3) == 'mpi' ) THEN
[4517]	2709
[4518]	2710	IF ( ALLOCATED( pt_2m_av ) ) CALL wrd_mpi_io( 'pt_2m_av', pt_2m_av )
[4757]	2711	IF ( ALLOCATED( rh_av ) ) CALL wrd_mpi_io( 'rh_av', rh_av )
[4861]	2712	IF ( ALLOCATED( ta_av ) ) CALL wrd_mpi_io( 'ta_av', ta_av )
[4518]	2713	IF ( ALLOCATED( ti_av ) ) CALL wrd_mpi_io( 'ti_av', ti_av )
	2714	IF ( ALLOCATED( u_center_av ) ) CALL wrd_mpi_io( 'u_center_av', u_center_av )
	2715	IF ( ALLOCATED( uu_av ) ) CALL wrd_mpi_io( 'uu_av', uu_av )
	2716	IF ( ALLOCATED( uv_10m_av ) ) CALL wrd_mpi_io( 'uv_10m_av', uv_10m_av )
	2717	IF ( ALLOCATED( vv_av ) ) CALL wrd_mpi_io( 'vv_av', vv_av )
	2718	IF ( ALLOCATED( v_center_av ) ) CALL wrd_mpi_io( 'v_center_av', v_center_av )
	2719	IF ( ALLOCATED( wtheta_av ) ) CALL wrd_mpi_io( 'wtheta_av', wtheta_av )
	2720	IF ( ALLOCATED( wq_av ) ) CALL wrd_mpi_io( 'wq_av', wq_av )
[4866]	2721	IF ( ALLOCATED( ws_av ) ) CALL wrd_mpi_io( 'ws_av', ws_av )
[4518]	2722	IF ( ALLOCATED( wspeed_av ) ) CALL wrd_mpi_io( 'wspeed_av', wspeed_av )
	2723	IF ( ALLOCATED( wu_av ) ) CALL wrd_mpi_io( 'wu_av', wu_av )
	2724	IF ( ALLOCATED( wv_av ) ) CALL wrd_mpi_io( 'wv_av', wv_av )
	2725	IF ( ALLOCATED( ww_av ) ) CALL wrd_mpi_io( 'ww_av', ww_av )
[4517]	2726
[4431]	2727	ENDIF
	2728
[4039]	2729	END SUBROUTINE doq_wrd_local
	2730
[4431]	2731
[3994]	2732	END MODULE diagnostic_output_quantities_mod

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |