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

Last change on this file since 2195 was 2191, checked in by raasch, 8 years ago

last commit documented

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