Home

Context Navigation

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

Last change on this file since 4716 was 4671, checked in by pavelkrc, 4 years ago
Radiative transfer model RTM version 4.1
Property svn:keywords set to `Id`
File size: 81.2 KB

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