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

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

last commit documented

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