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

Last change on this file since 2016 was 2001, checked in by knoop, 8 years ago

last commit documented

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