Home

Context Navigation

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

Last change on this file since 4526 was 4518, checked in by suehring, 5 years ago
Diagnostic output: Define arrays over ghost points in order to allow for standard mpi-io treatment. By this modularization of restart-data input is possible with the module interface. Move input of restart data to doq_rrd_local. Enable mpi-io for restart data. Bugfix: add missing restart input of wtheta_av, wq_av, wu_av, and wv_av.
Property svn:keywords set to `Id`
File size: 79.3 KB

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