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

Last change on this file since 1977 was 1977, checked in by maronga, 5 years ago

last commit documented

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