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

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

further modularization of land surface model (2D/3D output and restart data). Bugfix for restart runs without land surface model

  • Property svn:keywords set to Id
File size: 77.2 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 radiation quantities is now done directly in the respective module
22!
23! Former revisions:
24! -----------------
25! $Id: data_output_2d.f90 1976 2016-07-27 13:28:04Z maronga $
26!
27! 1972 2016-07-26 07:52:02Z maronga
28! Output of land surface quantities is now done directly in the respective
29! module
30!
31! 1960 2016-07-12 16:34:24Z suehring
32! Scalar surface flux added
33! Rename INTEGER variable s into s_ind, as s is already assigned to scalar
34!
35! 1849 2016-04-08 11:33:18Z hoffmann
36! precipitation_amount, precipitation_rate, prr moved to arrays_3d
37!
38! 1822 2016-04-07 07:49:42Z hoffmann
39! Output of bulk cloud physics simplified.
40!
41! 1788 2016-03-10 11:01:04Z maronga
42! Added output of z0q
43!
44! 1783 2016-03-06 18:36:17Z raasch
45! name change of netcdf routines and module + related changes
46!
47! 1745 2016-02-05 13:06:51Z gronemeier
48! Bugfix: test if time axis limit exceeds moved to point after call of check_open
49!
50! 1703 2015-11-02 12:38:44Z raasch
51! bugfix for output of single (*) xy-sections in case of parallel netcdf I/O
52!
53! 1701 2015-11-02 07:43:04Z maronga
54! Bugfix in output of RRTGM data
55!
56! 1691 2015-10-26 16:17:44Z maronga
57! Added output of Obukhov length (ol) and radiative heating rates  for RRTMG.
58! Formatting corrections.
59!
60! 1682 2015-10-07 23:56:08Z knoop
61! Code annotations made doxygen readable
62!
63! 1585 2015-04-30 07:05:52Z maronga
64! Added support for RRTMG
65!
66! 1555 2015-03-04 17:44:27Z maronga
67! Added output of r_a and r_s
68!
69! 1551 2015-03-03 14:18:16Z maronga
70! Added suppport for land surface model and radiation model output. In the course
71! of this action, the limits for vertical loops have been changed (from nzb and
72! nzt+1 to nzb_do and nzt_do, respectively in order to allow soil model output).
73! Moreover, a new vertical grid zs was introduced.
74!
75! 1359 2014-04-11 17:15:14Z hoffmann
76! New particle structure integrated.
77!
78! 1353 2014-04-08 15:21:23Z heinze
79! REAL constants provided with KIND-attribute
80!
81! 1327 2014-03-21 11:00:16Z raasch
82! parts concerning iso2d output removed,
83! -netcdf output queries
84!
85! 1320 2014-03-20 08:40:49Z raasch
86! ONLY-attribute added to USE-statements,
87! kind-parameters added to all INTEGER and REAL declaration statements,
88! kinds are defined in new module kinds,
89! revision history before 2012 removed,
90! comment fields (!:) to be used for variable explanations added to
91! all variable declaration statements
92!
93! 1318 2014-03-17 13:35:16Z raasch
94! barrier argument removed from cpu_log.
95! module interfaces removed
96!
97! 1311 2014-03-14 12:13:39Z heinze
98! bugfix: close #if defined( __netcdf )
99!
100! 1308 2014-03-13 14:58:42Z fricke
101! +local_2d_sections, local_2d_sections_l, ns
102! Check, if the limit of the time dimension is exceeded for parallel output
103! To increase the performance for parallel output, the following is done:
104! - Update of time axis is only done by PE0
105! - Cross sections are first stored on a local array and are written
106!   collectively to the output file by all PEs.
107!
108! 1115 2013-03-26 18:16:16Z hoffmann
109! ql is calculated by calc_liquid_water_content
110!
111! 1076 2012-12-05 08:30:18Z hoffmann
112! Bugfix in output of ql
113!
114! 1065 2012-11-22 17:42:36Z hoffmann
115! Bugfix: Output of cross sections of ql
116!
117! 1053 2012-11-13 17:11:03Z hoffmann
118! +qr, nr, qc and cross sections
119!
120! 1036 2012-10-22 13:43:42Z raasch
121! code put under GPL (PALM 3.9)
122!
123! 1031 2012-10-19 14:35:30Z raasch
124! netCDF4 without parallel file support implemented
125!
126! 1007 2012-09-19 14:30:36Z franke
127! Bugfix: missing calculation of ql_vp added
128!
129! 978 2012-08-09 08:28:32Z fricke
130! +z0h
131!
132! Revision 1.1  1997/08/11 06:24:09  raasch
133! Initial revision
134!
135!
136! Description:
137! ------------
138!> Data output of horizontal cross-sections in netCDF format or binary format
139!> compatible to old graphic software iso2d.
140!> Attention: The position of the sectional planes is still not always computed
141!> ---------  correctly. (zu is used always)!
142!------------------------------------------------------------------------------!
143 SUBROUTINE data_output_2d( mode, av )
144 
145
146    USE arrays_3d,                                                             &
147        ONLY:  dzw, e, nr, ol, p, pt, precipitation_amount, precipitation_rate,&
148               prr,q, qc, ql, ql_c, ql_v, ql_vp, qr, qsws, rho, s, sa, shf,    &
149               ssws, tend, ts, u, us, v, vpt, w, z0, z0h, z0q, zu, zw
150       
151    USE averaging
152       
153    USE cloud_parameters,                                                      &
154        ONLY:  hyrho, l_d_cp, pt_d_t
155               
156    USE control_parameters,                                                    &
157        ONLY:  cloud_physics, data_output_2d_on_each_pe, data_output_xy,       &
158               data_output_xz, data_output_yz, do2d,                           &
159               do2d_xy_last_time, do2d_xy_n, do2d_xy_time_count,               &
160               do2d_xz_last_time, do2d_xz_n, do2d_xz_time_count,               &
161               do2d_yz_last_time, do2d_yz_n, do2d_yz_time_count,               &
162               ibc_uv_b, io_blocks, io_group, message_string,                  &
163               ntdim_2d_xy, ntdim_2d_xz, ntdim_2d_yz,                          &
164               psolver, section, simulated_time, simulated_time_chr,           &
165               time_since_reference_point
166       
167    USE cpulog,                                                                &
168        ONLY:  cpu_log, log_point 
169       
170    USE grid_variables,                                                        &
171        ONLY:  dx, dy
172       
173    USE indices,                                                               &
174        ONLY:  nbgp, nx, nxl, nxlg, nxr, nxrg, ny, nyn, nyng, nys, nysg,       &
175               nz, nzb, nzt
176               
177    USE kinds
178   
179    USE land_surface_model_mod,                                                &
180        ONLY:  land_surface, lsm_data_output_2d, zs
181   
182#if defined( __netcdf )
183    USE NETCDF
184#endif
185
186    USE netcdf_interface,                                                      &
187        ONLY:  id_set_xy, id_set_xz, id_set_yz, id_var_do2d, id_var_time_xy,   &
188               id_var_time_xz, id_var_time_yz, nc_stat, netcdf_data_format,    &
189               netcdf_handle_error
190
191    USE particle_attributes,                                                   &
192        ONLY:  grid_particles, number_of_particles, particle_advection_start,  &
193               particles, prt_count
194   
195    USE pegrid
196
197    USE radiation_model_mod,                                                   &
198        ONLY:  radiation, radiation_data_output_2d
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
753
754             CASE ( 'rho_xy', 'rho_xz', 'rho_yz' )
755                IF ( av == 0 )  THEN
756                   to_be_resorted => rho
757                ELSE
758                   to_be_resorted => rho_av
759                ENDIF
760
761             CASE ( 's_xy', 's_xz', 's_yz' )
762                IF ( av == 0 )  THEN
763                   to_be_resorted => s
764                ELSE
765                   to_be_resorted => s_av
766                ENDIF
767
768             CASE ( 'sa_xy', 'sa_xz', 'sa_yz' )
769                IF ( av == 0 )  THEN
770                   to_be_resorted => sa
771                ELSE
772                   to_be_resorted => sa_av
773                ENDIF
774
775             CASE ( 'shf*_xy' )        ! 2d-array
776                IF ( av == 0 ) THEN
777                   DO  i = nxlg, nxrg
778                      DO  j = nysg, nyng
779                         local_pf(i,j,nzb+1) =  shf(j,i)
780                      ENDDO
781                   ENDDO
782                ELSE
783                   DO  i = nxlg, nxrg
784                      DO  j = nysg, nyng
785                         local_pf(i,j,nzb+1) =  shf_av(j,i)
786                      ENDDO
787                   ENDDO
788                ENDIF
789                resorted = .TRUE.
790                two_d = .TRUE.
791                level_z(nzb+1) = zu(nzb+1)
792               
793             CASE ( 'ssws*_xy' )        ! 2d-array
794                IF ( av == 0 ) THEN
795                   DO  i = nxlg, nxrg
796                      DO  j = nysg, nyng
797                         local_pf(i,j,nzb+1) =  ssws(j,i)
798                      ENDDO
799                   ENDDO
800                ELSE
801                   DO  i = nxlg, nxrg
802                      DO  j = nysg, nyng 
803                         local_pf(i,j,nzb+1) =  ssws_av(j,i)
804                      ENDDO
805                   ENDDO
806                ENDIF
807                resorted = .TRUE.
808                two_d = .TRUE.
809                level_z(nzb+1) = zu(nzb+1)               
810
811             CASE ( 't*_xy' )        ! 2d-array
812                IF ( av == 0 )  THEN
813                   DO  i = nxlg, nxrg
814                      DO  j = nysg, nyng
815                         local_pf(i,j,nzb+1) = ts(j,i)
816                      ENDDO
817                   ENDDO
818                ELSE
819                   DO  i = nxlg, nxrg
820                      DO  j = nysg, nyng
821                         local_pf(i,j,nzb+1) = ts_av(j,i)
822                      ENDDO
823                   ENDDO
824                ENDIF
825                resorted = .TRUE.
826                two_d = .TRUE.
827                level_z(nzb+1) = zu(nzb+1)
828
829             CASE ( 'u_xy', 'u_xz', 'u_yz' )
830                IF ( av == 0 )  THEN
831                   to_be_resorted => u
832                ELSE
833                   to_be_resorted => u_av
834                ENDIF
835                IF ( mode == 'xy' )  level_z = zu
836!
837!--             Substitute the values generated by "mirror" boundary condition
838!--             at the bottom boundary by the real surface values.
839                IF ( do2d(av,if) == 'u_xz'  .OR.  do2d(av,if) == 'u_yz' )  THEN
840                   IF ( ibc_uv_b == 0 )  local_pf(:,:,nzb) = 0.0_wp
841                ENDIF
842
843             CASE ( 'u*_xy' )        ! 2d-array
844                IF ( av == 0 )  THEN
845                   DO  i = nxlg, nxrg
846                      DO  j = nysg, nyng
847                         local_pf(i,j,nzb+1) = us(j,i)
848                      ENDDO
849                   ENDDO
850                ELSE
851                   DO  i = nxlg, nxrg
852                      DO  j = nysg, nyng
853                         local_pf(i,j,nzb+1) = us_av(j,i)
854                      ENDDO
855                   ENDDO
856                ENDIF
857                resorted = .TRUE.
858                two_d = .TRUE.
859                level_z(nzb+1) = zu(nzb+1)
860
861             CASE ( 'v_xy', 'v_xz', 'v_yz' )
862                IF ( av == 0 )  THEN
863                   to_be_resorted => v
864                ELSE
865                   to_be_resorted => v_av
866                ENDIF
867                IF ( mode == 'xy' )  level_z = zu
868!
869!--             Substitute the values generated by "mirror" boundary condition
870!--             at the bottom boundary by the real surface values.
871                IF ( do2d(av,if) == 'v_xz'  .OR.  do2d(av,if) == 'v_yz' )  THEN
872                   IF ( ibc_uv_b == 0 )  local_pf(:,:,nzb) = 0.0_wp
873                ENDIF
874
875             CASE ( 'vpt_xy', 'vpt_xz', 'vpt_yz' )
876                IF ( av == 0 )  THEN
877                   to_be_resorted => vpt
878                ELSE
879                   to_be_resorted => vpt_av
880                ENDIF
881                IF ( mode == 'xy' )  level_z = zu
882
883             CASE ( 'w_xy', 'w_xz', 'w_yz' )
884                IF ( av == 0 )  THEN
885                   to_be_resorted => w
886                ELSE
887                   to_be_resorted => w_av
888                ENDIF
889                IF ( mode == 'xy' )  level_z = zw
890
891             CASE ( 'z0*_xy' )        ! 2d-array
892                IF ( av == 0 ) THEN
893                   DO  i = nxlg, nxrg
894                      DO  j = nysg, nyng
895                         local_pf(i,j,nzb+1) =  z0(j,i)
896                      ENDDO
897                   ENDDO
898                ELSE
899                   DO  i = nxlg, nxrg
900                      DO  j = nysg, nyng
901                         local_pf(i,j,nzb+1) =  z0_av(j,i)
902                      ENDDO
903                   ENDDO
904                ENDIF
905                resorted = .TRUE.
906                two_d = .TRUE.
907                level_z(nzb+1) = zu(nzb+1)
908
909             CASE ( 'z0h*_xy' )        ! 2d-array
910                IF ( av == 0 ) THEN
911                   DO  i = nxlg, nxrg
912                      DO  j = nysg, nyng
913                         local_pf(i,j,nzb+1) =  z0h(j,i)
914                      ENDDO
915                   ENDDO
916                ELSE
917                   DO  i = nxlg, nxrg
918                      DO  j = nysg, nyng
919                         local_pf(i,j,nzb+1) =  z0h_av(j,i)
920                      ENDDO
921                   ENDDO
922                ENDIF
923                resorted = .TRUE.
924                two_d = .TRUE.
925                level_z(nzb+1) = zu(nzb+1)
926
927             CASE ( 'z0q*_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) =  z0q(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) =  z0q_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 DEFAULT
946
947!
948!--             Land surface model quantity
949                IF ( land_surface )  THEN
950                   CALL lsm_data_output_2d( av, do2d(av,if), found, grid, mode,&
951                                            local_pf, two_d, nzb_do, nzt_do )
952                ENDIF
953
954!
955!--             Radiation quantity
956                IF ( .NOT. found  .AND.  radiation )  THEN
957                   CALL radiation_data_output_2d( av, do2d(av,if), found, grid,&
958                                                  mode, local_pf, two_d  )
959                ENDIF
960
961!
962!--             User defined quantity
963                IF ( .NOT. found )  THEN
964                   CALL user_data_output_2d( av, do2d(av,if), found, grid,     &
965                                             local_pf, two_d, nzb_do, nzt_do )
966                ENDIF
967
968                resorted = .TRUE.
969
970                IF ( grid == 'zu' )  THEN
971                   IF ( mode == 'xy' )  level_z = zu
972                ELSEIF ( grid == 'zw' )  THEN
973                   IF ( mode == 'xy' )  level_z = zw
974                ELSEIF ( grid == 'zu1' ) THEN
975                   IF ( mode == 'xy' )  level_z(nzb+1) = zu(nzb+1)
976                ELSEIF ( grid == 'zs' ) THEN
977                   IF ( mode == 'xy' )  level_z = zs
978                ENDIF
979
980                IF ( .NOT. found )  THEN
981                   message_string = 'no output provided for: ' //              &
982                                    TRIM( do2d(av,if) )
983                   CALL message( 'data_output_2d', 'PA0181', 0, 0, 0, 6, 0 )
984                ENDIF
985
986          END SELECT
987
988!
989!--       Resort the array to be output, if not done above
990          IF ( .NOT. resorted )  THEN
991             DO  i = nxlg, nxrg
992                DO  j = nysg, nyng
993                   DO  k = nzb_do, nzt_do
994                      local_pf(i,j,k) = to_be_resorted(k,j,i)
995                   ENDDO
996                ENDDO
997             ENDDO
998          ENDIF
999
1000!
1001!--       Output of the individual cross-sections, depending on the cross-
1002!--       section mode chosen.
1003          is = 1
1004   loop1: DO WHILE ( section(is,s_ind) /= -9999  .OR.  two_d )
1005
1006             SELECT CASE ( mode )
1007
1008                CASE ( 'xy' )
1009!
1010!--                Determine the cross section index
1011                   IF ( two_d )  THEN
1012                      layer_xy = nzb+1
1013                   ELSE
1014                      layer_xy = section(is,s_ind)
1015                   ENDIF
1016
1017!
1018!--                Exit the loop for layers beyond the data output domain
1019!--                (used for soil model)
1020                   IF ( layer_xy > nzt_do )  THEN
1021                      EXIT loop1
1022                   ENDIF
1023
1024!
1025!--                Update the netCDF xy cross section time axis.
1026!--                In case of parallel output, this is only done by PE0
1027!--                to increase the performance.
1028                   IF ( simulated_time /= do2d_xy_last_time(av) )  THEN
1029                      do2d_xy_time_count(av) = do2d_xy_time_count(av) + 1
1030                      do2d_xy_last_time(av)  = simulated_time
1031                      IF ( myid == 0 )  THEN
1032                         IF ( .NOT. data_output_2d_on_each_pe  &
1033                              .OR.  netcdf_data_format > 4 )   &
1034                         THEN
1035#if defined( __netcdf )
1036                            nc_stat = NF90_PUT_VAR( id_set_xy(av),             &
1037                                                    id_var_time_xy(av),        &
1038                                             (/ time_since_reference_point /), &
1039                                         start = (/ do2d_xy_time_count(av) /), &
1040                                                    count = (/ 1 /) )
1041                            CALL netcdf_handle_error( 'data_output_2d', 53 )
1042#endif
1043                         ENDIF
1044                      ENDIF
1045                   ENDIF
1046!
1047!--                If required, carry out averaging along z
1048                   IF ( section(is,s_ind) == -1  .AND.  .NOT. two_d )  THEN
1049
1050                      local_2d = 0.0_wp
1051!
1052!--                   Carry out the averaging (all data are on the PE)
1053                      DO  k = nzb_do, nzt_do
1054                         DO  j = nysg, nyng
1055                            DO  i = nxlg, nxrg
1056                               local_2d(i,j) = local_2d(i,j) + local_pf(i,j,k)
1057                            ENDDO
1058                         ENDDO
1059                      ENDDO
1060
1061                      local_2d = local_2d / ( nzt_do - nzb_do + 1.0_wp)
1062
1063                   ELSE
1064!
1065!--                   Just store the respective section on the local array
1066                      local_2d = local_pf(:,:,layer_xy)
1067
1068                   ENDIF
1069
1070#if defined( __parallel )
1071                   IF ( netcdf_data_format > 4 )  THEN
1072!
1073!--                   Parallel output in netCDF4/HDF5 format.
1074                      IF ( two_d ) THEN
1075                         iis = 1
1076                      ELSE
1077                         iis = is
1078                      ENDIF
1079
1080#if defined( __netcdf )
1081!
1082!--                   For parallel output, all cross sections are first stored
1083!--                   here on a local array and will be written to the output
1084!--                   file afterwards to increase the performance.
1085                      DO  i = nxlg, nxrg
1086                         DO  j = nysg, nyng
1087                            local_2d_sections(i,j,iis) = local_2d(i,j)
1088                         ENDDO
1089                      ENDDO
1090#endif
1091                   ELSE
1092
1093                      IF ( data_output_2d_on_each_pe )  THEN
1094!
1095!--                      Output of partial arrays on each PE
1096#if defined( __netcdf )
1097                         IF ( myid == 0 )  THEN
1098                            WRITE ( 21 )  time_since_reference_point,          &
1099                                          do2d_xy_time_count(av), av
1100                         ENDIF
1101#endif
1102                         DO  i = 0, io_blocks-1
1103                            IF ( i == io_group )  THEN
1104                               WRITE ( 21 )  nxlg, nxrg, nysg, nyng, nysg, nyng
1105                               WRITE ( 21 )  local_2d
1106                            ENDIF
1107#if defined( __parallel )
1108                            CALL MPI_BARRIER( comm2d, ierr )
1109#endif
1110                         ENDDO
1111
1112                      ELSE
1113!
1114!--                      PE0 receives partial arrays from all processors and
1115!--                      then outputs them. Here a barrier has to be set,
1116!--                      because otherwise "-MPI- FATAL: Remote protocol queue
1117!--                      full" may occur.
1118                         CALL MPI_BARRIER( comm2d, ierr )
1119
1120                         ngp = ( nxrg-nxlg+1 ) * ( nyng-nysg+1 )
1121                         IF ( myid == 0 )  THEN
1122!
1123!--                         Local array can be relocated directly.
1124                            total_2d(nxlg:nxrg,nysg:nyng) = local_2d
1125!
1126!--                         Receive data from all other PEs.
1127                            DO  n = 1, numprocs-1
1128!
1129!--                            Receive index limits first, then array.
1130!--                            Index limits are received in arbitrary order from
1131!--                            the PEs.
1132                               CALL MPI_RECV( ind(1), 4, MPI_INTEGER,          &
1133                                              MPI_ANY_SOURCE, 0, comm2d,       &
1134                                              status, ierr )
1135                               sender = status(MPI_SOURCE)
1136                               DEALLOCATE( local_2d )
1137                               ALLOCATE( local_2d(ind(1):ind(2),ind(3):ind(4)) )
1138                               CALL MPI_RECV( local_2d(ind(1),ind(3)), ngp,    &
1139                                              MPI_REAL, sender, 1, comm2d,     &
1140                                              status, ierr )
1141                               total_2d(ind(1):ind(2),ind(3):ind(4)) = local_2d
1142                            ENDDO
1143!
1144!--                         Relocate the local array for the next loop increment
1145                            DEALLOCATE( local_2d )
1146                            ALLOCATE( local_2d(nxlg:nxrg,nysg:nyng) )
1147
1148#if defined( __netcdf )
1149                            IF ( two_d ) THEN
1150                               nc_stat = NF90_PUT_VAR( id_set_xy(av),       &
1151                                                       id_var_do2d(av,if),  &
1152                                                   total_2d(0:nx+1,0:ny+1), &
1153                             start = (/ 1, 1, 1, do2d_xy_time_count(av) /), &
1154                                             count = (/ nx+2, ny+2, 1, 1 /) )
1155                            ELSE
1156                               nc_stat = NF90_PUT_VAR( id_set_xy(av),       &
1157                                                       id_var_do2d(av,if),  &
1158                                                   total_2d(0:nx+1,0:ny+1), &
1159                            start = (/ 1, 1, is, do2d_xy_time_count(av) /), &
1160                                             count = (/ nx+2, ny+2, 1, 1 /) )
1161                            ENDIF
1162                            CALL netcdf_handle_error( 'data_output_2d', 54 )
1163#endif
1164
1165                         ELSE
1166!
1167!--                         First send the local index limits to PE0
1168                            ind(1) = nxlg; ind(2) = nxrg
1169                            ind(3) = nysg; ind(4) = nyng
1170                            CALL MPI_SEND( ind(1), 4, MPI_INTEGER, 0, 0,       &
1171                                           comm2d, ierr )
1172!
1173!--                         Send data to PE0
1174                            CALL MPI_SEND( local_2d(nxlg,nysg), ngp,           &
1175                                           MPI_REAL, 0, 1, comm2d, ierr )
1176                         ENDIF
1177!
1178!--                      A barrier has to be set, because otherwise some PEs may
1179!--                      proceed too fast so that PE0 may receive wrong data on
1180!--                      tag 0
1181                         CALL MPI_BARRIER( comm2d, ierr )
1182                      ENDIF
1183
1184                   ENDIF
1185#else
1186#if defined( __netcdf )
1187                   IF ( two_d ) THEN
1188                      nc_stat = NF90_PUT_VAR( id_set_xy(av),                &
1189                                              id_var_do2d(av,if),           &
1190                                             local_2d(nxl:nxr+1,nys:nyn+1), &
1191                             start = (/ 1, 1, 1, do2d_xy_time_count(av) /), &
1192                                           count = (/ nx+2, ny+2, 1, 1 /) )
1193                   ELSE
1194                      nc_stat = NF90_PUT_VAR( id_set_xy(av),                &
1195                                              id_var_do2d(av,if),           &
1196                                             local_2d(nxl:nxr+1,nys:nyn+1), &
1197                            start = (/ 1, 1, is, do2d_xy_time_count(av) /), &
1198                                           count = (/ nx+2, ny+2, 1, 1 /) )
1199                   ENDIF
1200                   CALL netcdf_handle_error( 'data_output_2d', 447 )
1201#endif
1202#endif
1203                   do2d_xy_n = do2d_xy_n + 1
1204!
1205!--                For 2D-arrays (e.g. u*) only one cross-section is available.
1206!--                Hence exit loop of output levels.
1207                   IF ( two_d )  THEN
1208                      IF ( netcdf_data_format < 5 )  two_d = .FALSE.
1209                      EXIT loop1
1210                   ENDIF
1211
1212                CASE ( 'xz' )
1213!
1214!--                Update the netCDF xz cross section time axis.
1215!--                In case of parallel output, this is only done by PE0
1216!--                to increase the performance.
1217                   IF ( simulated_time /= do2d_xz_last_time(av) )  THEN
1218                      do2d_xz_time_count(av) = do2d_xz_time_count(av) + 1
1219                      do2d_xz_last_time(av)  = simulated_time
1220                      IF ( myid == 0 )  THEN
1221                         IF ( .NOT. data_output_2d_on_each_pe  &
1222                              .OR.  netcdf_data_format > 4 )   &
1223                         THEN
1224#if defined( __netcdf )
1225                            nc_stat = NF90_PUT_VAR( id_set_xz(av),             &
1226                                                    id_var_time_xz(av),        &
1227                                             (/ time_since_reference_point /), &
1228                                         start = (/ do2d_xz_time_count(av) /), &
1229                                                    count = (/ 1 /) )
1230                            CALL netcdf_handle_error( 'data_output_2d', 56 )
1231#endif
1232                         ENDIF
1233                      ENDIF
1234                   ENDIF
1235
1236!
1237!--                If required, carry out averaging along y
1238                   IF ( section(is,s_ind) == -1 )  THEN
1239
1240                      ALLOCATE( local_2d_l(nxlg:nxrg,nzb_do:nzt_do) )
1241                      local_2d_l = 0.0_wp
1242                      ngp = ( nxrg-nxlg + 1 ) * ( nzt_do-nzb_do + 1 )
1243!
1244!--                   First local averaging on the PE
1245                      DO  k = nzb_do, nzt_do
1246                         DO  j = nys, nyn
1247                            DO  i = nxlg, nxrg
1248                               local_2d_l(i,k) = local_2d_l(i,k) +             &
1249                                                 local_pf(i,j,k)
1250                            ENDDO
1251                         ENDDO
1252                      ENDDO
1253#if defined( __parallel )
1254!
1255!--                   Now do the averaging over all PEs along y
1256                      IF ( collective_wait )  CALL MPI_BARRIER( comm2d, ierr )
1257                      CALL MPI_ALLREDUCE( local_2d_l(nxlg,nzb_do),                &
1258                                          local_2d(nxlg,nzb_do), ngp, MPI_REAL,   &
1259                                          MPI_SUM, comm1dy, ierr )
1260#else
1261                      local_2d = local_2d_l
1262#endif
1263                      local_2d = local_2d / ( ny + 1.0_wp )
1264
1265                      DEALLOCATE( local_2d_l )
1266
1267                   ELSE
1268!
1269!--                   Just store the respective section on the local array
1270!--                   (but only if it is available on this PE!)
1271                      IF ( section(is,s_ind) >= nys  .AND.  section(is,s_ind) <= nyn ) &
1272                      THEN
1273                         local_2d = local_pf(:,section(is,s_ind),nzb_do:nzt_do)
1274                      ENDIF
1275
1276                   ENDIF
1277
1278#if defined( __parallel )
1279                   IF ( netcdf_data_format > 4 )  THEN
1280!
1281!--                   Output in netCDF4/HDF5 format.
1282!--                   Output only on those PEs where the respective cross
1283!--                   sections reside. Cross sections averaged along y are
1284!--                   output on the respective first PE along y (myidy=0).
1285                      IF ( ( section(is,s_ind) >= nys  .AND.                   &
1286                             section(is,s_ind) <= nyn )  .OR.                  &
1287                           ( section(is,s_ind) == -1  .AND.  myidy == 0 ) )  THEN
1288#if defined( __netcdf )
1289!
1290!--                      For parallel output, all cross sections are first
1291!--                      stored here on a local array and will be written to the
1292!--                      output file afterwards to increase the performance.
1293                         DO  i = nxlg, nxrg
1294                            DO  k = nzb_do, nzt_do
1295                               local_2d_sections_l(i,is,k) = local_2d(i,k)
1296                            ENDDO
1297                         ENDDO
1298#endif
1299                      ENDIF
1300
1301                   ELSE
1302
1303                      IF ( data_output_2d_on_each_pe )  THEN
1304!
1305!--                      Output of partial arrays on each PE. If the cross
1306!--                      section does not reside on the PE, output special
1307!--                      index values.
1308#if defined( __netcdf )
1309                         IF ( myid == 0 )  THEN
1310                            WRITE ( 22 )  time_since_reference_point,          &
1311                                          do2d_xz_time_count(av), av
1312                         ENDIF
1313#endif
1314                         DO  i = 0, io_blocks-1
1315                            IF ( i == io_group )  THEN
1316                               IF ( ( section(is,s_ind) >= nys  .AND.          &
1317                                      section(is,s_ind) <= nyn )  .OR.         &
1318                                    ( section(is,s_ind) == -1  .AND.           &
1319                                      nys-1 == -1 ) )                          &
1320                               THEN
1321                                  WRITE (22)  nxlg, nxrg, nzb_do, nzt_do, nzb, nzt+1
1322                                  WRITE (22)  local_2d
1323                               ELSE
1324                                  WRITE (22)  -1, -1, -1, -1, -1, -1
1325                               ENDIF
1326                            ENDIF
1327#if defined( __parallel )
1328                            CALL MPI_BARRIER( comm2d, ierr )
1329#endif
1330                         ENDDO
1331
1332                      ELSE
1333!
1334!--                      PE0 receives partial arrays from all processors of the
1335!--                      respective cross section and outputs them. Here a
1336!--                      barrier has to be set, because otherwise
1337!--                      "-MPI- FATAL: Remote protocol queue full" may occur.
1338                         CALL MPI_BARRIER( comm2d, ierr )
1339
1340                         ngp = ( nxrg-nxlg + 1 ) * ( nzt_do-nzb_do + 1 )
1341                         IF ( myid == 0 )  THEN
1342!
1343!--                         Local array can be relocated directly.
1344                            IF ( ( section(is,s_ind) >= nys  .AND.              &
1345                                   section(is,s_ind) <= nyn )  .OR.             &
1346                                 ( section(is,s_ind) == -1  .AND.               &
1347                                   nys-1 == -1 ) )  THEN
1348                               total_2d(nxlg:nxrg,nzb_do:nzt_do) = local_2d
1349                            ENDIF
1350!
1351!--                         Receive data from all other PEs.
1352                            DO  n = 1, numprocs-1
1353!
1354!--                            Receive index limits first, then array.
1355!--                            Index limits are received in arbitrary order from
1356!--                            the PEs.
1357                               CALL MPI_RECV( ind(1), 4, MPI_INTEGER,          &
1358                                              MPI_ANY_SOURCE, 0, comm2d,       &
1359                                              status, ierr )
1360!
1361!--                            Not all PEs have data for XZ-cross-section.
1362                               IF ( ind(1) /= -9999 )  THEN
1363                                  sender = status(MPI_SOURCE)
1364                                  DEALLOCATE( local_2d )
1365                                  ALLOCATE( local_2d(ind(1):ind(2),            &
1366                                                     ind(3):ind(4)) )
1367                                  CALL MPI_RECV( local_2d(ind(1),ind(3)), ngp, &
1368                                                 MPI_REAL, sender, 1, comm2d,  &
1369                                                 status, ierr )
1370                                  total_2d(ind(1):ind(2),ind(3):ind(4)) =      &
1371                                                                        local_2d
1372                               ENDIF
1373                            ENDDO
1374!
1375!--                         Relocate the local array for the next loop increment
1376                            DEALLOCATE( local_2d )
1377                            ALLOCATE( local_2d(nxlg:nxrg,nzb_do:nzt_do) )
1378
1379#if defined( __netcdf )
1380                            nc_stat = NF90_PUT_VAR( id_set_xz(av),          &
1381                                                 id_var_do2d(av,if),        &
1382                                                 total_2d(0:nx+1,nzb_do:nzt_do),&
1383                            start = (/ 1, is, 1, do2d_xz_time_count(av) /), &
1384                                             count = (/ nx+2, 1, nzt_do-nzb_do+1, 1 /) )
1385                            CALL netcdf_handle_error( 'data_output_2d', 58 )
1386#endif
1387
1388                         ELSE
1389!
1390!--                         If the cross section resides on the PE, send the
1391!--                         local index limits, otherwise send -9999 to PE0.
1392                            IF ( ( section(is,s_ind) >= nys  .AND.              &
1393                                   section(is,s_ind) <= nyn )  .OR.             &
1394                                 ( section(is,s_ind) == -1  .AND.  nys-1 == -1 ) ) &
1395                            THEN
1396                               ind(1) = nxlg; ind(2) = nxrg
1397                               ind(3) = nzb_do;   ind(4) = nzt_do
1398                            ELSE
1399                               ind(1) = -9999; ind(2) = -9999
1400                               ind(3) = -9999; ind(4) = -9999
1401                            ENDIF
1402                            CALL MPI_SEND( ind(1), 4, MPI_INTEGER, 0, 0,       &
1403                                           comm2d, ierr )
1404!
1405!--                         If applicable, send data to PE0.
1406                            IF ( ind(1) /= -9999 )  THEN
1407                               CALL MPI_SEND( local_2d(nxlg,nzb_do), ngp,         &
1408                                              MPI_REAL, 0, 1, comm2d, ierr )
1409                            ENDIF
1410                         ENDIF
1411!
1412!--                      A barrier has to be set, because otherwise some PEs may
1413!--                      proceed too fast so that PE0 may receive wrong data on
1414!--                      tag 0
1415                         CALL MPI_BARRIER( comm2d, ierr )
1416                      ENDIF
1417
1418                   ENDIF
1419#else
1420#if defined( __netcdf )
1421                   nc_stat = NF90_PUT_VAR( id_set_xz(av),                   &
1422                                           id_var_do2d(av,if),              &
1423                                           local_2d(nxl:nxr+1,nzb_do:nzt_do),   &
1424                            start = (/ 1, is, 1, do2d_xz_time_count(av) /), &
1425                                           count = (/ nx+2, 1, nzt_do-nzb_do+1, 1 /) )
1426                   CALL netcdf_handle_error( 'data_output_2d', 451 )
1427#endif
1428#endif
1429                   do2d_xz_n = do2d_xz_n + 1
1430
1431                CASE ( 'yz' )
1432!
1433!--                Update the netCDF yz cross section time axis.
1434!--                In case of parallel output, this is only done by PE0
1435!--                to increase the performance.
1436                   IF ( simulated_time /= do2d_yz_last_time(av) )  THEN
1437                      do2d_yz_time_count(av) = do2d_yz_time_count(av) + 1
1438                      do2d_yz_last_time(av)  = simulated_time
1439                      IF ( myid == 0 )  THEN
1440                         IF ( .NOT. data_output_2d_on_each_pe  &
1441                              .OR.  netcdf_data_format > 4 )   &
1442                         THEN
1443#if defined( __netcdf )
1444                            nc_stat = NF90_PUT_VAR( id_set_yz(av),             &
1445                                                    id_var_time_yz(av),        &
1446                                             (/ time_since_reference_point /), &
1447                                         start = (/ do2d_yz_time_count(av) /), &
1448                                                    count = (/ 1 /) )
1449                            CALL netcdf_handle_error( 'data_output_2d', 59 )
1450#endif
1451                         ENDIF
1452                      ENDIF
1453                   ENDIF
1454
1455!
1456!--                If required, carry out averaging along x
1457                   IF ( section(is,s_ind) == -1 )  THEN
1458
1459                      ALLOCATE( local_2d_l(nysg:nyng,nzb_do:nzt_do) )
1460                      local_2d_l = 0.0_wp
1461                      ngp = ( nyng-nysg+1 ) * ( nzt_do-nzb_do+1 )
1462!
1463!--                   First local averaging on the PE
1464                      DO  k = nzb_do, nzt_do
1465                         DO  j = nysg, nyng
1466                            DO  i = nxl, nxr
1467                               local_2d_l(j,k) = local_2d_l(j,k) +             &
1468                                                 local_pf(i,j,k)
1469                            ENDDO
1470                         ENDDO
1471                      ENDDO
1472#if defined( __parallel )
1473!
1474!--                   Now do the averaging over all PEs along x
1475                      IF ( collective_wait )  CALL MPI_BARRIER( comm2d, ierr )
1476                      CALL MPI_ALLREDUCE( local_2d_l(nysg,nzb_do),                &
1477                                          local_2d(nysg,nzb_do), ngp, MPI_REAL,   &
1478                                          MPI_SUM, comm1dx, ierr )
1479#else
1480                      local_2d = local_2d_l
1481#endif
1482                      local_2d = local_2d / ( nx + 1.0_wp )
1483
1484                      DEALLOCATE( local_2d_l )
1485
1486                   ELSE
1487!
1488!--                   Just store the respective section on the local array
1489!--                   (but only if it is available on this PE!)
1490                      IF ( section(is,s_ind) >= nxl  .AND.  section(is,s_ind) <= nxr ) &
1491                      THEN
1492                         local_2d = local_pf(section(is,s_ind),:,nzb_do:nzt_do)
1493                      ENDIF
1494
1495                   ENDIF
1496
1497#if defined( __parallel )
1498                   IF ( netcdf_data_format > 4 )  THEN
1499!
1500!--                   Output in netCDF4/HDF5 format.
1501!--                   Output only on those PEs where the respective cross
1502!--                   sections reside. Cross sections averaged along x are
1503!--                   output on the respective first PE along x (myidx=0).
1504                      IF ( ( section(is,s_ind) >= nxl  .AND.                       &
1505                             section(is,s_ind) <= nxr )  .OR.                      &
1506                           ( section(is,s_ind) == -1  .AND.  myidx == 0 ) )  THEN
1507#if defined( __netcdf )
1508!
1509!--                      For parallel output, all cross sections are first
1510!--                      stored here on a local array and will be written to the
1511!--                      output file afterwards to increase the performance.
1512                         DO  j = nysg, nyng
1513                            DO  k = nzb_do, nzt_do
1514                               local_2d_sections_l(is,j,k) = local_2d(j,k)
1515                            ENDDO
1516                         ENDDO
1517#endif
1518                      ENDIF
1519
1520                   ELSE
1521
1522                      IF ( data_output_2d_on_each_pe )  THEN
1523!
1524!--                      Output of partial arrays on each PE. If the cross
1525!--                      section does not reside on the PE, output special
1526!--                      index values.
1527#if defined( __netcdf )
1528                         IF ( myid == 0 )  THEN
1529                            WRITE ( 23 )  time_since_reference_point,          &
1530                                          do2d_yz_time_count(av), av
1531                         ENDIF
1532#endif
1533                         DO  i = 0, io_blocks-1
1534                            IF ( i == io_group )  THEN
1535                               IF ( ( section(is,s_ind) >= nxl  .AND.          &
1536                                      section(is,s_ind) <= nxr )  .OR.         &
1537                                    ( section(is,s_ind) == -1  .AND.           &
1538                                      nxl-1 == -1 ) )                          &
1539                               THEN
1540                                  WRITE (23)  nysg, nyng, nzb_do, nzt_do, nzb, nzt+1
1541                                  WRITE (23)  local_2d
1542                               ELSE
1543                                  WRITE (23)  -1, -1, -1, -1, -1, -1
1544                               ENDIF
1545                            ENDIF
1546#if defined( __parallel )
1547                            CALL MPI_BARRIER( comm2d, ierr )
1548#endif
1549                         ENDDO
1550
1551                      ELSE
1552!
1553!--                      PE0 receives partial arrays from all processors of the
1554!--                      respective cross section and outputs them. Here a
1555!--                      barrier has to be set, because otherwise
1556!--                      "-MPI- FATAL: Remote protocol queue full" may occur.
1557                         CALL MPI_BARRIER( comm2d, ierr )
1558
1559                         ngp = ( nyng-nysg+1 ) * ( nzt_do-nzb_do+1 )
1560                         IF ( myid == 0 )  THEN
1561!
1562!--                         Local array can be relocated directly.
1563                            IF ( ( section(is,s_ind) >= nxl  .AND.             &
1564                                   section(is,s_ind) <= nxr )   .OR.           &
1565                                 ( section(is,s_ind) == -1  .AND.  nxl-1 == -1 ) ) &
1566                            THEN
1567                               total_2d(nysg:nyng,nzb_do:nzt_do) = local_2d
1568                            ENDIF
1569!
1570!--                         Receive data from all other PEs.
1571                            DO  n = 1, numprocs-1
1572!
1573!--                            Receive index limits first, then array.
1574!--                            Index limits are received in arbitrary order from
1575!--                            the PEs.
1576                               CALL MPI_RECV( ind(1), 4, MPI_INTEGER,          &
1577                                              MPI_ANY_SOURCE, 0, comm2d,       &
1578                                              status, ierr )
1579!
1580!--                            Not all PEs have data for YZ-cross-section.
1581                               IF ( ind(1) /= -9999 )  THEN
1582                                  sender = status(MPI_SOURCE)
1583                                  DEALLOCATE( local_2d )
1584                                  ALLOCATE( local_2d(ind(1):ind(2),            &
1585                                                     ind(3):ind(4)) )
1586                                  CALL MPI_RECV( local_2d(ind(1),ind(3)), ngp, &
1587                                                 MPI_REAL, sender, 1, comm2d,  &
1588                                                 status, ierr )
1589                                  total_2d(ind(1):ind(2),ind(3):ind(4)) =      &
1590                                                                        local_2d
1591                               ENDIF
1592                            ENDDO
1593!
1594!--                         Relocate the local array for the next loop increment
1595                            DEALLOCATE( local_2d )
1596                            ALLOCATE( local_2d(nysg:nyng,nzb_do:nzt_do) )
1597
1598#if defined( __netcdf )
1599                            nc_stat = NF90_PUT_VAR( id_set_yz(av),          &
1600                                                 id_var_do2d(av,if),        &
1601                                                 total_2d(0:ny+1,nzb_do:nzt_do),&
1602                            start = (/ is, 1, 1, do2d_yz_time_count(av) /), &
1603                                             count = (/ 1, ny+2, nzt_do-nzb_do+1, 1 /) )
1604                            CALL netcdf_handle_error( 'data_output_2d', 61 )
1605#endif
1606
1607                         ELSE
1608!
1609!--                         If the cross section resides on the PE, send the
1610!--                         local index limits, otherwise send -9999 to PE0.
1611                            IF ( ( section(is,s_ind) >= nxl  .AND.              &
1612                                   section(is,s_ind) <= nxr )  .OR.             &
1613                                 ( section(is,s_ind) == -1  .AND.  nxl-1 == -1 ) ) &
1614                            THEN
1615                               ind(1) = nysg; ind(2) = nyng
1616                               ind(3) = nzb_do;   ind(4) = nzt_do
1617                            ELSE
1618                               ind(1) = -9999; ind(2) = -9999
1619                               ind(3) = -9999; ind(4) = -9999
1620                            ENDIF
1621                            CALL MPI_SEND( ind(1), 4, MPI_INTEGER, 0, 0,       &
1622                                           comm2d, ierr )
1623!
1624!--                         If applicable, send data to PE0.
1625                            IF ( ind(1) /= -9999 )  THEN
1626                               CALL MPI_SEND( local_2d(nysg,nzb_do), ngp,         &
1627                                              MPI_REAL, 0, 1, comm2d, ierr )
1628                            ENDIF
1629                         ENDIF
1630!
1631!--                      A barrier has to be set, because otherwise some PEs may
1632!--                      proceed too fast so that PE0 may receive wrong data on
1633!--                      tag 0
1634                         CALL MPI_BARRIER( comm2d, ierr )
1635                      ENDIF
1636
1637                   ENDIF
1638#else
1639#if defined( __netcdf )
1640                   nc_stat = NF90_PUT_VAR( id_set_yz(av),                   &
1641                                           id_var_do2d(av,if),              &
1642                                           local_2d(nys:nyn+1,nzb_do:nzt_do),   &
1643                            start = (/ is, 1, 1, do2d_xz_time_count(av) /), &
1644                                           count = (/ 1, ny+2, nzt_do-nzb_do+1, 1 /) )
1645                   CALL netcdf_handle_error( 'data_output_2d', 452 )
1646#endif
1647#endif
1648                   do2d_yz_n = do2d_yz_n + 1
1649
1650             END SELECT
1651
1652             is = is + 1
1653          ENDDO loop1
1654
1655!
1656!--       For parallel output, all data were collected before on a local array
1657!--       and are written now to the netcdf file. This must be done to increase
1658!--       the performance of the parallel output.
1659#if defined( __netcdf )
1660          IF ( netcdf_data_format > 4 )  THEN
1661
1662                SELECT CASE ( mode )
1663
1664                   CASE ( 'xy' )
1665                      IF ( two_d ) THEN
1666                         nis = 1
1667                         two_d = .FALSE.
1668                      ELSE
1669                         nis = ns
1670                      ENDIF
1671!
1672!--                   Do not output redundant ghost point data except for the
1673!--                   boundaries of the total domain.
1674                      IF ( nxr == nx  .AND.  nyn /= ny )  THEN
1675                         nc_stat = NF90_PUT_VAR( id_set_xy(av),                &
1676                                                 id_var_do2d(av,if),           &
1677                                                 local_2d_sections(nxl:nxr+1,  &
1678                                                    nys:nyn,1:nis),            &
1679                                                 start = (/ nxl+1, nys+1, 1,   &
1680                                                    do2d_xy_time_count(av) /), &
1681                                                 count = (/ nxr-nxl+2,         &
1682                                                            nyn-nys+1, nis, 1  &
1683                                                          /) )
1684                      ELSEIF ( nxr /= nx  .AND.  nyn == ny )  THEN
1685                         nc_stat = NF90_PUT_VAR( id_set_xy(av),                &
1686                                                 id_var_do2d(av,if),           &
1687                                                 local_2d_sections(nxl:nxr,    &
1688                                                    nys:nyn+1,1:nis),          &
1689                                                 start = (/ nxl+1, nys+1, 1,   &
1690                                                    do2d_xy_time_count(av) /), &
1691                                                 count = (/ nxr-nxl+1,         &
1692                                                            nyn-nys+2, nis, 1  &
1693                                                          /) )
1694                      ELSEIF ( nxr == nx  .AND.  nyn == ny )  THEN
1695                         nc_stat = NF90_PUT_VAR( id_set_xy(av),                &
1696                                                 id_var_do2d(av,if),           &
1697                                                 local_2d_sections(nxl:nxr+1,  &
1698                                                    nys:nyn+1,1:nis),          &
1699                                                 start = (/ nxl+1, nys+1, 1,   &
1700                                                    do2d_xy_time_count(av) /), &
1701                                                 count = (/ nxr-nxl+2,         &
1702                                                            nyn-nys+2, nis, 1  &
1703                                                          /) )
1704                      ELSE
1705                         nc_stat = NF90_PUT_VAR( id_set_xy(av),                &
1706                                                 id_var_do2d(av,if),           &
1707                                                 local_2d_sections(nxl:nxr,    &
1708                                                    nys:nyn,1:nis),            &
1709                                                 start = (/ nxl+1, nys+1, 1,   &
1710                                                    do2d_xy_time_count(av) /), &
1711                                                 count = (/ nxr-nxl+1,         &
1712                                                            nyn-nys+1, nis, 1  &
1713                                                          /) )
1714                      ENDIF   
1715
1716                      CALL netcdf_handle_error( 'data_output_2d', 55 )
1717
1718                   CASE ( 'xz' )
1719!
1720!--                   First, all PEs get the information of all cross-sections.
1721!--                   Then the data are written to the output file by all PEs
1722!--                   while NF90_COLLECTIVE is set in subroutine
1723!--                   define_netcdf_header. Although redundant information are
1724!--                   written to the output file in that case, the performance
1725!--                   is significantly better compared to the case where only
1726!--                   the first row of PEs in x-direction (myidx = 0) is given
1727!--                   the output while NF90_INDEPENDENT is set.
1728                      IF ( npey /= 1 )  THEN
1729                         
1730#if defined( __parallel )
1731!
1732!--                      Distribute data over all PEs along y
1733                         ngp = ( nxrg-nxlg+1 ) * ( nzt_do-nzb_do+1 ) * ns
1734                         IF ( collective_wait ) CALL MPI_BARRIER( comm2d, ierr )
1735                         CALL MPI_ALLREDUCE( local_2d_sections_l(nxlg,1,nzb_do),  &
1736                                             local_2d_sections(nxlg,1,nzb_do),    &
1737                                             ngp, MPI_REAL, MPI_SUM, comm1dy,  &
1738                                             ierr )
1739#else
1740                         local_2d_sections = local_2d_sections_l
1741#endif
1742                      ENDIF
1743!
1744!--                   Do not output redundant ghost point data except for the
1745!--                   boundaries of the total domain.
1746                      IF ( nxr == nx )  THEN
1747                         nc_stat = NF90_PUT_VAR( id_set_xz(av),                &
1748                                             id_var_do2d(av,if),               & 
1749                                             local_2d_sections(nxl:nxr+1,1:ns, &
1750                                                nzb_do:nzt_do),                &
1751                                             start = (/ nxl+1, 1, 1,           &
1752                                                do2d_xz_time_count(av) /),     &
1753                                             count = (/ nxr-nxl+2, ns, nzt_do-nzb_do+1,  &
1754                                                        1 /) )
1755                      ELSE
1756                         nc_stat = NF90_PUT_VAR( id_set_xz(av),                &
1757                                             id_var_do2d(av,if),               &
1758                                             local_2d_sections(nxl:nxr,1:ns,   &
1759                                                nzb_do:nzt_do),                &
1760                                             start = (/ nxl+1, 1, 1,           &
1761                                                do2d_xz_time_count(av) /),     &
1762                                             count = (/ nxr-nxl+1, ns, nzt_do-nzb_do+1,  &
1763                                                1 /) )
1764                      ENDIF
1765
1766                      CALL netcdf_handle_error( 'data_output_2d', 57 )
1767
1768                   CASE ( 'yz' )
1769!
1770!--                   First, all PEs get the information of all cross-sections.
1771!--                   Then the data are written to the output file by all PEs
1772!--                   while NF90_COLLECTIVE is set in subroutine
1773!--                   define_netcdf_header. Although redundant information are
1774!--                   written to the output file in that case, the performance
1775!--                   is significantly better compared to the case where only
1776!--                   the first row of PEs in y-direction (myidy = 0) is given
1777!--                   the output while NF90_INDEPENDENT is set.
1778                      IF ( npex /= 1 )  THEN
1779
1780#if defined( __parallel )
1781!
1782!--                      Distribute data over all PEs along x
1783                         ngp = ( nyng-nysg+1 ) * ( nzt-nzb + 2 ) * ns
1784                         IF ( collective_wait ) CALL MPI_BARRIER( comm2d, ierr )
1785                         CALL MPI_ALLREDUCE( local_2d_sections_l(1,nysg,nzb_do),  &
1786                                             local_2d_sections(1,nysg,nzb_do),    &
1787                                             ngp, MPI_REAL, MPI_SUM, comm1dx,  &
1788                                             ierr )
1789#else
1790                         local_2d_sections = local_2d_sections_l
1791#endif
1792                      ENDIF
1793!
1794!--                   Do not output redundant ghost point data except for the
1795!--                   boundaries of the total domain.
1796                      IF ( nyn == ny )  THEN
1797                         nc_stat = NF90_PUT_VAR( id_set_yz(av),                &
1798                                             id_var_do2d(av,if),               &
1799                                             local_2d_sections(1:ns,           &
1800                                                nys:nyn+1,nzb_do:nzt_do),      &
1801                                             start = (/ 1, nys+1, 1,           &
1802                                                do2d_yz_time_count(av) /),     &
1803                                             count = (/ ns, nyn-nys+2,         &
1804                                                        nzt_do-nzb_do+1, 1 /) )
1805                      ELSE
1806                         nc_stat = NF90_PUT_VAR( id_set_yz(av),                &
1807                                             id_var_do2d(av,if),               &
1808                                             local_2d_sections(1:ns,nys:nyn,   &
1809                                                nzb_do:nzt_do),                &
1810                                             start = (/ 1, nys+1, 1,           &
1811                                                do2d_yz_time_count(av) /),     &
1812                                             count = (/ ns, nyn-nys+1,         &
1813                                                        nzt_do-nzb_do+1, 1 /) )
1814                      ENDIF
1815
1816                      CALL netcdf_handle_error( 'data_output_2d', 60 )
1817
1818                   CASE DEFAULT
1819                      message_string = 'unknown cross-section: ' // TRIM( mode )
1820                      CALL message( 'data_output_2d', 'PA0180', 1, 2, 0, 6, 0 )
1821
1822                END SELECT                     
1823
1824          ENDIF
1825#endif
1826       ENDIF
1827
1828       if = if + 1
1829       l = MAX( 2, LEN_TRIM( do2d(av,if) ) )
1830       do2d_mode = do2d(av,if)(l-1:l)
1831
1832    ENDDO
1833
1834!
1835!-- Deallocate temporary arrays.
1836    IF ( ALLOCATED( level_z ) )  DEALLOCATE( level_z )
1837    IF ( netcdf_data_format > 4 )  THEN
1838       DEALLOCATE( local_pf, local_2d, local_2d_sections )
1839       IF( mode == 'xz' .OR. mode == 'yz' ) DEALLOCATE( local_2d_sections_l )
1840    ENDIF
1841#if defined( __parallel )
1842    IF ( .NOT.  data_output_2d_on_each_pe  .AND.  myid == 0 )  THEN
1843       DEALLOCATE( total_2d )
1844    ENDIF
1845#endif
1846
1847!
1848!-- Close plot output file.
1849    file_id = 20 + s_ind
1850
1851    IF ( data_output_2d_on_each_pe )  THEN
1852       DO  i = 0, io_blocks-1
1853          IF ( i == io_group )  THEN
1854             CALL close_file( file_id )
1855          ENDIF
1856#if defined( __parallel )
1857          CALL MPI_BARRIER( comm2d, ierr )
1858#endif
1859       ENDDO
1860    ELSE
1861       IF ( myid == 0 )  CALL close_file( file_id )
1862    ENDIF
1863
1864    CALL cpu_log( log_point(3), 'data_output_2d', 'stop' )
1865
1866 END SUBROUTINE data_output_2d
Note: See TracBrowser for help on using the repository browser.