Home

Context Navigation

source: palm/trunk/SOURCE/data_output_2d.f90 @ 4337

Last change on this file since 4337 was 4331, checked in by suehring, 5 years ago
New diagnostic output for 10-m wind speed; Diagnostic output of 2-m potential temperature moved to diagnostic output
Property svn:keywords set to `Id`
File size: 88.4 KB

Line
1	!> @file data_output_2d.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	!
17	! Copyright 1997-2019 Leibniz Universitaet Hannover
18	!------------------------------------------------------------------------------!
19	!
20	! Current revisions:
21	! ------------------
22	!
23	!
24	! Former revisions:
25	! -----------------
26	! $Id: data_output_2d.f90 4331 2019-12-10 18:25:02Z Giersch $
27	! Move 2-m potential temperature output to diagnostic_output_quantities
28	!
29	! 4329 2019-12-10 15:46:36Z motisi
30	! Renamed wall_flags_0 to wall_flags_static_0
31	!
32	! 4182 2019-08-22 15:20:23Z scharf
33	! Corrected "Former revisions" section
34	!
35	! 4048 2019-06-21 21:00:21Z knoop
36	! Removed turbulence_closure_mod dependency
37	!
38	! 4039 2019-06-18 10:32:41Z suehring
39	! modularize diagnostic output
40	!
41	! 3994 2019-05-22 18:08:09Z suehring
42	! output of turbulence intensity added
43	!
44	! 3987 2019-05-22 09:52:13Z kanani
45	! Introduce alternative switch for debug output during timestepping
46	!
47	! 3943 2019-05-02 09:50:41Z maronga
48	! Added output of qsws for green roofs.
49	!
50	! 3885 2019-04-11 11:29:34Z kanani
51	! Changes related to global restructuring of location messages and introduction
52	! of additional debug messages
53	!
54	! 3766 2019-02-26 16:23:41Z raasch
55	! unused variables removed
56	!
57	! 3655 2019-01-07 16:51:22Z knoop
58	! Bugfix: use time_since_reference_point instead of simulated_time (relevant
59	! when using wall/soil spinup)
60	!
61	! Revision 1.1 1997/08/11 06:24:09 raasch
62	! Initial revision
63	!
64	!
65	! Description:
66	! ------------
67	!> Data output of cross-sections in netCDF format or binary format
68	!> to be later converted to NetCDF by helper routine combine_plot_fields.
69	!> Attention: The position of the sectional planes is still not always computed
70	!> --------- correctly. (zu is used always)!
71	!------------------------------------------------------------------------------!
72	SUBROUTINE data_output_2d( mode, av )
73
74
75	USE arrays_3d, &
76	ONLY: dzw, d_exner, e, heatflux_output_conversion, p, pt, q, ql, ql_c, ql_v, s, tend, u, &
77	v, vpt, w, waterflux_output_conversion, zu, zw
78
79	USE averaging
80
81	USE basic_constants_and_equations_mod, &
82	ONLY: c_p, lv_d_cp, l_v
83
84	USE bulk_cloud_model_mod, &
85	ONLY: bulk_cloud_model
86
87	USE control_parameters, &
88	ONLY: data_output_2d_on_each_pe, &
89	data_output_xy, data_output_xz, data_output_yz, &
90	debug_output_timestep, &
91	do2d, &
92	do2d_xy_last_time, do2d_xy_time_count, &
93	do2d_xz_last_time, do2d_xz_time_count, &
94	do2d_yz_last_time, do2d_yz_time_count, &
95	ibc_uv_b, io_blocks, io_group, message_string, &
96	ntdim_2d_xy, ntdim_2d_xz, ntdim_2d_yz, &
97	psolver, section, &
98	time_since_reference_point
99
100	USE cpulog, &
101	ONLY: cpu_log, log_point
102
103	USE indices, &
104	ONLY: nbgp, nx, nxl, nxlg, nxr, nxrg, ny, nyn, nyng, nys, nysg, &
105	nzb, nzt, wall_flags_static_0
106
107	USE kinds
108
109	USE land_surface_model_mod, &
110	ONLY: zs
111
112	USE module_interface, &
113	ONLY: module_interface_data_output_2d
114
115	#if defined( __netcdf )
116	USE NETCDF
117	#endif
118
119	USE netcdf_interface, &
120	ONLY: fill_value, id_set_xy, id_set_xz, id_set_yz, id_var_do2d, &
121	id_var_time_xy, id_var_time_xz, id_var_time_yz, nc_stat, &
122	netcdf_data_format, netcdf_handle_error
123
124	USE particle_attributes, &
125	ONLY: grid_particles, number_of_particles, particle_advection_start, &
126	particles, prt_count
127
128	USE pegrid
129
130	USE surface_mod, &
131	ONLY: ind_pav_green, ind_veg_wall, ind_wat_win, surf_def_h, &
132	surf_lsm_h, surf_usm_h
133
134
135	IMPLICIT NONE
136
137	CHARACTER (LEN=2) :: do2d_mode !< output mode of variable ('xy', 'xz', 'yz')
138	CHARACTER (LEN=2) :: mode !< mode with which the routine is called ('xy', 'xz', 'yz')
139	CHARACTER (LEN=4) :: grid !< string defining the vertical grid
140
141	INTEGER(iwp) :: av !< flag for (non-)average output
142	INTEGER(iwp) :: ngp !< number of grid points of an output slice
143	INTEGER(iwp) :: file_id !< id of output files
144	INTEGER(iwp) :: flag_nr !< number of masking flag
145	INTEGER(iwp) :: i !< loop index
146	INTEGER(iwp) :: iis !< vertical index of a xy slice in array 'local_2d_sections'
147	INTEGER(iwp) :: is !< slice index
148	INTEGER(iwp) :: ivar !< variable index
149	INTEGER(iwp) :: j !< loop index
150	INTEGER(iwp) :: k !< loop index
151	INTEGER(iwp) :: l !< loop index
152	INTEGER(iwp) :: layer_xy !< vertical index of a xy slice in array 'local_pf'
153	INTEGER(iwp) :: m !< loop index
154	INTEGER(iwp) :: n !< loop index
155	INTEGER(iwp) :: nis !< number of vertical slices to be written via parallel NetCDF output
156	INTEGER(iwp) :: ns !< number of output slices
157	INTEGER(iwp) :: nzb_do !< lower limit of the data field (usually nzb)
158	INTEGER(iwp) :: nzt_do !< upper limit of the data field (usually nzt+1)
159	INTEGER(iwp) :: s_ind !< index of slice types (xy=1, xz=2, yz=3)
160	INTEGER(iwp) :: sender !< PE id of sending PE
161	INTEGER(iwp) :: ind(4) !< index limits (lower/upper bounds) of array 'local_2d'
162
163	LOGICAL :: found !< true if output variable was found
164	LOGICAL :: resorted !< true if variable is resorted
165	LOGICAL :: two_d !< true if variable is only two dimensional
166
167	REAL(wp) :: mean_r !< mean particle radius
168	REAL(wp) :: s_r2 !< sum( particle-radius**2 )
169	REAL(wp) :: s_r3 !< sum( particle-radius**3 )
170
171	REAL(wp), DIMENSION(:), ALLOCATABLE :: level_z !< z levels for output array
172	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: local_2d !< local 2-dimensional array containing output values
173	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: local_2d_l !< local 2-dimensional array containing output values
174
175	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: local_pf !< output array
176	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: local_2d_sections !< local array containing values at all slices
177	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: local_2d_sections_l !< local array containing values at all slices
178
179	#if defined( __parallel )
180	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: total_2d !< same as local_2d
181	#endif
182	REAL(wp), DIMENSION(:,:,:), POINTER :: to_be_resorted !< points to array which shall be output
183
184
185	IF ( debug_output_timestep ) CALL debug_message( 'data_output_2d', 'start' )
186	!
187	!-- Immediate return, if no output is requested (no respective sections
188	!-- found in parameter data_output)
189	IF ( mode == 'xy' .AND. .NOT. data_output_xy(av) ) RETURN
190	IF ( mode == 'xz' .AND. .NOT. data_output_xz(av) ) RETURN
191	IF ( mode == 'yz' .AND. .NOT. data_output_yz(av) ) RETURN
192
193	CALL cpu_log (log_point(3),'data_output_2d','start')
194
195	two_d = .FALSE. ! local variable to distinguish between output of pure 2D
196	! arrays and cross-sections of 3D arrays.
197
198	!
199	!-- Depending on the orientation of the cross-section, the respective output
200	!-- files have to be opened.
201	SELECT CASE ( mode )
202
203	CASE ( 'xy' )
204	s_ind = 1
205	ALLOCATE( level_z(nzb:nzt+1), local_2d(nxl:nxr,nys:nyn) )
206
207	IF ( netcdf_data_format > 4 ) THEN
208	ns = 1
209	DO WHILE ( section(ns,s_ind) /= -9999 .AND. ns <= 100 )
210	ns = ns + 1
211	ENDDO
212	ns = ns - 1
213	ALLOCATE( local_2d_sections(nxl:nxr,nys:nyn,1:ns) )
214	local_2d_sections = 0.0_wp
215	ENDIF
216
217	!
218	!-- Parallel netCDF4/HDF5 output is done on all PEs, all other on PE0 only
219	IF ( myid == 0 .OR. netcdf_data_format > 4 ) THEN
220	CALL check_open( 101+av*10 )
221	ENDIF
222	IF ( data_output_2d_on_each_pe .AND. netcdf_data_format < 5 ) THEN
223	CALL check_open( 21 )
224	ELSE
225	IF ( myid == 0 ) THEN
226	#if defined( __parallel )
227	ALLOCATE( total_2d(0:nx,0:ny) )
228	#endif
229	ENDIF
230	ENDIF
231
232	CASE ( 'xz' )
233	s_ind = 2
234	ALLOCATE( local_2d(nxl:nxr,nzb:nzt+1) )
235
236	IF ( netcdf_data_format > 4 ) THEN
237	ns = 1
238	DO WHILE ( section(ns,s_ind) /= -9999 .AND. ns <= 100 )
239	ns = ns + 1
240	ENDDO
241	ns = ns - 1
242	ALLOCATE( local_2d_sections(nxl:nxr,1:ns,nzb:nzt+1) )
243	ALLOCATE( local_2d_sections_l(nxl:nxr,1:ns,nzb:nzt+1) )
244	local_2d_sections = 0.0_wp; local_2d_sections_l = 0.0_wp
245	ENDIF
246
247	!
248	!-- Parallel netCDF4/HDF5 output is done on all PEs, all other on PE0 only
249	IF ( myid == 0 .OR. netcdf_data_format > 4 ) THEN
250	CALL check_open( 102+av*10 )
251	ENDIF
252
253	IF ( data_output_2d_on_each_pe .AND. netcdf_data_format < 5 ) THEN
254	CALL check_open( 22 )
255	ELSE
256	IF ( myid == 0 ) THEN
257	#if defined( __parallel )
258	ALLOCATE( total_2d(0:nx,nzb:nzt+1) )
259	#endif
260	ENDIF
261	ENDIF
262
263	CASE ( 'yz' )
264	s_ind = 3
265	ALLOCATE( local_2d(nys:nyn,nzb:nzt+1) )
266
267	IF ( netcdf_data_format > 4 ) THEN
268	ns = 1
269	DO WHILE ( section(ns,s_ind) /= -9999 .AND. ns <= 100 )
270	ns = ns + 1
271	ENDDO
272	ns = ns - 1
273	ALLOCATE( local_2d_sections(1:ns,nys:nyn,nzb:nzt+1) )
274	ALLOCATE( local_2d_sections_l(1:ns,nys:nyn,nzb:nzt+1) )
275	local_2d_sections = 0.0_wp; local_2d_sections_l = 0.0_wp
276	ENDIF
277
278	!
279	!-- Parallel netCDF4/HDF5 output is done on all PEs, all other on PE0 only
280	IF ( myid == 0 .OR. netcdf_data_format > 4 ) THEN
281	CALL check_open( 103+av*10 )
282	ENDIF
283
284	IF ( data_output_2d_on_each_pe .AND. netcdf_data_format < 5 ) THEN
285	CALL check_open( 23 )
286	ELSE
287	IF ( myid == 0 ) THEN
288	#if defined( __parallel )
289	ALLOCATE( total_2d(0:ny,nzb:nzt+1) )
290	#endif
291	ENDIF
292	ENDIF
293
294	CASE DEFAULT
295	message_string = 'unknown cross-section: ' // TRIM( mode )
296	CALL message( 'data_output_2d', 'PA0180', 1, 2, 0, 6, 0 )
297
298	END SELECT
299
300	!
301	!-- For parallel netcdf output the time axis must be limited. Return, if this
302	!-- limit is exceeded. This could be the case, if the simulated time exceeds
303	!-- the given end time by the length of the given output interval.
304	IF ( netcdf_data_format > 4 ) THEN
305	IF ( mode == 'xy' .AND. do2d_xy_time_count(av) + 1 > &
306	ntdim_2d_xy(av) ) THEN
307	WRITE ( message_string, * ) 'Output of xy cross-sections is not ', &
308	'given at t=', time_since_reference_point, 's because the', &
309	' maximum number of output time levels is exceeded.'
310	CALL message( 'data_output_2d', 'PA0384', 0, 1, 0, 6, 0 )
311	CALL cpu_log( log_point(3), 'data_output_2d', 'stop' )
312	RETURN
313	ENDIF
314	IF ( mode == 'xz' .AND. do2d_xz_time_count(av) + 1 > &
315	ntdim_2d_xz(av) ) THEN
316	WRITE ( message_string, * ) 'Output of xz cross-sections is not ', &
317	'given at t=', time_since_reference_point, 's because the', &
318	' maximum number of output time levels is exceeded.'
319	CALL message( 'data_output_2d', 'PA0385', 0, 1, 0, 6, 0 )
320	CALL cpu_log( log_point(3), 'data_output_2d', 'stop' )
321	RETURN
322	ENDIF
323	IF ( mode == 'yz' .AND. do2d_yz_time_count(av) + 1 > &
324	ntdim_2d_yz(av) ) THEN
325	WRITE ( message_string, * ) 'Output of yz cross-sections is not ', &
326	'given at t=', time_since_reference_point, 's because the', &
327	' maximum number of output time levels is exceeded.'
328	CALL message( 'data_output_2d', 'PA0386', 0, 1, 0, 6, 0 )
329	CALL cpu_log( log_point(3), 'data_output_2d', 'stop' )
330	RETURN
331	ENDIF
332	ENDIF
333
334	!
335	!-- Allocate a temporary array for resorting (kji -> ijk).
336	ALLOCATE( local_pf(nxl:nxr,nys:nyn,nzb:nzt+1) )
337	local_pf = 0.0
338
339	!
340	!-- Loop of all variables to be written.
341	!-- Output dimensions chosen
342	ivar = 1
343	l = MAX( 2, LEN_TRIM( do2d(av,ivar) ) )
344	do2d_mode = do2d(av,ivar)(l-1:l)
345
346	DO WHILE ( do2d(av,ivar)(1:1) /= ' ' )
347
348	IF ( do2d_mode == mode ) THEN
349	!
350	!-- Set flag to steer output of radiation, land-surface, or user-defined
351	!-- quantities
352	found = .FALSE.
353
354	nzb_do = nzb
355	nzt_do = nzt+1
356	!
357	!-- Before each output, set array local_pf to fill value
358	local_pf = fill_value
359	!
360	!-- Set masking flag for topography for not resorted arrays
361	flag_nr = 0
362
363	!
364	!-- Store the array chosen on the temporary array.
365	resorted = .FALSE.
366	SELECT CASE ( TRIM( do2d(av,ivar) ) )
367	CASE ( 'e_xy', 'e_xz', 'e_yz' )
368	IF ( av == 0 ) THEN
369	to_be_resorted => e
370	ELSE
371	IF ( .NOT. ALLOCATED( e_av ) ) THEN
372	ALLOCATE( e_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
373	e_av = REAL( fill_value, KIND = wp )
374	ENDIF
375	to_be_resorted => e_av
376	ENDIF
377	IF ( mode == 'xy' ) level_z = zu
378
379	CASE ( 'thetal_xy', 'thetal_xz', 'thetal_yz' )
380	IF ( av == 0 ) THEN
381	to_be_resorted => pt
382	ELSE
383	IF ( .NOT. ALLOCATED( lpt_av ) ) THEN
384	ALLOCATE( lpt_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
385	lpt_av = REAL( fill_value, KIND = wp )
386	ENDIF
387	to_be_resorted => lpt_av
388	ENDIF
389	IF ( mode == 'xy' ) level_z = zu
390
391	CASE ( 'lwp*_xy' ) ! 2d-array
392	IF ( av == 0 ) THEN
393	DO i = nxl, nxr
394	DO j = nys, nyn
395	local_pf(i,j,nzb+1) = SUM( ql(nzb:nzt,j,i) * &
396	dzw(1:nzt+1) )
397	ENDDO
398	ENDDO
399	ELSE
400	IF ( .NOT. ALLOCATED( lwp_av ) ) THEN
401	ALLOCATE( lwp_av(nysg:nyng,nxlg:nxrg) )
402	lwp_av = REAL( fill_value, KIND = wp )
403	ENDIF
404	DO i = nxl, nxr
405	DO j = nys, nyn
406	local_pf(i,j,nzb+1) = lwp_av(j,i)
407	ENDDO
408	ENDDO
409	ENDIF
410	resorted = .TRUE.
411	two_d = .TRUE.
412	level_z(nzb+1) = zu(nzb+1)
413
414	CASE ( 'ghf*_xy' ) ! 2d-array
415	IF ( av == 0 ) THEN
416	DO m = 1, surf_lsm_h%ns
417	i = surf_lsm_h%i(m)
418	j = surf_lsm_h%j(m)
419	local_pf(i,j,nzb+1) = surf_lsm_h%ghf(m)
420	ENDDO
421	DO m = 1, surf_usm_h%ns
422	i = surf_usm_h%i(m)
423	j = surf_usm_h%j(m)
424	local_pf(i,j,nzb+1) = surf_usm_h%frac(ind_veg_wall,m) * &
425	surf_usm_h%wghf_eb(m) + &
426	surf_usm_h%frac(ind_pav_green,m) * &
427	surf_usm_h%wghf_eb_green(m) + &
428	surf_usm_h%frac(ind_wat_win,m) * &
429	surf_usm_h%wghf_eb_window(m)
430	ENDDO
431	ELSE
432	IF ( .NOT. ALLOCATED( ghf_av ) ) THEN
433	ALLOCATE( ghf_av(nysg:nyng,nxlg:nxrg) )
434	ghf_av = REAL( fill_value, KIND = wp )
435	ENDIF
436	DO i = nxl, nxr
437	DO j = nys, nyn
438	local_pf(i,j,nzb+1) = ghf_av(j,i)
439	ENDDO
440	ENDDO
441	ENDIF
442
443	resorted = .TRUE.
444	two_d = .TRUE.
445	level_z(nzb+1) = zu(nzb+1)
446
447	CASE ( 'ol*_xy' ) ! 2d-array
448	IF ( av == 0 ) THEN
449	DO m = 1, surf_def_h(0)%ns
450	i = surf_def_h(0)%i(m)
451	j = surf_def_h(0)%j(m)
452	local_pf(i,j,nzb+1) = surf_def_h(0)%ol(m)
453	ENDDO
454	DO m = 1, surf_lsm_h%ns
455	i = surf_lsm_h%i(m)
456	j = surf_lsm_h%j(m)
457	local_pf(i,j,nzb+1) = surf_lsm_h%ol(m)
458	ENDDO
459	DO m = 1, surf_usm_h%ns
460	i = surf_usm_h%i(m)
461	j = surf_usm_h%j(m)
462	local_pf(i,j,nzb+1) = surf_usm_h%ol(m)
463	ENDDO
464	ELSE
465	IF ( .NOT. ALLOCATED( ol_av ) ) THEN
466	ALLOCATE( ol_av(nysg:nyng,nxlg:nxrg) )
467	ol_av = REAL( fill_value, KIND = wp )
468	ENDIF
469	DO i = nxl, nxr
470	DO j = nys, nyn
471	local_pf(i,j,nzb+1) = ol_av(j,i)
472	ENDDO
473	ENDDO
474	ENDIF
475	resorted = .TRUE.
476	two_d = .TRUE.
477	level_z(nzb+1) = zu(nzb+1)
478
479	CASE ( 'p_xy', 'p_xz', 'p_yz' )
480	IF ( av == 0 ) THEN
481	IF ( psolver /= 'sor' ) CALL exchange_horiz( p, nbgp )
482	to_be_resorted => p
483	ELSE
484	IF ( .NOT. ALLOCATED( p_av ) ) THEN
485	ALLOCATE( p_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
486	p_av = REAL( fill_value, KIND = wp )
487	ENDIF
488	IF ( psolver /= 'sor' ) CALL exchange_horiz( p_av, nbgp )
489	to_be_resorted => p_av
490	ENDIF
491	IF ( mode == 'xy' ) level_z = zu
492
493	CASE ( 'pc_xy', 'pc_xz', 'pc_yz' ) ! particle concentration
494	IF ( av == 0 ) THEN
495	IF ( time_since_reference_point >= particle_advection_start ) THEN
496	tend = prt_count
497	! CALL exchange_horiz( tend, nbgp )
498	ELSE
499	tend = 0.0_wp
500	ENDIF
501	DO i = nxl, nxr
502	DO j = nys, nyn
503	DO k = nzb, nzt+1
504	local_pf(i,j,k) = tend(k,j,i)
505	ENDDO
506	ENDDO
507	ENDDO
508	resorted = .TRUE.
509	ELSE
510	IF ( .NOT. ALLOCATED( pc_av ) ) THEN
511	ALLOCATE( pc_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
512	pc_av = REAL( fill_value, KIND = wp )
513	ENDIF
514	! CALL exchange_horiz( pc_av, nbgp )
515	to_be_resorted => pc_av
516	ENDIF
517
518	CASE ( 'pr_xy', 'pr_xz', 'pr_yz' ) ! mean particle radius (effective radius)
519	IF ( av == 0 ) THEN
520	IF ( time_since_reference_point >= particle_advection_start ) THEN
521	DO i = nxl, nxr
522	DO j = nys, nyn
523	DO k = nzb, nzt+1
524	number_of_particles = prt_count(k,j,i)
525	IF (number_of_particles <= 0) CYCLE
526	particles => grid_particles(k,j,i)%particles(1:number_of_particles)
527	s_r2 = 0.0_wp
528	s_r3 = 0.0_wp
529	DO n = 1, number_of_particles
530	IF ( particles(n)%particle_mask ) THEN
531	s_r2 = s_r2 + particles(n)%radius*2 &
532	particles(n)%weight_factor
533	s_r3 = s_r3 + particles(n)%radius*3 &
534	particles(n)%weight_factor
535	ENDIF
536	ENDDO
537	IF ( s_r2 > 0.0_wp ) THEN
538	mean_r = s_r3 / s_r2
539	ELSE
540	mean_r = 0.0_wp
541	ENDIF
542	tend(k,j,i) = mean_r
543	ENDDO
544	ENDDO
545	ENDDO
546	! CALL exchange_horiz( tend, nbgp )
547	ELSE
548	tend = 0.0_wp
549	ENDIF
550	DO i = nxl, nxr
551	DO j = nys, nyn
552	DO k = nzb, nzt+1
553	local_pf(i,j,k) = tend(k,j,i)
554	ENDDO
555	ENDDO
556	ENDDO
557	resorted = .TRUE.
558	ELSE
559	IF ( .NOT. ALLOCATED( pr_av ) ) THEN
560	ALLOCATE( pr_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
561	pr_av = REAL( fill_value, KIND = wp )
562	ENDIF
563	! CALL exchange_horiz( pr_av, nbgp )
564	to_be_resorted => pr_av
565	ENDIF
566
567	CASE ( 'theta_xy', 'theta_xz', 'theta_yz' )
568	IF ( av == 0 ) THEN
569	IF ( .NOT. bulk_cloud_model ) THEN
570	to_be_resorted => pt
571	ELSE
572	DO i = nxl, nxr
573	DO j = nys, nyn
574	DO k = nzb, nzt+1
575	local_pf(i,j,k) = pt(k,j,i) + lv_d_cp * &
576	d_exner(k) * &
577	ql(k,j,i)
578	ENDDO
579	ENDDO
580	ENDDO
581	resorted = .TRUE.
582	ENDIF
583	ELSE
584	IF ( .NOT. ALLOCATED( pt_av ) ) THEN
585	ALLOCATE( pt_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
586	pt_av = REAL( fill_value, KIND = wp )
587	ENDIF
588	to_be_resorted => pt_av
589	ENDIF
590	IF ( mode == 'xy' ) level_z = zu
591
592	CASE ( 'q_xy', 'q_xz', 'q_yz' )
593	IF ( av == 0 ) THEN
594	to_be_resorted => q
595	ELSE
596	IF ( .NOT. ALLOCATED( q_av ) ) THEN
597	ALLOCATE( q_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
598	q_av = REAL( fill_value, KIND = wp )
599	ENDIF
600	to_be_resorted => q_av
601	ENDIF
602	IF ( mode == 'xy' ) level_z = zu
603
604	CASE ( 'ql_xy', 'ql_xz', 'ql_yz' )
605	IF ( av == 0 ) THEN
606	to_be_resorted => ql
607	ELSE
608	IF ( .NOT. ALLOCATED( ql_av ) ) THEN
609	ALLOCATE( ql_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
610	ql_av = REAL( fill_value, KIND = wp )
611	ENDIF
612	to_be_resorted => ql_av
613	ENDIF
614	IF ( mode == 'xy' ) level_z = zu
615
616	CASE ( 'ql_c_xy', 'ql_c_xz', 'ql_c_yz' )
617	IF ( av == 0 ) THEN
618	to_be_resorted => ql_c
619	ELSE
620	IF ( .NOT. ALLOCATED( ql_c_av ) ) THEN
621	ALLOCATE( ql_c_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
622	ql_c_av = REAL( fill_value, KIND = wp )
623	ENDIF
624	to_be_resorted => ql_c_av
625	ENDIF
626	IF ( mode == 'xy' ) level_z = zu
627
628	CASE ( 'ql_v_xy', 'ql_v_xz', 'ql_v_yz' )
629	IF ( av == 0 ) THEN
630	to_be_resorted => ql_v
631	ELSE
632	IF ( .NOT. ALLOCATED( ql_v_av ) ) THEN
633	ALLOCATE( ql_v_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
634	ql_v_av = REAL( fill_value, KIND = wp )
635	ENDIF
636	to_be_resorted => ql_v_av
637	ENDIF
638	IF ( mode == 'xy' ) level_z = zu
639
640	CASE ( 'ql_vp_xy', 'ql_vp_xz', 'ql_vp_yz' )
641	IF ( av == 0 ) THEN
642	IF ( time_since_reference_point >= particle_advection_start ) THEN
643	DO i = nxl, nxr
644	DO j = nys, nyn
645	DO k = nzb, nzt+1
646	number_of_particles = prt_count(k,j,i)
647	IF (number_of_particles <= 0) CYCLE
648	particles => grid_particles(k,j,i)%particles(1:number_of_particles)
649	DO n = 1, number_of_particles
650	IF ( particles(n)%particle_mask ) THEN
651	tend(k,j,i) = tend(k,j,i) + &
652	particles(n)%weight_factor / &
653	prt_count(k,j,i)
654	ENDIF
655	ENDDO
656	ENDDO
657	ENDDO
658	ENDDO
659	! CALL exchange_horiz( tend, nbgp )
660	ELSE
661	tend = 0.0_wp
662	ENDIF
663	DO i = nxl, nxr
664	DO j = nys, nyn
665	DO k = nzb, nzt+1
666	local_pf(i,j,k) = tend(k,j,i)
667	ENDDO
668	ENDDO
669	ENDDO
670	resorted = .TRUE.
671	ELSE
672	IF ( .NOT. ALLOCATED( ql_vp_av ) ) THEN
673	ALLOCATE( ql_vp_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
674	ql_vp_av = REAL( fill_value, KIND = wp )
675	ENDIF
676	! CALL exchange_horiz( ql_vp_av, nbgp )
677	to_be_resorted => ql_vp_av
678	ENDIF
679	IF ( mode == 'xy' ) level_z = zu
680
681	CASE ( 'qsws*_xy' ) ! 2d-array
682	IF ( av == 0 ) THEN
683	local_pf(:,:,nzb+1) = REAL( fill_value, KIND = wp )
684	!
685	!-- In case of default surfaces, clean-up flux by density.
686	!-- In case of land-surfaces, convert fluxes into
687	!-- dynamic units
688	DO m = 1, surf_def_h(0)%ns
689	i = surf_def_h(0)%i(m)
690	j = surf_def_h(0)%j(m)
691	k = surf_def_h(0)%k(m)
692	local_pf(i,j,nzb+1) = surf_def_h(0)%qsws(m) * &
693	waterflux_output_conversion(k)
694	ENDDO
695	DO m = 1, surf_lsm_h%ns
696	i = surf_lsm_h%i(m)
697	j = surf_lsm_h%j(m)
698	k = surf_lsm_h%k(m)
699	local_pf(i,j,nzb+1) = surf_lsm_h%qsws(m) * l_v
700	ENDDO
701	DO m = 1, surf_usm_h%ns
702	i = surf_usm_h%i(m)
703	j = surf_usm_h%j(m)
704	k = surf_usm_h%k(m)
705	local_pf(i,j,nzb+1) = surf_usm_h%qsws(m) * l_v
706	ENDDO
707	ELSE
708	IF ( .NOT. ALLOCATED( qsws_av ) ) THEN
709	ALLOCATE( qsws_av(nysg:nyng,nxlg:nxrg) )
710	qsws_av = REAL( fill_value, KIND = wp )
711	ENDIF
712	DO i = nxl, nxr
713	DO j = nys, nyn
714	local_pf(i,j,nzb+1) = qsws_av(j,i)
715	ENDDO
716	ENDDO
717	ENDIF
718	resorted = .TRUE.
719	two_d = .TRUE.
720	level_z(nzb+1) = zu(nzb+1)
721
722	CASE ( 'qv_xy', 'qv_xz', 'qv_yz' )
723	IF ( av == 0 ) THEN
724	DO i = nxl, nxr
725	DO j = nys, nyn
726	DO k = nzb, nzt+1
727	local_pf(i,j,k) = q(k,j,i) - ql(k,j,i)
728	ENDDO
729	ENDDO
730	ENDDO
731	resorted = .TRUE.
732	ELSE
733	IF ( .NOT. ALLOCATED( qv_av ) ) THEN
734	ALLOCATE( qv_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
735	qv_av = REAL( fill_value, KIND = wp )
736	ENDIF
737	to_be_resorted => qv_av
738	ENDIF
739	IF ( mode == 'xy' ) level_z = zu
740
741	CASE ( 'r_a*_xy' ) ! 2d-array
742	IF ( av == 0 ) THEN
743	DO m = 1, surf_lsm_h%ns
744	i = surf_lsm_h%i(m)
745	j = surf_lsm_h%j(m)
746	local_pf(i,j,nzb+1) = surf_lsm_h%r_a(m)
747	ENDDO
748
749	DO m = 1, surf_usm_h%ns
750	i = surf_usm_h%i(m)
751	j = surf_usm_h%j(m)
752	local_pf(i,j,nzb+1) = &
753	( surf_usm_h%frac(ind_veg_wall,m) * &
754	surf_usm_h%r_a(m) + &
755	surf_usm_h%frac(ind_pav_green,m) * &
756	surf_usm_h%r_a_green(m) + &
757	surf_usm_h%frac(ind_wat_win,m) * &
758	surf_usm_h%r_a_window(m) )
759	ENDDO
760	ELSE
761	IF ( .NOT. ALLOCATED( r_a_av ) ) THEN
762	ALLOCATE( r_a_av(nysg:nyng,nxlg:nxrg) )
763	r_a_av = REAL( fill_value, KIND = wp )
764	ENDIF
765	DO i = nxl, nxr
766	DO j = nys, nyn
767	local_pf(i,j,nzb+1) = r_a_av(j,i)
768	ENDDO
769	ENDDO
770	ENDIF
771	resorted = .TRUE.
772	two_d = .TRUE.
773	level_z(nzb+1) = zu(nzb+1)
774
775	CASE ( 's_xy', 's_xz', 's_yz' )
776	IF ( av == 0 ) THEN
777	to_be_resorted => s
778	ELSE
779	IF ( .NOT. ALLOCATED( s_av ) ) THEN
780	ALLOCATE( s_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
781	s_av = REAL( fill_value, KIND = wp )
782	ENDIF
783	to_be_resorted => s_av
784	ENDIF
785
786	CASE ( 'shf*_xy' ) ! 2d-array
787	IF ( av == 0 ) THEN
788	!
789	!-- In case of default surfaces, clean-up flux by density.
790	!-- In case of land- and urban-surfaces, convert fluxes into
791	!-- dynamic units.
792	DO m = 1, surf_def_h(0)%ns
793	i = surf_def_h(0)%i(m)
794	j = surf_def_h(0)%j(m)
795	k = surf_def_h(0)%k(m)
796	local_pf(i,j,nzb+1) = surf_def_h(0)%shf(m) * &
797	heatflux_output_conversion(k)
798	ENDDO
799	DO m = 1, surf_lsm_h%ns
800	i = surf_lsm_h%i(m)
801	j = surf_lsm_h%j(m)
802	k = surf_lsm_h%k(m)
803	local_pf(i,j,nzb+1) = surf_lsm_h%shf(m) * c_p
804	ENDDO
805	DO m = 1, surf_usm_h%ns
806	i = surf_usm_h%i(m)
807	j = surf_usm_h%j(m)
808	k = surf_usm_h%k(m)
809	local_pf(i,j,nzb+1) = surf_usm_h%shf(m) * c_p
810	ENDDO
811	ELSE
812	IF ( .NOT. ALLOCATED( shf_av ) ) THEN
813	ALLOCATE( shf_av(nysg:nyng,nxlg:nxrg) )
814	shf_av = REAL( fill_value, KIND = wp )
815	ENDIF
816	DO i = nxl, nxr
817	DO j = nys, nyn
818	local_pf(i,j,nzb+1) = shf_av(j,i)
819	ENDDO
820	ENDDO
821	ENDIF
822	resorted = .TRUE.
823	two_d = .TRUE.
824	level_z(nzb+1) = zu(nzb+1)
825
826	CASE ( 'ssws*_xy' ) ! 2d-array
827	IF ( av == 0 ) THEN
828	DO m = 1, surf_def_h(0)%ns
829	i = surf_def_h(0)%i(m)
830	j = surf_def_h(0)%j(m)
831	local_pf(i,j,nzb+1) = surf_def_h(0)%ssws(m)
832	ENDDO
833	DO m = 1, surf_lsm_h%ns
834	i = surf_lsm_h%i(m)
835	j = surf_lsm_h%j(m)
836	local_pf(i,j,nzb+1) = surf_lsm_h%ssws(m)
837	ENDDO
838	DO m = 1, surf_usm_h%ns
839	i = surf_usm_h%i(m)
840	j = surf_usm_h%j(m)
841	local_pf(i,j,nzb+1) = surf_usm_h%ssws(m)
842	ENDDO
843	ELSE
844	IF ( .NOT. ALLOCATED( ssws_av ) ) THEN
845	ALLOCATE( ssws_av(nysg:nyng,nxlg:nxrg) )
846	ssws_av = REAL( fill_value, KIND = wp )
847	ENDIF
848	DO i = nxl, nxr
849	DO j = nys, nyn
850	local_pf(i,j,nzb+1) = ssws_av(j,i)
851	ENDDO
852	ENDDO
853	ENDIF
854	resorted = .TRUE.
855	two_d = .TRUE.
856	level_z(nzb+1) = zu(nzb+1)
857
858	CASE ( 't*_xy' ) ! 2d-array
859	IF ( av == 0 ) THEN
860	DO m = 1, surf_def_h(0)%ns
861	i = surf_def_h(0)%i(m)
862	j = surf_def_h(0)%j(m)
863	local_pf(i,j,nzb+1) = surf_def_h(0)%ts(m)
864	ENDDO
865	DO m = 1, surf_lsm_h%ns
866	i = surf_lsm_h%i(m)
867	j = surf_lsm_h%j(m)
868	local_pf(i,j,nzb+1) = surf_lsm_h%ts(m)
869	ENDDO
870	DO m = 1, surf_usm_h%ns
871	i = surf_usm_h%i(m)
872	j = surf_usm_h%j(m)
873	local_pf(i,j,nzb+1) = surf_usm_h%ts(m)
874	ENDDO
875	ELSE
876	IF ( .NOT. ALLOCATED( ts_av ) ) THEN
877	ALLOCATE( ts_av(nysg:nyng,nxlg:nxrg) )
878	ts_av = REAL( fill_value, KIND = wp )
879	ENDIF
880	DO i = nxl, nxr
881	DO j = nys, nyn
882	local_pf(i,j,nzb+1) = ts_av(j,i)
883	ENDDO
884	ENDDO
885	ENDIF
886	resorted = .TRUE.
887	two_d = .TRUE.
888	level_z(nzb+1) = zu(nzb+1)
889
890	CASE ( 'tsurf*_xy' ) ! 2d-array
891	IF ( av == 0 ) THEN
892	DO m = 1, surf_def_h(0)%ns
893	i = surf_def_h(0)%i(m)
894	j = surf_def_h(0)%j(m)
895	local_pf(i,j,nzb+1) = surf_def_h(0)%pt_surface(m)
896	ENDDO
897
898	DO m = 1, surf_lsm_h%ns
899	i = surf_lsm_h%i(m)
900	j = surf_lsm_h%j(m)
901	local_pf(i,j,nzb+1) = surf_lsm_h%pt_surface(m)
902	ENDDO
903
904	DO m = 1, surf_usm_h%ns
905	i = surf_usm_h%i(m)
906	j = surf_usm_h%j(m)
907	local_pf(i,j,nzb+1) = surf_usm_h%pt_surface(m)
908	ENDDO
909
910	ELSE
911	IF ( .NOT. ALLOCATED( tsurf_av ) ) THEN
912	ALLOCATE( tsurf_av(nysg:nyng,nxlg:nxrg) )
913	tsurf_av = REAL( fill_value, KIND = wp )
914	ENDIF
915	DO i = nxl, nxr
916	DO j = nys, nyn
917	local_pf(i,j,nzb+1) = tsurf_av(j,i)
918	ENDDO
919	ENDDO
920	ENDIF
921	resorted = .TRUE.
922	two_d = .TRUE.
923	level_z(nzb+1) = zu(nzb+1)
924
925	CASE ( 'u_xy', 'u_xz', 'u_yz' )
926	flag_nr = 1
927	IF ( av == 0 ) THEN
928	to_be_resorted => u
929	ELSE
930	IF ( .NOT. ALLOCATED( u_av ) ) THEN
931	ALLOCATE( u_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
932	u_av = REAL( fill_value, KIND = wp )
933	ENDIF
934	to_be_resorted => u_av
935	ENDIF
936	IF ( mode == 'xy' ) level_z = zu
937	!
938	!-- Substitute the values generated by "mirror" boundary condition
939	!-- at the bottom boundary by the real surface values.
940	IF ( do2d(av,ivar) == 'u_xz' .OR. do2d(av,ivar) == 'u_yz' ) THEN
941	IF ( ibc_uv_b == 0 ) local_pf(:,:,nzb) = 0.0_wp
942	ENDIF
943
944	CASE ( 'us*_xy' ) ! 2d-array
945	IF ( av == 0 ) THEN
946	DO m = 1, surf_def_h(0)%ns
947	i = surf_def_h(0)%i(m)
948	j = surf_def_h(0)%j(m)
949	local_pf(i,j,nzb+1) = surf_def_h(0)%us(m)
950	ENDDO
951	DO m = 1, surf_lsm_h%ns
952	i = surf_lsm_h%i(m)
953	j = surf_lsm_h%j(m)
954	local_pf(i,j,nzb+1) = surf_lsm_h%us(m)
955	ENDDO
956	DO m = 1, surf_usm_h%ns
957	i = surf_usm_h%i(m)
958	j = surf_usm_h%j(m)
959	local_pf(i,j,nzb+1) = surf_usm_h%us(m)
960	ENDDO
961	ELSE
962	IF ( .NOT. ALLOCATED( us_av ) ) THEN
963	ALLOCATE( us_av(nysg:nyng,nxlg:nxrg) )
964	us_av = REAL( fill_value, KIND = wp )
965	ENDIF
966	DO i = nxl, nxr
967	DO j = nys, nyn
968	local_pf(i,j,nzb+1) = us_av(j,i)
969	ENDDO
970	ENDDO
971	ENDIF
972	resorted = .TRUE.
973	two_d = .TRUE.
974	level_z(nzb+1) = zu(nzb+1)
975
976	CASE ( 'v_xy', 'v_xz', 'v_yz' )
977	flag_nr = 2
978	IF ( av == 0 ) THEN
979	to_be_resorted => v
980	ELSE
981	IF ( .NOT. ALLOCATED( v_av ) ) THEN
982	ALLOCATE( v_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
983	v_av = REAL( fill_value, KIND = wp )
984	ENDIF
985	to_be_resorted => v_av
986	ENDIF
987	IF ( mode == 'xy' ) level_z = zu
988	!
989	!-- Substitute the values generated by "mirror" boundary condition
990	!-- at the bottom boundary by the real surface values.
991	IF ( do2d(av,ivar) == 'v_xz' .OR. do2d(av,ivar) == 'v_yz' ) THEN
992	IF ( ibc_uv_b == 0 ) local_pf(:,:,nzb) = 0.0_wp
993	ENDIF
994
995	CASE ( 'thetav_xy', 'thetav_xz', 'thetav_yz' )
996	IF ( av == 0 ) THEN
997	to_be_resorted => vpt
998	ELSE
999	IF ( .NOT. ALLOCATED( vpt_av ) ) THEN
1000	ALLOCATE( vpt_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
1001	vpt_av = REAL( fill_value, KIND = wp )
1002	ENDIF
1003	to_be_resorted => vpt_av
1004	ENDIF
1005	IF ( mode == 'xy' ) level_z = zu
1006
1007	CASE ( 'w_xy', 'w_xz', 'w_yz' )
1008	flag_nr = 3
1009	IF ( av == 0 ) THEN
1010	to_be_resorted => w
1011	ELSE
1012	IF ( .NOT. ALLOCATED( w_av ) ) THEN
1013	ALLOCATE( w_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
1014	w_av = REAL( fill_value, KIND = wp )
1015	ENDIF
1016	to_be_resorted => w_av
1017	ENDIF
1018	IF ( mode == 'xy' ) level_z = zw
1019
1020	CASE ( 'z0*_xy' ) ! 2d-array
1021	IF ( av == 0 ) THEN
1022	DO m = 1, surf_def_h(0)%ns
1023	i = surf_def_h(0)%i(m)
1024	j = surf_def_h(0)%j(m)
1025	local_pf(i,j,nzb+1) = surf_def_h(0)%z0(m)
1026	ENDDO
1027	DO m = 1, surf_lsm_h%ns
1028	i = surf_lsm_h%i(m)
1029	j = surf_lsm_h%j(m)
1030	local_pf(i,j,nzb+1) = surf_lsm_h%z0(m)
1031	ENDDO
1032	DO m = 1, surf_usm_h%ns
1033	i = surf_usm_h%i(m)
1034	j = surf_usm_h%j(m)
1035	local_pf(i,j,nzb+1) = surf_usm_h%z0(m)
1036	ENDDO
1037	ELSE
1038	IF ( .NOT. ALLOCATED( z0_av ) ) THEN
1039	ALLOCATE( z0_av(nysg:nyng,nxlg:nxrg) )
1040	z0_av = REAL( fill_value, KIND = wp )
1041	ENDIF
1042	DO i = nxl, nxr
1043	DO j = nys, nyn
1044	local_pf(i,j,nzb+1) = z0_av(j,i)
1045	ENDDO
1046	ENDDO
1047	ENDIF
1048	resorted = .TRUE.
1049	two_d = .TRUE.
1050	level_z(nzb+1) = zu(nzb+1)
1051
1052	CASE ( 'z0h*_xy' ) ! 2d-array
1053	IF ( av == 0 ) THEN
1054	DO m = 1, surf_def_h(0)%ns
1055	i = surf_def_h(0)%i(m)
1056	j = surf_def_h(0)%j(m)
1057	local_pf(i,j,nzb+1) = surf_def_h(0)%z0h(m)
1058	ENDDO
1059	DO m = 1, surf_lsm_h%ns
1060	i = surf_lsm_h%i(m)
1061	j = surf_lsm_h%j(m)
1062	local_pf(i,j,nzb+1) = surf_lsm_h%z0h(m)
1063	ENDDO
1064	DO m = 1, surf_usm_h%ns
1065	i = surf_usm_h%i(m)
1066	j = surf_usm_h%j(m)
1067	local_pf(i,j,nzb+1) = surf_usm_h%z0h(m)
1068	ENDDO
1069	ELSE
1070	IF ( .NOT. ALLOCATED( z0h_av ) ) THEN
1071	ALLOCATE( z0h_av(nysg:nyng,nxlg:nxrg) )
1072	z0h_av = REAL( fill_value, KIND = wp )
1073	ENDIF
1074	DO i = nxl, nxr
1075	DO j = nys, nyn
1076	local_pf(i,j,nzb+1) = z0h_av(j,i)
1077	ENDDO
1078	ENDDO
1079	ENDIF
1080	resorted = .TRUE.
1081	two_d = .TRUE.
1082	level_z(nzb+1) = zu(nzb+1)
1083
1084	CASE ( 'z0q*_xy' ) ! 2d-array
1085	IF ( av == 0 ) THEN
1086	DO m = 1, surf_def_h(0)%ns
1087	i = surf_def_h(0)%i(m)
1088	j = surf_def_h(0)%j(m)
1089	local_pf(i,j,nzb+1) = surf_def_h(0)%z0q(m)
1090	ENDDO
1091	DO m = 1, surf_lsm_h%ns
1092	i = surf_lsm_h%i(m)
1093	j = surf_lsm_h%j(m)
1094	local_pf(i,j,nzb+1) = surf_lsm_h%z0q(m)
1095	ENDDO
1096	DO m = 1, surf_usm_h%ns
1097	i = surf_usm_h%i(m)
1098	j = surf_usm_h%j(m)
1099	local_pf(i,j,nzb+1) = surf_usm_h%z0q(m)
1100	ENDDO
1101	ELSE
1102	IF ( .NOT. ALLOCATED( z0q_av ) ) THEN
1103	ALLOCATE( z0q_av(nysg:nyng,nxlg:nxrg) )
1104	z0q_av = REAL( fill_value, KIND = wp )
1105	ENDIF
1106	DO i = nxl, nxr
1107	DO j = nys, nyn
1108	local_pf(i,j,nzb+1) = z0q_av(j,i)
1109	ENDDO
1110	ENDDO
1111	ENDIF
1112	resorted = .TRUE.
1113	two_d = .TRUE.
1114	level_z(nzb+1) = zu(nzb+1)
1115
1116	CASE DEFAULT
1117
1118	!
1119	!-- Quantities of other modules
1120	IF ( .NOT. found ) THEN
1121	CALL module_interface_data_output_2d( &
1122	av, do2d(av,ivar), found, grid, mode, &
1123	local_pf, two_d, nzb_do, nzt_do, &
1124	fill_value &
1125	)
1126	ENDIF
1127
1128	resorted = .TRUE.
1129
1130	IF ( grid == 'zu' ) THEN
1131	IF ( mode == 'xy' ) level_z = zu
1132	ELSEIF ( grid == 'zw' ) THEN
1133	IF ( mode == 'xy' ) level_z = zw
1134	ELSEIF ( grid == 'zu1' ) THEN
1135	IF ( mode == 'xy' ) level_z(nzb+1) = zu(nzb+1)
1136	ELSEIF ( grid == 'zs' ) THEN
1137	IF ( mode == 'xy' ) level_z = zs
1138	ENDIF
1139
1140	IF ( .NOT. found ) THEN
1141	message_string = 'no output provided for: ' // &
1142	TRIM( do2d(av,ivar) )
1143	CALL message( 'data_output_2d', 'PA0181', 0, 0, 0, 6, 0 )
1144	ENDIF
1145
1146	END SELECT
1147
1148	!
1149	!-- Resort the array to be output, if not done above. Flag topography
1150	!-- grid points with fill values, using the corresponding maksing flag.
1151	IF ( .NOT. resorted ) THEN
1152	DO i = nxl, nxr
1153	DO j = nys, nyn
1154	DO k = nzb_do, nzt_do
1155	local_pf(i,j,k) = MERGE( to_be_resorted(k,j,i), &
1156	REAL( fill_value, KIND = wp ), &
1157	BTEST( wall_flags_static_0(k,j,i), &
1158	flag_nr ) )
1159	ENDDO
1160	ENDDO
1161	ENDDO
1162	ENDIF
1163
1164	!
1165	!-- Output of the individual cross-sections, depending on the cross-
1166	!-- section mode chosen.
1167	is = 1
1168	loop1: DO WHILE ( section(is,s_ind) /= -9999 .OR. two_d )
1169
1170	SELECT CASE ( mode )
1171
1172	CASE ( 'xy' )
1173	!
1174	!-- Determine the cross section index
1175	IF ( two_d ) THEN
1176	layer_xy = nzb+1
1177	ELSE
1178	layer_xy = section(is,s_ind)
1179	ENDIF
1180
1181	!
1182	!-- Exit the loop for layers beyond the data output domain
1183	!-- (used for soil model)
1184	IF ( layer_xy > nzt_do ) THEN
1185	EXIT loop1
1186	ENDIF
1187
1188	!
1189	!-- Update the netCDF xy cross section time axis.
1190	!-- In case of parallel output, this is only done by PE0
1191	!-- to increase the performance.
1192	IF ( time_since_reference_point /= do2d_xy_last_time(av) ) THEN
1193	do2d_xy_time_count(av) = do2d_xy_time_count(av) + 1
1194	do2d_xy_last_time(av) = time_since_reference_point
1195	IF ( myid == 0 ) THEN
1196	IF ( .NOT. data_output_2d_on_each_pe &
1197	.OR. netcdf_data_format > 4 ) &
1198	THEN
1199	#if defined( __netcdf )
1200	nc_stat = NF90_PUT_VAR( id_set_xy(av), &
1201	id_var_time_xy(av), &
1202	(/ time_since_reference_point /), &
1203	start = (/ do2d_xy_time_count(av) /), &
1204	count = (/ 1 /) )
1205	CALL netcdf_handle_error( 'data_output_2d', 53 )
1206	#endif
1207	ENDIF
1208	ENDIF
1209	ENDIF
1210	!
1211	!-- If required, carry out averaging along z
1212	IF ( section(is,s_ind) == -1 .AND. .NOT. two_d ) THEN
1213
1214	local_2d = 0.0_wp
1215	!
1216	!-- Carry out the averaging (all data are on the PE)
1217	DO k = nzb_do, nzt_do
1218	DO j = nys, nyn
1219	DO i = nxl, nxr
1220	local_2d(i,j) = local_2d(i,j) + local_pf(i,j,k)
1221	ENDDO
1222	ENDDO
1223	ENDDO
1224
1225	local_2d = local_2d / ( nzt_do - nzb_do + 1.0_wp)
1226
1227	ELSE
1228	!
1229	!-- Just store the respective section on the local array
1230	local_2d = local_pf(:,:,layer_xy)
1231
1232	ENDIF
1233
1234	#if defined( __parallel )
1235	IF ( netcdf_data_format > 4 ) THEN
1236	!
1237	!-- Parallel output in netCDF4/HDF5 format.
1238	IF ( two_d ) THEN
1239	iis = 1
1240	ELSE
1241	iis = is
1242	ENDIF
1243
1244	#if defined( __netcdf )
1245	!
1246	!-- For parallel output, all cross sections are first stored
1247	!-- here on a local array and will be written to the output
1248	!-- file afterwards to increase the performance.
1249	DO i = nxl, nxr
1250	DO j = nys, nyn
1251	local_2d_sections(i,j,iis) = local_2d(i,j)
1252	ENDDO
1253	ENDDO
1254	#endif
1255	ELSE
1256
1257	IF ( data_output_2d_on_each_pe ) THEN
1258	!
1259	!-- Output of partial arrays on each PE
1260	#if defined( __netcdf )
1261	IF ( myid == 0 ) THEN
1262	WRITE ( 21 ) time_since_reference_point, &
1263	do2d_xy_time_count(av), av
1264	ENDIF
1265	#endif
1266	DO i = 0, io_blocks-1
1267	IF ( i == io_group ) THEN
1268	WRITE ( 21 ) nxl, nxr, nys, nyn, nys, nyn
1269	WRITE ( 21 ) local_2d
1270	ENDIF
1271	#if defined( __parallel )
1272	CALL MPI_BARRIER( comm2d, ierr )
1273	#endif
1274	ENDDO
1275
1276	ELSE
1277	!
1278	!-- PE0 receives partial arrays from all processors and
1279	!-- then outputs them. Here a barrier has to be set,
1280	!-- because otherwise "-MPI- FATAL: Remote protocol queue
1281	!-- full" may occur.
1282	CALL MPI_BARRIER( comm2d, ierr )
1283
1284	ngp = ( nxr-nxl+1 ) * ( nyn-nys+1 )
1285	IF ( myid == 0 ) THEN
1286	!
1287	!-- Local array can be relocated directly.
1288	total_2d(nxl:nxr,nys:nyn) = local_2d
1289	!
1290	!-- Receive data from all other PEs.
1291	DO n = 1, numprocs-1
1292	!
1293	!-- Receive index limits first, then array.
1294	!-- Index limits are received in arbitrary order from
1295	!-- the PEs.
1296	CALL MPI_RECV( ind(1), 4, MPI_INTEGER, &
1297	MPI_ANY_SOURCE, 0, comm2d, &
1298	status, ierr )
1299	sender = status(MPI_SOURCE)
1300	DEALLOCATE( local_2d )
1301	ALLOCATE( local_2d(ind(1):ind(2),ind(3):ind(4)) )
1302	CALL MPI_RECV( local_2d(ind(1),ind(3)), ngp, &
1303	MPI_REAL, sender, 1, comm2d, &
1304	status, ierr )
1305	total_2d(ind(1):ind(2),ind(3):ind(4)) = local_2d
1306	ENDDO
1307	!
1308	!-- Relocate the local array for the next loop increment
1309	DEALLOCATE( local_2d )
1310	ALLOCATE( local_2d(nxl:nxr,nys:nyn) )
1311
1312	#if defined( __netcdf )
1313	IF ( two_d ) THEN
1314	nc_stat = NF90_PUT_VAR( id_set_xy(av), &
1315	id_var_do2d(av,ivar), &
1316	total_2d(0:nx,0:ny), &
1317	start = (/ 1, 1, 1, do2d_xy_time_count(av) /), &
1318	count = (/ nx+1, ny+1, 1, 1 /) )
1319	ELSE
1320	nc_stat = NF90_PUT_VAR( id_set_xy(av), &
1321	id_var_do2d(av,ivar), &
1322	total_2d(0:nx,0:ny), &
1323	start = (/ 1, 1, is, do2d_xy_time_count(av) /), &
1324	count = (/ nx+1, ny+1, 1, 1 /) )
1325	ENDIF
1326	CALL netcdf_handle_error( 'data_output_2d', 54 )
1327	#endif
1328
1329	ELSE
1330	!
1331	!-- First send the local index limits to PE0
1332	ind(1) = nxl; ind(2) = nxr
1333	ind(3) = nys; ind(4) = nyn
1334	CALL MPI_SEND( ind(1), 4, MPI_INTEGER, 0, 0, &
1335	comm2d, ierr )
1336	!
1337	!-- Send data to PE0
1338	CALL MPI_SEND( local_2d(nxl,nys), ngp, &
1339	MPI_REAL, 0, 1, comm2d, ierr )
1340	ENDIF
1341	!
1342	!-- A barrier has to be set, because otherwise some PEs may
1343	!-- proceed too fast so that PE0 may receive wrong data on
1344	!-- tag 0
1345	CALL MPI_BARRIER( comm2d, ierr )
1346	ENDIF
1347
1348	ENDIF
1349	#else
1350	#if defined( __netcdf )
1351	IF ( two_d ) THEN
1352	nc_stat = NF90_PUT_VAR( id_set_xy(av), &
1353	id_var_do2d(av,ivar), &
1354	local_2d(nxl:nxr,nys:nyn), &
1355	start = (/ 1, 1, 1, do2d_xy_time_count(av) /), &
1356	count = (/ nx+1, ny+1, 1, 1 /) )
1357	ELSE
1358	nc_stat = NF90_PUT_VAR( id_set_xy(av), &
1359	id_var_do2d(av,ivar), &
1360	local_2d(nxl:nxr,nys:nyn), &
1361	start = (/ 1, 1, is, do2d_xy_time_count(av) /), &
1362	count = (/ nx+1, ny+1, 1, 1 /) )
1363	ENDIF
1364	CALL netcdf_handle_error( 'data_output_2d', 447 )
1365	#endif
1366	#endif
1367
1368	!
1369	!-- For 2D-arrays (e.g. u*) only one cross-section is available.
1370	!-- Hence exit loop of output levels.
1371	IF ( two_d ) THEN
1372	IF ( netcdf_data_format < 5 ) two_d = .FALSE.
1373	EXIT loop1
1374	ENDIF
1375
1376	CASE ( 'xz' )
1377	!
1378	!-- Update the netCDF xz cross section time axis.
1379	!-- In case of parallel output, this is only done by PE0
1380	!-- to increase the performance.
1381	IF ( time_since_reference_point /= do2d_xz_last_time(av) ) THEN
1382	do2d_xz_time_count(av) = do2d_xz_time_count(av) + 1
1383	do2d_xz_last_time(av) = time_since_reference_point
1384	IF ( myid == 0 ) THEN
1385	IF ( .NOT. data_output_2d_on_each_pe &
1386	.OR. netcdf_data_format > 4 ) &
1387	THEN
1388	#if defined( __netcdf )
1389	nc_stat = NF90_PUT_VAR( id_set_xz(av), &
1390	id_var_time_xz(av), &
1391	(/ time_since_reference_point /), &
1392	start = (/ do2d_xz_time_count(av) /), &
1393	count = (/ 1 /) )
1394	CALL netcdf_handle_error( 'data_output_2d', 56 )
1395	#endif
1396	ENDIF
1397	ENDIF
1398	ENDIF
1399
1400	!
1401	!-- If required, carry out averaging along y
1402	IF ( section(is,s_ind) == -1 ) THEN
1403
1404	ALLOCATE( local_2d_l(nxl:nxr,nzb_do:nzt_do) )
1405	local_2d_l = 0.0_wp
1406	ngp = ( nxr-nxl + 1 ) * ( nzt_do-nzb_do + 1 )
1407	!
1408	!-- First local averaging on the PE
1409	DO k = nzb_do, nzt_do
1410	DO j = nys, nyn
1411	DO i = nxl, nxr
1412	local_2d_l(i,k) = local_2d_l(i,k) + &
1413	local_pf(i,j,k)
1414	ENDDO
1415	ENDDO
1416	ENDDO
1417	#if defined( __parallel )
1418	!
1419	!-- Now do the averaging over all PEs along y
1420	IF ( collective_wait ) CALL MPI_BARRIER( comm2d, ierr )
1421	CALL MPI_ALLREDUCE( local_2d_l(nxl,nzb_do), &
1422	local_2d(nxl,nzb_do), ngp, MPI_REAL, &
1423	MPI_SUM, comm1dy, ierr )
1424	#else
1425	local_2d = local_2d_l
1426	#endif
1427	local_2d = local_2d / ( ny + 1.0_wp )
1428
1429	DEALLOCATE( local_2d_l )
1430
1431	ELSE
1432	!
1433	!-- Just store the respective section on the local array
1434	!-- (but only if it is available on this PE!)
1435	IF ( section(is,s_ind) >= nys .AND. section(is,s_ind) <= nyn ) &
1436	THEN
1437	local_2d = local_pf(:,section(is,s_ind),nzb_do:nzt_do)
1438	ENDIF
1439
1440	ENDIF
1441
1442	#if defined( __parallel )
1443	IF ( netcdf_data_format > 4 ) THEN
1444	!
1445	!-- Output in netCDF4/HDF5 format.
1446	!-- Output only on those PEs where the respective cross
1447	!-- sections reside. Cross sections averaged along y are
1448	!-- output on the respective first PE along y (myidy=0).
1449	IF ( ( section(is,s_ind) >= nys .AND. &
1450	section(is,s_ind) <= nyn ) .OR. &
1451	( section(is,s_ind) == -1 .AND. myidy == 0 ) ) THEN
1452	#if defined( __netcdf )
1453	!
1454	!-- For parallel output, all cross sections are first
1455	!-- stored here on a local array and will be written to the
1456	!-- output file afterwards to increase the performance.
1457	DO i = nxl, nxr
1458	DO k = nzb_do, nzt_do
1459	local_2d_sections_l(i,is,k) = local_2d(i,k)
1460	ENDDO
1461	ENDDO
1462	#endif
1463	ENDIF
1464
1465	ELSE
1466
1467	IF ( data_output_2d_on_each_pe ) THEN
1468	!
1469	!-- Output of partial arrays on each PE. If the cross
1470	!-- section does not reside on the PE, output special
1471	!-- index values.
1472	#if defined( __netcdf )
1473	IF ( myid == 0 ) THEN
1474	WRITE ( 22 ) time_since_reference_point, &
1475	do2d_xz_time_count(av), av
1476	ENDIF
1477	#endif
1478	DO i = 0, io_blocks-1
1479	IF ( i == io_group ) THEN
1480	IF ( ( section(is,s_ind) >= nys .AND. &
1481	section(is,s_ind) <= nyn ) .OR. &
1482	( section(is,s_ind) == -1 .AND. &
1483	nys-1 == -1 ) ) &
1484	THEN
1485	WRITE (22) nxl, nxr, nzb_do, nzt_do, nzb, nzt+1
1486	WRITE (22) local_2d
1487	ELSE
1488	WRITE (22) -1, -1, -1, -1, -1, -1
1489	ENDIF
1490	ENDIF
1491	#if defined( __parallel )
1492	CALL MPI_BARRIER( comm2d, ierr )
1493	#endif
1494	ENDDO
1495
1496	ELSE
1497	!
1498	!-- PE0 receives partial arrays from all processors of the
1499	!-- respective cross section and outputs them. Here a
1500	!-- barrier has to be set, because otherwise
1501	!-- "-MPI- FATAL: Remote protocol queue full" may occur.
1502	CALL MPI_BARRIER( comm2d, ierr )
1503
1504	ngp = ( nxr-nxl + 1 ) * ( nzt_do-nzb_do + 1 )
1505	IF ( myid == 0 ) THEN
1506	!
1507	!-- Local array can be relocated directly.
1508	IF ( ( section(is,s_ind) >= nys .AND. &
1509	section(is,s_ind) <= nyn ) .OR. &
1510	( section(is,s_ind) == -1 .AND. &
1511	nys-1 == -1 ) ) THEN
1512	total_2d(nxl:nxr,nzb_do:nzt_do) = local_2d
1513	ENDIF
1514	!
1515	!-- Receive data from all other PEs.
1516	DO n = 1, numprocs-1
1517	!
1518	!-- Receive index limits first, then array.
1519	!-- Index limits are received in arbitrary order from
1520	!-- the PEs.
1521	CALL MPI_RECV( ind(1), 4, MPI_INTEGER, &
1522	MPI_ANY_SOURCE, 0, comm2d, &
1523	status, ierr )
1524	!
1525	!-- Not all PEs have data for XZ-cross-section.
1526	IF ( ind(1) /= -9999 ) THEN
1527	sender = status(MPI_SOURCE)
1528	DEALLOCATE( local_2d )
1529	ALLOCATE( local_2d(ind(1):ind(2), &
1530	ind(3):ind(4)) )
1531	CALL MPI_RECV( local_2d(ind(1),ind(3)), ngp, &
1532	MPI_REAL, sender, 1, comm2d, &
1533	status, ierr )
1534	total_2d(ind(1):ind(2),ind(3):ind(4)) = &
1535	local_2d
1536	ENDIF
1537	ENDDO
1538	!
1539	!-- Relocate the local array for the next loop increment
1540	DEALLOCATE( local_2d )
1541	ALLOCATE( local_2d(nxl:nxr,nzb_do:nzt_do) )
1542
1543	#if defined( __netcdf )
1544	nc_stat = NF90_PUT_VAR( id_set_xz(av), &
1545	id_var_do2d(av,ivar), &
1546	total_2d(0:nx,nzb_do:nzt_do), &
1547	start = (/ 1, is, 1, do2d_xz_time_count(av) /), &
1548	count = (/ nx+1, 1, nzt_do-nzb_do+1, 1 /) )
1549	CALL netcdf_handle_error( 'data_output_2d', 58 )
1550	#endif
1551
1552	ELSE
1553	!
1554	!-- If the cross section resides on the PE, send the
1555	!-- local index limits, otherwise send -9999 to PE0.
1556	IF ( ( section(is,s_ind) >= nys .AND. &
1557	section(is,s_ind) <= nyn ) .OR. &
1558	( section(is,s_ind) == -1 .AND. nys-1 == -1 ) ) &
1559	THEN
1560	ind(1) = nxl; ind(2) = nxr
1561	ind(3) = nzb_do; ind(4) = nzt_do
1562	ELSE
1563	ind(1) = -9999; ind(2) = -9999
1564	ind(3) = -9999; ind(4) = -9999
1565	ENDIF
1566	CALL MPI_SEND( ind(1), 4, MPI_INTEGER, 0, 0, &
1567	comm2d, ierr )
1568	!
1569	!-- If applicable, send data to PE0.
1570	IF ( ind(1) /= -9999 ) THEN
1571	CALL MPI_SEND( local_2d(nxl,nzb_do), ngp, &
1572	MPI_REAL, 0, 1, comm2d, ierr )
1573	ENDIF
1574	ENDIF
1575	!
1576	!-- A barrier has to be set, because otherwise some PEs may
1577	!-- proceed too fast so that PE0 may receive wrong data on
1578	!-- tag 0
1579	CALL MPI_BARRIER( comm2d, ierr )
1580	ENDIF
1581
1582	ENDIF
1583	#else
1584	#if defined( __netcdf )
1585	nc_stat = NF90_PUT_VAR( id_set_xz(av), &
1586	id_var_do2d(av,ivar), &
1587	local_2d(nxl:nxr,nzb_do:nzt_do), &
1588	start = (/ 1, is, 1, do2d_xz_time_count(av) /), &
1589	count = (/ nx+1, 1, nzt_do-nzb_do+1, 1 /) )
1590	CALL netcdf_handle_error( 'data_output_2d', 451 )
1591	#endif
1592	#endif
1593
1594	CASE ( 'yz' )
1595	!
1596	!-- Update the netCDF yz cross section time axis.
1597	!-- In case of parallel output, this is only done by PE0
1598	!-- to increase the performance.
1599	IF ( time_since_reference_point /= do2d_yz_last_time(av) ) THEN
1600	do2d_yz_time_count(av) = do2d_yz_time_count(av) + 1
1601	do2d_yz_last_time(av) = time_since_reference_point
1602	IF ( myid == 0 ) THEN
1603	IF ( .NOT. data_output_2d_on_each_pe &
1604	.OR. netcdf_data_format > 4 ) &
1605	THEN
1606	#if defined( __netcdf )
1607	nc_stat = NF90_PUT_VAR( id_set_yz(av), &
1608	id_var_time_yz(av), &
1609	(/ time_since_reference_point /), &
1610	start = (/ do2d_yz_time_count(av) /), &
1611	count = (/ 1 /) )
1612	CALL netcdf_handle_error( 'data_output_2d', 59 )
1613	#endif
1614	ENDIF
1615	ENDIF
1616	ENDIF
1617
1618	!
1619	!-- If required, carry out averaging along x
1620	IF ( section(is,s_ind) == -1 ) THEN
1621
1622	ALLOCATE( local_2d_l(nys:nyn,nzb_do:nzt_do) )
1623	local_2d_l = 0.0_wp
1624	ngp = ( nyn-nys+1 ) * ( nzt_do-nzb_do+1 )
1625	!
1626	!-- First local averaging on the PE
1627	DO k = nzb_do, nzt_do
1628	DO j = nys, nyn
1629	DO i = nxl, nxr
1630	local_2d_l(j,k) = local_2d_l(j,k) + &
1631	local_pf(i,j,k)
1632	ENDDO
1633	ENDDO
1634	ENDDO
1635	#if defined( __parallel )
1636	!
1637	!-- Now do the averaging over all PEs along x
1638	IF ( collective_wait ) CALL MPI_BARRIER( comm2d, ierr )
1639	CALL MPI_ALLREDUCE( local_2d_l(nys,nzb_do), &
1640	local_2d(nys,nzb_do), ngp, MPI_REAL, &
1641	MPI_SUM, comm1dx, ierr )
1642	#else
1643	local_2d = local_2d_l
1644	#endif
1645	local_2d = local_2d / ( nx + 1.0_wp )
1646
1647	DEALLOCATE( local_2d_l )
1648
1649	ELSE
1650	!
1651	!-- Just store the respective section on the local array
1652	!-- (but only if it is available on this PE!)
1653	IF ( section(is,s_ind) >= nxl .AND. section(is,s_ind) <= nxr ) &
1654	THEN
1655	local_2d = local_pf(section(is,s_ind),:,nzb_do:nzt_do)
1656	ENDIF
1657
1658	ENDIF
1659
1660	#if defined( __parallel )
1661	IF ( netcdf_data_format > 4 ) THEN
1662	!
1663	!-- Output in netCDF4/HDF5 format.
1664	!-- Output only on those PEs where the respective cross
1665	!-- sections reside. Cross sections averaged along x are
1666	!-- output on the respective first PE along x (myidx=0).
1667	IF ( ( section(is,s_ind) >= nxl .AND. &
1668	section(is,s_ind) <= nxr ) .OR. &
1669	( section(is,s_ind) == -1 .AND. myidx == 0 ) ) THEN
1670	#if defined( __netcdf )
1671	!
1672	!-- For parallel output, all cross sections are first
1673	!-- stored here on a local array and will be written to the
1674	!-- output file afterwards to increase the performance.
1675	DO j = nys, nyn
1676	DO k = nzb_do, nzt_do
1677	local_2d_sections_l(is,j,k) = local_2d(j,k)
1678	ENDDO
1679	ENDDO
1680	#endif
1681	ENDIF
1682
1683	ELSE
1684
1685	IF ( data_output_2d_on_each_pe ) THEN
1686	!
1687	!-- Output of partial arrays on each PE. If the cross
1688	!-- section does not reside on the PE, output special
1689	!-- index values.
1690	#if defined( __netcdf )
1691	IF ( myid == 0 ) THEN
1692	WRITE ( 23 ) time_since_reference_point, &
1693	do2d_yz_time_count(av), av
1694	ENDIF
1695	#endif
1696	DO i = 0, io_blocks-1
1697	IF ( i == io_group ) THEN
1698	IF ( ( section(is,s_ind) >= nxl .AND. &
1699	section(is,s_ind) <= nxr ) .OR. &
1700	( section(is,s_ind) == -1 .AND. &
1701	nxl-1 == -1 ) ) &
1702	THEN
1703	WRITE (23) nys, nyn, nzb_do, nzt_do, nzb, nzt+1
1704	WRITE (23) local_2d
1705	ELSE
1706	WRITE (23) -1, -1, -1, -1, -1, -1
1707	ENDIF
1708	ENDIF
1709	#if defined( __parallel )
1710	CALL MPI_BARRIER( comm2d, ierr )
1711	#endif
1712	ENDDO
1713
1714	ELSE
1715	!
1716	!-- PE0 receives partial arrays from all processors of the
1717	!-- respective cross section and outputs them. Here a
1718	!-- barrier has to be set, because otherwise
1719	!-- "-MPI- FATAL: Remote protocol queue full" may occur.
1720	CALL MPI_BARRIER( comm2d, ierr )
1721
1722	ngp = ( nyn-nys+1 ) * ( nzt_do-nzb_do+1 )
1723	IF ( myid == 0 ) THEN
1724	!
1725	!-- Local array can be relocated directly.
1726	IF ( ( section(is,s_ind) >= nxl .AND. &
1727	section(is,s_ind) <= nxr ) .OR. &
1728	( section(is,s_ind) == -1 .AND. nxl-1 == -1 ) ) &
1729	THEN
1730	total_2d(nys:nyn,nzb_do:nzt_do) = local_2d
1731	ENDIF
1732	!
1733	!-- Receive data from all other PEs.
1734	DO n = 1, numprocs-1
1735	!
1736	!-- Receive index limits first, then array.
1737	!-- Index limits are received in arbitrary order from
1738	!-- the PEs.
1739	CALL MPI_RECV( ind(1), 4, MPI_INTEGER, &
1740	MPI_ANY_SOURCE, 0, comm2d, &
1741	status, ierr )
1742	!
1743	!-- Not all PEs have data for YZ-cross-section.
1744	IF ( ind(1) /= -9999 ) THEN
1745	sender = status(MPI_SOURCE)
1746	DEALLOCATE( local_2d )
1747	ALLOCATE( local_2d(ind(1):ind(2), &
1748	ind(3):ind(4)) )
1749	CALL MPI_RECV( local_2d(ind(1),ind(3)), ngp, &
1750	MPI_REAL, sender, 1, comm2d, &
1751	status, ierr )
1752	total_2d(ind(1):ind(2),ind(3):ind(4)) = &
1753	local_2d
1754	ENDIF
1755	ENDDO
1756	!
1757	!-- Relocate the local array for the next loop increment
1758	DEALLOCATE( local_2d )
1759	ALLOCATE( local_2d(nys:nyn,nzb_do:nzt_do) )
1760
1761	#if defined( __netcdf )
1762	nc_stat = NF90_PUT_VAR( id_set_yz(av), &
1763	id_var_do2d(av,ivar), &
1764	total_2d(0:ny,nzb_do:nzt_do), &
1765	start = (/ is, 1, 1, do2d_yz_time_count(av) /), &
1766	count = (/ 1, ny+1, nzt_do-nzb_do+1, 1 /) )
1767	CALL netcdf_handle_error( 'data_output_2d', 61 )
1768	#endif
1769
1770	ELSE
1771	!
1772	!-- If the cross section resides on the PE, send the
1773	!-- local index limits, otherwise send -9999 to PE0.
1774	IF ( ( section(is,s_ind) >= nxl .AND. &
1775	section(is,s_ind) <= nxr ) .OR. &
1776	( section(is,s_ind) == -1 .AND. nxl-1 == -1 ) ) &
1777	THEN
1778	ind(1) = nys; ind(2) = nyn
1779	ind(3) = nzb_do; ind(4) = nzt_do
1780	ELSE
1781	ind(1) = -9999; ind(2) = -9999
1782	ind(3) = -9999; ind(4) = -9999
1783	ENDIF
1784	CALL MPI_SEND( ind(1), 4, MPI_INTEGER, 0, 0, &
1785	comm2d, ierr )
1786	!
1787	!-- If applicable, send data to PE0.
1788	IF ( ind(1) /= -9999 ) THEN
1789	CALL MPI_SEND( local_2d(nys,nzb_do), ngp, &
1790	MPI_REAL, 0, 1, comm2d, ierr )
1791	ENDIF
1792	ENDIF
1793	!
1794	!-- A barrier has to be set, because otherwise some PEs may
1795	!-- proceed too fast so that PE0 may receive wrong data on
1796	!-- tag 0
1797	CALL MPI_BARRIER( comm2d, ierr )
1798	ENDIF
1799
1800	ENDIF
1801	#else
1802	#if defined( __netcdf )
1803	nc_stat = NF90_PUT_VAR( id_set_yz(av), &
1804	id_var_do2d(av,ivar), &
1805	local_2d(nys:nyn,nzb_do:nzt_do), &
1806	start = (/ is, 1, 1, do2d_xz_time_count(av) /), &
1807	count = (/ 1, ny+1, nzt_do-nzb_do+1, 1 /) )
1808	CALL netcdf_handle_error( 'data_output_2d', 452 )
1809	#endif
1810	#endif
1811
1812	END SELECT
1813
1814	is = is + 1
1815	ENDDO loop1
1816
1817	!
1818	!-- For parallel output, all data were collected before on a local array
1819	!-- and are written now to the netcdf file. This must be done to increase
1820	!-- the performance of the parallel output.
1821	#if defined( __netcdf )
1822	IF ( netcdf_data_format > 4 ) THEN
1823
1824	SELECT CASE ( mode )
1825
1826	CASE ( 'xy' )
1827	IF ( two_d ) THEN
1828	nis = 1
1829	two_d = .FALSE.
1830	ELSE
1831	nis = ns
1832	ENDIF
1833	!
1834	!-- Do not output redundant ghost point data except for the
1835	!-- boundaries of the total domain.
1836	! IF ( nxr == nx .AND. nyn /= ny ) THEN
1837	! nc_stat = NF90_PUT_VAR( id_set_xy(av), &
1838	! id_var_do2d(av,ivar), &
1839	! local_2d_sections(nxl:nxr+1, &
1840	! nys:nyn,1:nis), &
1841	! start = (/ nxl+1, nys+1, 1, &
1842	! do2d_xy_time_count(av) /), &
1843	! count = (/ nxr-nxl+2, &
1844	! nyn-nys+1, nis, 1 &
1845	! /) )
1846	! ELSEIF ( nxr /= nx .AND. nyn == ny ) THEN
1847	! nc_stat = NF90_PUT_VAR( id_set_xy(av), &
1848	! id_var_do2d(av,ivar), &
1849	! local_2d_sections(nxl:nxr, &
1850	! nys:nyn+1,1:nis), &
1851	! start = (/ nxl+1, nys+1, 1, &
1852	! do2d_xy_time_count(av) /), &
1853	! count = (/ nxr-nxl+1, &
1854	! nyn-nys+2, nis, 1 &
1855	! /) )
1856	! ELSEIF ( nxr == nx .AND. nyn == ny ) THEN
1857	! nc_stat = NF90_PUT_VAR( id_set_xy(av), &
1858	! id_var_do2d(av,ivar), &
1859	! local_2d_sections(nxl:nxr+1, &
1860	! nys:nyn+1,1:nis), &
1861	! start = (/ nxl+1, nys+1, 1, &
1862	! do2d_xy_time_count(av) /), &
1863	! count = (/ nxr-nxl+2, &
1864	! nyn-nys+2, nis, 1 &
1865	! /) )
1866	! ELSE
1867	nc_stat = NF90_PUT_VAR( id_set_xy(av), &
1868	id_var_do2d(av,ivar), &
1869	local_2d_sections(nxl:nxr, &
1870	nys:nyn,1:nis), &
1871	start = (/ nxl+1, nys+1, 1, &
1872	do2d_xy_time_count(av) /), &
1873	count = (/ nxr-nxl+1, &
1874	nyn-nys+1, nis, 1 &
1875	/) )
1876	! ENDIF
1877
1878	CALL netcdf_handle_error( 'data_output_2d', 55 )
1879
1880	CASE ( 'xz' )
1881	!
1882	!-- First, all PEs get the information of all cross-sections.
1883	!-- Then the data are written to the output file by all PEs
1884	!-- while NF90_COLLECTIVE is set in subroutine
1885	!-- define_netcdf_header. Although redundant information are
1886	!-- written to the output file in that case, the performance
1887	!-- is significantly better compared to the case where only
1888	!-- the first row of PEs in x-direction (myidx = 0) is given
1889	!-- the output while NF90_INDEPENDENT is set.
1890	IF ( npey /= 1 ) THEN
1891
1892	#if defined( __parallel )
1893	!
1894	!-- Distribute data over all PEs along y
1895	ngp = ( nxr-nxl+1 ) * ( nzt_do-nzb_do+1 ) * ns
1896	IF ( collective_wait ) CALL MPI_BARRIER( comm2d, ierr )
1897	CALL MPI_ALLREDUCE( local_2d_sections_l(nxl,1,nzb_do), &
1898	local_2d_sections(nxl,1,nzb_do), &
1899	ngp, MPI_REAL, MPI_SUM, comm1dy, &
1900	ierr )
1901	#else
1902	local_2d_sections = local_2d_sections_l
1903	#endif
1904	ENDIF
1905	!
1906	!-- Do not output redundant ghost point data except for the
1907	!-- boundaries of the total domain.
1908	! IF ( nxr == nx ) THEN
1909	! nc_stat = NF90_PUT_VAR( id_set_xz(av), &
1910	! id_var_do2d(av,ivar), &
1911	! local_2d_sections(nxl:nxr+1,1:ns, &
1912	! nzb_do:nzt_do), &
1913	! start = (/ nxl+1, 1, 1, &
1914	! do2d_xz_time_count(av) /), &
1915	! count = (/ nxr-nxl+2, ns, nzt_do-nzb_do+1, &
1916	! 1 /) )
1917	! ELSE
1918	nc_stat = NF90_PUT_VAR( id_set_xz(av), &
1919	id_var_do2d(av,ivar), &
1920	local_2d_sections(nxl:nxr,1:ns, &
1921	nzb_do:nzt_do), &
1922	start = (/ nxl+1, 1, 1, &
1923	do2d_xz_time_count(av) /), &
1924	count = (/ nxr-nxl+1, ns, nzt_do-nzb_do+1, &
1925	1 /) )
1926	! ENDIF
1927
1928	CALL netcdf_handle_error( 'data_output_2d', 57 )
1929
1930	CASE ( 'yz' )
1931	!
1932	!-- First, all PEs get the information of all cross-sections.
1933	!-- Then the data are written to the output file by all PEs
1934	!-- while NF90_COLLECTIVE is set in subroutine
1935	!-- define_netcdf_header. Although redundant information are
1936	!-- written to the output file in that case, the performance
1937	!-- is significantly better compared to the case where only
1938	!-- the first row of PEs in y-direction (myidy = 0) is given
1939	!-- the output while NF90_INDEPENDENT is set.
1940	IF ( npex /= 1 ) THEN
1941
1942	#if defined( __parallel )
1943	!
1944	!-- Distribute data over all PEs along x
1945	ngp = ( nyn-nys+1 ) * ( nzt-nzb + 2 ) * ns
1946	IF ( collective_wait ) CALL MPI_BARRIER( comm2d, ierr )
1947	CALL MPI_ALLREDUCE( local_2d_sections_l(1,nys,nzb_do), &
1948	local_2d_sections(1,nys,nzb_do), &
1949	ngp, MPI_REAL, MPI_SUM, comm1dx, &
1950	ierr )
1951	#else
1952	local_2d_sections = local_2d_sections_l
1953	#endif
1954	ENDIF
1955	!
1956	!-- Do not output redundant ghost point data except for the
1957	!-- boundaries of the total domain.
1958	! IF ( nyn == ny ) THEN
1959	! nc_stat = NF90_PUT_VAR( id_set_yz(av), &
1960	! id_var_do2d(av,ivar), &
1961	! local_2d_sections(1:ns, &
1962	! nys:nyn+1,nzb_do:nzt_do), &
1963	! start = (/ 1, nys+1, 1, &
1964	! do2d_yz_time_count(av) /), &
1965	! count = (/ ns, nyn-nys+2, &
1966	! nzt_do-nzb_do+1, 1 /) )
1967	! ELSE
1968	nc_stat = NF90_PUT_VAR( id_set_yz(av), &
1969	id_var_do2d(av,ivar), &
1970	local_2d_sections(1:ns,nys:nyn, &
1971	nzb_do:nzt_do), &
1972	start = (/ 1, nys+1, 1, &
1973	do2d_yz_time_count(av) /), &
1974	count = (/ ns, nyn-nys+1, &
1975	nzt_do-nzb_do+1, 1 /) )
1976	! ENDIF
1977
1978	CALL netcdf_handle_error( 'data_output_2d', 60 )
1979
1980	CASE DEFAULT
1981	message_string = 'unknown cross-section: ' // TRIM( mode )
1982	CALL message( 'data_output_2d', 'PA0180', 1, 2, 0, 6, 0 )
1983
1984	END SELECT
1985
1986	ENDIF
1987	#endif
1988	ENDIF
1989
1990	ivar = ivar + 1
1991	l = MAX( 2, LEN_TRIM( do2d(av,ivar) ) )
1992	do2d_mode = do2d(av,ivar)(l-1:l)
1993
1994	ENDDO
1995
1996	!
1997	!-- Deallocate temporary arrays.
1998	IF ( ALLOCATED( level_z ) ) DEALLOCATE( level_z )
1999	IF ( netcdf_data_format > 4 ) THEN
2000	DEALLOCATE( local_pf, local_2d, local_2d_sections )
2001	IF( mode == 'xz' .OR. mode == 'yz' ) DEALLOCATE( local_2d_sections_l )
2002	ENDIF
2003	#if defined( __parallel )
2004	IF ( .NOT. data_output_2d_on_each_pe .AND. myid == 0 ) THEN
2005	DEALLOCATE( total_2d )
2006	ENDIF
2007	#endif
2008
2009	!
2010	!-- Close plot output file.
2011	file_id = 20 + s_ind
2012
2013	IF ( data_output_2d_on_each_pe ) THEN
2014	DO i = 0, io_blocks-1
2015	IF ( i == io_group ) THEN
2016	CALL close_file( file_id )
2017	ENDIF
2018	#if defined( __parallel )
2019	CALL MPI_BARRIER( comm2d, ierr )
2020	#endif
2021	ENDDO
2022	ELSE
2023	IF ( myid == 0 ) CALL close_file( file_id )
2024	ENDIF
2025
2026	CALL cpu_log( log_point(3), 'data_output_2d', 'stop' )
2027
2028	IF ( debug_output_timestep ) CALL debug_message( 'data_output_2d', 'end' )
2029
2030
2031	END SUBROUTINE data_output_2d

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |