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

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

last commit documented

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