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

Last change on this file since 1788 was 1788, checked in by maronga, 6 years ago

added support for water and paved surfaced in land surface model / minor changes

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