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

Last change on this file since 1972 was 1972, checked in by maronga, 8 years ago

further modularization of land surface model (2D/3D output and restart data)

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