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

Last change on this file since 2059 was 2032, checked in by knoop, 8 years ago

last commit documented

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