Home

Context Navigation

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

Last change on this file since 1806 was 1789, checked in by maronga, 9 years ago
last commit documented
Property svn:keywords set to `Id`
File size: 88.5 KB

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

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |