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

Last change on this file since 1784 was 1784, checked in by raasch, 6 years ago

last commit documented

  • Property svn:keywords set to Id
File size: 87.8 KB
Line 
1!> @file data_output_2d.f90
2!--------------------------------------------------------------------------------!
3! This file is part of PALM.
4!
5! PALM is free software: you can redistribute it and/or modify it under the terms
6! of the GNU General Public License as published by the Free Software Foundation,
7! either version 3 of the License, or (at your option) any later version.
8!
9! PALM is distributed in the hope that it will be useful, but WITHOUT ANY
10! WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
11! A PARTICULAR PURPOSE.  See the GNU General Public License for more details.
12!
13! You should have received a copy of the GNU General Public License along with
14! PALM. If not, see <http://www.gnu.org/licenses/>.
15!
16! Copyright 1997-2015 Leibniz Universitaet Hannover
17!--------------------------------------------------------------------------------!
18!
19! Current revisions:
20! -----------------
21!
22!
23! Former revisions:
24! -----------------
25! $Id: data_output_2d.f90 1784 2016-03-06 19:14:40Z raasch $
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, 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 DEFAULT
1257!
1258!--             User defined quantity
1259                CALL user_data_output_2d( av, do2d(av,if), found, grid,        &
1260                                          local_pf, two_d, nzb_do, nzt_do )
1261                resorted = .TRUE.
1262
1263                IF ( grid == 'zu' )  THEN
1264                   IF ( mode == 'xy' )  level_z = zu
1265                ELSEIF ( grid == 'zw' )  THEN
1266                   IF ( mode == 'xy' )  level_z = zw
1267                ELSEIF ( grid == 'zu1' ) THEN
1268                   IF ( mode == 'xy' )  level_z(nzb+1) = zu(nzb+1)
1269                ELSEIF ( grid == 'zs' ) THEN
1270                   IF ( mode == 'xy' )  level_z = zs
1271                ENDIF
1272
1273                IF ( .NOT. found )  THEN
1274                   message_string = 'no output provided for: ' //              &
1275                                    TRIM( do2d(av,if) )
1276                   CALL message( 'data_output_2d', 'PA0181', 0, 0, 0, 6, 0 )
1277                ENDIF
1278
1279          END SELECT
1280
1281!
1282!--       Resort the array to be output, if not done above
1283          IF ( .NOT. resorted )  THEN
1284             DO  i = nxlg, nxrg
1285                DO  j = nysg, nyng
1286                   DO  k = nzb_do, nzt_do
1287                      local_pf(i,j,k) = to_be_resorted(k,j,i)
1288                   ENDDO
1289                ENDDO
1290             ENDDO
1291          ENDIF
1292
1293!
1294!--       Output of the individual cross-sections, depending on the cross-
1295!--       section mode chosen.
1296          is = 1
1297   loop1: DO WHILE ( section(is,s) /= -9999  .OR.  two_d )
1298
1299             SELECT CASE ( mode )
1300
1301                CASE ( 'xy' )
1302!
1303!--                Determine the cross section index
1304                   IF ( two_d )  THEN
1305                      layer_xy = nzb+1
1306                   ELSE
1307                      layer_xy = section(is,s)
1308                   ENDIF
1309
1310!
1311!--                Exit the loop for layers beyond the data output domain
1312!--                (used for soil model)
1313                   IF ( layer_xy > nzt_do )  THEN
1314                      EXIT loop1
1315                   ENDIF
1316
1317!
1318!--                Update the netCDF xy cross section time axis.
1319!--                In case of parallel output, this is only done by PE0
1320!--                to increase the performance.
1321                   IF ( simulated_time /= do2d_xy_last_time(av) )  THEN
1322                      do2d_xy_time_count(av) = do2d_xy_time_count(av) + 1
1323                      do2d_xy_last_time(av)  = simulated_time
1324                      IF ( myid == 0 )  THEN
1325                         IF ( .NOT. data_output_2d_on_each_pe  &
1326                              .OR.  netcdf_data_format > 4 )   &
1327                         THEN
1328#if defined( __netcdf )
1329                            nc_stat = NF90_PUT_VAR( id_set_xy(av),             &
1330                                                    id_var_time_xy(av),        &
1331                                             (/ time_since_reference_point /), &
1332                                         start = (/ do2d_xy_time_count(av) /), &
1333                                                    count = (/ 1 /) )
1334                            CALL netcdf_handle_error( 'data_output_2d', 53 )
1335#endif
1336                         ENDIF
1337                      ENDIF
1338                   ENDIF
1339!
1340!--                If required, carry out averaging along z
1341                   IF ( section(is,s) == -1  .AND.  .NOT. two_d )  THEN
1342
1343                      local_2d = 0.0_wp
1344!
1345!--                   Carry out the averaging (all data are on the PE)
1346                      DO  k = nzb_do, nzt_do
1347                         DO  j = nysg, nyng
1348                            DO  i = nxlg, nxrg
1349                               local_2d(i,j) = local_2d(i,j) + local_pf(i,j,k)
1350                            ENDDO
1351                         ENDDO
1352                      ENDDO
1353
1354                      local_2d = local_2d / ( nzt_do - nzb_do + 1.0_wp)
1355
1356                   ELSE
1357!
1358!--                   Just store the respective section on the local array
1359                      local_2d = local_pf(:,:,layer_xy)
1360
1361                   ENDIF
1362
1363#if defined( __parallel )
1364                   IF ( netcdf_data_format > 4 )  THEN
1365!
1366!--                   Parallel output in netCDF4/HDF5 format.
1367                      IF ( two_d ) THEN
1368                         iis = 1
1369                      ELSE
1370                         iis = is
1371                      ENDIF
1372
1373#if defined( __netcdf )
1374!
1375!--                   For parallel output, all cross sections are first stored
1376!--                   here on a local array and will be written to the output
1377!--                   file afterwards to increase the performance.
1378                      DO  i = nxlg, nxrg
1379                         DO  j = nysg, nyng
1380                            local_2d_sections(i,j,iis) = local_2d(i,j)
1381                         ENDDO
1382                      ENDDO
1383#endif
1384                   ELSE
1385
1386                      IF ( data_output_2d_on_each_pe )  THEN
1387!
1388!--                      Output of partial arrays on each PE
1389#if defined( __netcdf )
1390                         IF ( myid == 0 )  THEN
1391                            WRITE ( 21 )  time_since_reference_point,          &
1392                                          do2d_xy_time_count(av), av
1393                         ENDIF
1394#endif
1395                         DO  i = 0, io_blocks-1
1396                            IF ( i == io_group )  THEN
1397                               WRITE ( 21 )  nxlg, nxrg, nysg, nyng, nysg, nyng
1398                               WRITE ( 21 )  local_2d
1399                            ENDIF
1400#if defined( __parallel )
1401                            CALL MPI_BARRIER( comm2d, ierr )
1402#endif
1403                         ENDDO
1404
1405                      ELSE
1406!
1407!--                      PE0 receives partial arrays from all processors and
1408!--                      then outputs them. Here a barrier has to be set,
1409!--                      because otherwise "-MPI- FATAL: Remote protocol queue
1410!--                      full" may occur.
1411                         CALL MPI_BARRIER( comm2d, ierr )
1412
1413                         ngp = ( nxrg-nxlg+1 ) * ( nyng-nysg+1 )
1414                         IF ( myid == 0 )  THEN
1415!
1416!--                         Local array can be relocated directly.
1417                            total_2d(nxlg:nxrg,nysg:nyng) = local_2d
1418!
1419!--                         Receive data from all other PEs.
1420                            DO  n = 1, numprocs-1
1421!
1422!--                            Receive index limits first, then array.
1423!--                            Index limits are received in arbitrary order from
1424!--                            the PEs.
1425                               CALL MPI_RECV( ind(1), 4, MPI_INTEGER,          &
1426                                              MPI_ANY_SOURCE, 0, comm2d,       &
1427                                              status, ierr )
1428                               sender = status(MPI_SOURCE)
1429                               DEALLOCATE( local_2d )
1430                               ALLOCATE( local_2d(ind(1):ind(2),ind(3):ind(4)) )
1431                               CALL MPI_RECV( local_2d(ind(1),ind(3)), ngp,    &
1432                                              MPI_REAL, sender, 1, comm2d,     &
1433                                              status, ierr )
1434                               total_2d(ind(1):ind(2),ind(3):ind(4)) = local_2d
1435                            ENDDO
1436!
1437!--                         Relocate the local array for the next loop increment
1438                            DEALLOCATE( local_2d )
1439                            ALLOCATE( local_2d(nxlg:nxrg,nysg:nyng) )
1440
1441#if defined( __netcdf )
1442                            IF ( two_d ) THEN
1443                               nc_stat = NF90_PUT_VAR( id_set_xy(av),       &
1444                                                       id_var_do2d(av,if),  &
1445                                                   total_2d(0:nx+1,0:ny+1), &
1446                             start = (/ 1, 1, 1, do2d_xy_time_count(av) /), &
1447                                             count = (/ nx+2, ny+2, 1, 1 /) )
1448                            ELSE
1449                               nc_stat = NF90_PUT_VAR( id_set_xy(av),       &
1450                                                       id_var_do2d(av,if),  &
1451                                                   total_2d(0:nx+1,0:ny+1), &
1452                            start = (/ 1, 1, is, do2d_xy_time_count(av) /), &
1453                                             count = (/ nx+2, ny+2, 1, 1 /) )
1454                            ENDIF
1455                            CALL netcdf_handle_error( 'data_output_2d', 54 )
1456#endif
1457
1458                         ELSE
1459!
1460!--                         First send the local index limits to PE0
1461                            ind(1) = nxlg; ind(2) = nxrg
1462                            ind(3) = nysg; ind(4) = nyng
1463                            CALL MPI_SEND( ind(1), 4, MPI_INTEGER, 0, 0,       &
1464                                           comm2d, ierr )
1465!
1466!--                         Send data to PE0
1467                            CALL MPI_SEND( local_2d(nxlg,nysg), ngp,           &
1468                                           MPI_REAL, 0, 1, comm2d, ierr )
1469                         ENDIF
1470!
1471!--                      A barrier has to be set, because otherwise some PEs may
1472!--                      proceed too fast so that PE0 may receive wrong data on
1473!--                      tag 0
1474                         CALL MPI_BARRIER( comm2d, ierr )
1475                      ENDIF
1476
1477                   ENDIF
1478#else
1479#if defined( __netcdf )
1480                   IF ( two_d ) THEN
1481                      nc_stat = NF90_PUT_VAR( id_set_xy(av),                &
1482                                              id_var_do2d(av,if),           &
1483                                             local_2d(nxl:nxr+1,nys:nyn+1), &
1484                             start = (/ 1, 1, 1, do2d_xy_time_count(av) /), &
1485                                           count = (/ nx+2, ny+2, 1, 1 /) )
1486                   ELSE
1487                      nc_stat = NF90_PUT_VAR( id_set_xy(av),                &
1488                                              id_var_do2d(av,if),           &
1489                                             local_2d(nxl:nxr+1,nys:nyn+1), &
1490                            start = (/ 1, 1, is, do2d_xy_time_count(av) /), &
1491                                           count = (/ nx+2, ny+2, 1, 1 /) )
1492                   ENDIF
1493                   CALL netcdf_handle_error( 'data_output_2d', 447 )
1494#endif
1495#endif
1496                   do2d_xy_n = do2d_xy_n + 1
1497!
1498!--                For 2D-arrays (e.g. u*) only one cross-section is available.
1499!--                Hence exit loop of output levels.
1500                   IF ( two_d )  THEN
1501                      IF ( netcdf_data_format < 5 )  two_d = .FALSE.
1502                      EXIT loop1
1503                   ENDIF
1504
1505                CASE ( 'xz' )
1506!
1507!--                Update the netCDF xz cross section time axis.
1508!--                In case of parallel output, this is only done by PE0
1509!--                to increase the performance.
1510                   IF ( simulated_time /= do2d_xz_last_time(av) )  THEN
1511                      do2d_xz_time_count(av) = do2d_xz_time_count(av) + 1
1512                      do2d_xz_last_time(av)  = simulated_time
1513                      IF ( myid == 0 )  THEN
1514                         IF ( .NOT. data_output_2d_on_each_pe  &
1515                              .OR.  netcdf_data_format > 4 )   &
1516                         THEN
1517#if defined( __netcdf )
1518                            nc_stat = NF90_PUT_VAR( id_set_xz(av),             &
1519                                                    id_var_time_xz(av),        &
1520                                             (/ time_since_reference_point /), &
1521                                         start = (/ do2d_xz_time_count(av) /), &
1522                                                    count = (/ 1 /) )
1523                            CALL netcdf_handle_error( 'data_output_2d', 56 )
1524#endif
1525                         ENDIF
1526                      ENDIF
1527                   ENDIF
1528
1529!
1530!--                If required, carry out averaging along y
1531                   IF ( section(is,s) == -1 )  THEN
1532
1533                      ALLOCATE( local_2d_l(nxlg:nxrg,nzb_do:nzt_do) )
1534                      local_2d_l = 0.0_wp
1535                      ngp = ( nxrg-nxlg + 1 ) * ( nzt_do-nzb_do + 1 )
1536!
1537!--                   First local averaging on the PE
1538                      DO  k = nzb_do, nzt_do
1539                         DO  j = nys, nyn
1540                            DO  i = nxlg, nxrg
1541                               local_2d_l(i,k) = local_2d_l(i,k) +             &
1542                                                 local_pf(i,j,k)
1543                            ENDDO
1544                         ENDDO
1545                      ENDDO
1546#if defined( __parallel )
1547!
1548!--                   Now do the averaging over all PEs along y
1549                      IF ( collective_wait )  CALL MPI_BARRIER( comm2d, ierr )
1550                      CALL MPI_ALLREDUCE( local_2d_l(nxlg,nzb_do),                &
1551                                          local_2d(nxlg,nzb_do), ngp, MPI_REAL,   &
1552                                          MPI_SUM, comm1dy, ierr )
1553#else
1554                      local_2d = local_2d_l
1555#endif
1556                      local_2d = local_2d / ( ny + 1.0_wp )
1557
1558                      DEALLOCATE( local_2d_l )
1559
1560                   ELSE
1561!
1562!--                   Just store the respective section on the local array
1563!--                   (but only if it is available on this PE!)
1564                      IF ( section(is,s) >= nys  .AND.  section(is,s) <= nyn ) &
1565                      THEN
1566                         local_2d = local_pf(:,section(is,s),nzb_do:nzt_do)
1567                      ENDIF
1568
1569                   ENDIF
1570
1571#if defined( __parallel )
1572                   IF ( netcdf_data_format > 4 )  THEN
1573!
1574!--                   Output in netCDF4/HDF5 format.
1575!--                   Output only on those PEs where the respective cross
1576!--                   sections reside. Cross sections averaged along y are
1577!--                   output on the respective first PE along y (myidy=0).
1578                      IF ( ( section(is,s) >= nys  .AND.                       &
1579                             section(is,s) <= nyn )  .OR.                      &
1580                           ( section(is,s) == -1  .AND.  myidy == 0 ) )  THEN
1581#if defined( __netcdf )
1582!
1583!--                      For parallel output, all cross sections are first
1584!--                      stored here on a local array and will be written to the
1585!--                      output file afterwards to increase the performance.
1586                         DO  i = nxlg, nxrg
1587                            DO  k = nzb_do, nzt_do
1588                               local_2d_sections_l(i,is,k) = local_2d(i,k)
1589                            ENDDO
1590                         ENDDO
1591#endif
1592                      ENDIF
1593
1594                   ELSE
1595
1596                      IF ( data_output_2d_on_each_pe )  THEN
1597!
1598!--                      Output of partial arrays on each PE. If the cross
1599!--                      section does not reside on the PE, output special
1600!--                      index values.
1601#if defined( __netcdf )
1602                         IF ( myid == 0 )  THEN
1603                            WRITE ( 22 )  time_since_reference_point,          &
1604                                          do2d_xz_time_count(av), av
1605                         ENDIF
1606#endif
1607                         DO  i = 0, io_blocks-1
1608                            IF ( i == io_group )  THEN
1609                               IF ( ( section(is,s) >= nys  .AND.              &
1610                                      section(is,s) <= nyn )  .OR.             &
1611                                    ( section(is,s) == -1  .AND.               &
1612                                      nys-1 == -1 ) )                          &
1613                               THEN
1614                                  WRITE (22)  nxlg, nxrg, nzb_do, nzt_do, nzb, nzt+1
1615                                  WRITE (22)  local_2d
1616                               ELSE
1617                                  WRITE (22)  -1, -1, -1, -1, -1, -1
1618                               ENDIF
1619                            ENDIF
1620#if defined( __parallel )
1621                            CALL MPI_BARRIER( comm2d, ierr )
1622#endif
1623                         ENDDO
1624
1625                      ELSE
1626!
1627!--                      PE0 receives partial arrays from all processors of the
1628!--                      respective cross section and outputs them. Here a
1629!--                      barrier has to be set, because otherwise
1630!--                      "-MPI- FATAL: Remote protocol queue full" may occur.
1631                         CALL MPI_BARRIER( comm2d, ierr )
1632
1633                         ngp = ( nxrg-nxlg + 1 ) * ( nzt_do-nzb_do + 1 )
1634                         IF ( myid == 0 )  THEN
1635!
1636!--                         Local array can be relocated directly.
1637                            IF ( ( section(is,s) >= nys  .AND.                 &
1638                                   section(is,s) <= nyn )  .OR.                &
1639                                 ( section(is,s) == -1  .AND.  nys-1 == -1 ) ) &
1640                            THEN
1641                               total_2d(nxlg:nxrg,nzb_do:nzt_do) = local_2d
1642                            ENDIF
1643!
1644!--                         Receive data from all other PEs.
1645                            DO  n = 1, numprocs-1
1646!
1647!--                            Receive index limits first, then array.
1648!--                            Index limits are received in arbitrary order from
1649!--                            the PEs.
1650                               CALL MPI_RECV( ind(1), 4, MPI_INTEGER,          &
1651                                              MPI_ANY_SOURCE, 0, comm2d,       &
1652                                              status, ierr )
1653!
1654!--                            Not all PEs have data for XZ-cross-section.
1655                               IF ( ind(1) /= -9999 )  THEN
1656                                  sender = status(MPI_SOURCE)
1657                                  DEALLOCATE( local_2d )
1658                                  ALLOCATE( local_2d(ind(1):ind(2),            &
1659                                                     ind(3):ind(4)) )
1660                                  CALL MPI_RECV( local_2d(ind(1),ind(3)), ngp, &
1661                                                 MPI_REAL, sender, 1, comm2d,  &
1662                                                 status, ierr )
1663                                  total_2d(ind(1):ind(2),ind(3):ind(4)) =      &
1664                                                                        local_2d
1665                               ENDIF
1666                            ENDDO
1667!
1668!--                         Relocate the local array for the next loop increment
1669                            DEALLOCATE( local_2d )
1670                            ALLOCATE( local_2d(nxlg:nxrg,nzb_do:nzt_do) )
1671
1672#if defined( __netcdf )
1673                            nc_stat = NF90_PUT_VAR( id_set_xz(av),          &
1674                                                 id_var_do2d(av,if),        &
1675                                                 total_2d(0:nx+1,nzb_do:nzt_do),&
1676                            start = (/ 1, is, 1, do2d_xz_time_count(av) /), &
1677                                             count = (/ nx+2, 1, nzt_do-nzb_do+1, 1 /) )
1678                            CALL netcdf_handle_error( 'data_output_2d', 58 )
1679#endif
1680
1681                         ELSE
1682!
1683!--                         If the cross section resides on the PE, send the
1684!--                         local index limits, otherwise send -9999 to PE0.
1685                            IF ( ( section(is,s) >= nys  .AND.                 &
1686                                   section(is,s) <= nyn )  .OR.                &
1687                                 ( section(is,s) == -1  .AND.  nys-1 == -1 ) ) &
1688                            THEN
1689                               ind(1) = nxlg; ind(2) = nxrg
1690                               ind(3) = nzb_do;   ind(4) = nzt_do
1691                            ELSE
1692                               ind(1) = -9999; ind(2) = -9999
1693                               ind(3) = -9999; ind(4) = -9999
1694                            ENDIF
1695                            CALL MPI_SEND( ind(1), 4, MPI_INTEGER, 0, 0,       &
1696                                           comm2d, ierr )
1697!
1698!--                         If applicable, send data to PE0.
1699                            IF ( ind(1) /= -9999 )  THEN
1700                               CALL MPI_SEND( local_2d(nxlg,nzb_do), ngp,         &
1701                                              MPI_REAL, 0, 1, comm2d, ierr )
1702                            ENDIF
1703                         ENDIF
1704!
1705!--                      A barrier has to be set, because otherwise some PEs may
1706!--                      proceed too fast so that PE0 may receive wrong data on
1707!--                      tag 0
1708                         CALL MPI_BARRIER( comm2d, ierr )
1709                      ENDIF
1710
1711                   ENDIF
1712#else
1713#if defined( __netcdf )
1714                   nc_stat = NF90_PUT_VAR( id_set_xz(av),                   &
1715                                           id_var_do2d(av,if),              &
1716                                           local_2d(nxl:nxr+1,nzb_do:nzt_do),   &
1717                            start = (/ 1, is, 1, do2d_xz_time_count(av) /), &
1718                                           count = (/ nx+2, 1, nzt_do-nzb_do+1, 1 /) )
1719                   CALL netcdf_handle_error( 'data_output_2d', 451 )
1720#endif
1721#endif
1722                   do2d_xz_n = do2d_xz_n + 1
1723
1724                CASE ( 'yz' )
1725!
1726!--                Update the netCDF yz cross section time axis.
1727!--                In case of parallel output, this is only done by PE0
1728!--                to increase the performance.
1729                   IF ( simulated_time /= do2d_yz_last_time(av) )  THEN
1730                      do2d_yz_time_count(av) = do2d_yz_time_count(av) + 1
1731                      do2d_yz_last_time(av)  = simulated_time
1732                      IF ( myid == 0 )  THEN
1733                         IF ( .NOT. data_output_2d_on_each_pe  &
1734                              .OR.  netcdf_data_format > 4 )   &
1735                         THEN
1736#if defined( __netcdf )
1737                            nc_stat = NF90_PUT_VAR( id_set_yz(av),             &
1738                                                    id_var_time_yz(av),        &
1739                                             (/ time_since_reference_point /), &
1740                                         start = (/ do2d_yz_time_count(av) /), &
1741                                                    count = (/ 1 /) )
1742                            CALL netcdf_handle_error( 'data_output_2d', 59 )
1743#endif
1744                         ENDIF
1745                      ENDIF
1746                   ENDIF
1747
1748!
1749!--                If required, carry out averaging along x
1750                   IF ( section(is,s) == -1 )  THEN
1751
1752                      ALLOCATE( local_2d_l(nysg:nyng,nzb_do:nzt_do) )
1753                      local_2d_l = 0.0_wp
1754                      ngp = ( nyng-nysg+1 ) * ( nzt_do-nzb_do+1 )
1755!
1756!--                   First local averaging on the PE
1757                      DO  k = nzb_do, nzt_do
1758                         DO  j = nysg, nyng
1759                            DO  i = nxl, nxr
1760                               local_2d_l(j,k) = local_2d_l(j,k) +             &
1761                                                 local_pf(i,j,k)
1762                            ENDDO
1763                         ENDDO
1764                      ENDDO
1765#if defined( __parallel )
1766!
1767!--                   Now do the averaging over all PEs along x
1768                      IF ( collective_wait )  CALL MPI_BARRIER( comm2d, ierr )
1769                      CALL MPI_ALLREDUCE( local_2d_l(nysg,nzb_do),                &
1770                                          local_2d(nysg,nzb_do), ngp, MPI_REAL,   &
1771                                          MPI_SUM, comm1dx, ierr )
1772#else
1773                      local_2d = local_2d_l
1774#endif
1775                      local_2d = local_2d / ( nx + 1.0_wp )
1776
1777                      DEALLOCATE( local_2d_l )
1778
1779                   ELSE
1780!
1781!--                   Just store the respective section on the local array
1782!--                   (but only if it is available on this PE!)
1783                      IF ( section(is,s) >= nxl  .AND.  section(is,s) <= nxr ) &
1784                      THEN
1785                         local_2d = local_pf(section(is,s),:,nzb_do:nzt_do)
1786                      ENDIF
1787
1788                   ENDIF
1789
1790#if defined( __parallel )
1791                   IF ( netcdf_data_format > 4 )  THEN
1792!
1793!--                   Output in netCDF4/HDF5 format.
1794!--                   Output only on those PEs where the respective cross
1795!--                   sections reside. Cross sections averaged along x are
1796!--                   output on the respective first PE along x (myidx=0).
1797                      IF ( ( section(is,s) >= nxl  .AND.                       &
1798                             section(is,s) <= nxr )  .OR.                      &
1799                           ( section(is,s) == -1  .AND.  myidx == 0 ) )  THEN
1800#if defined( __netcdf )
1801!
1802!--                      For parallel output, all cross sections are first
1803!--                      stored here on a local array and will be written to the
1804!--                      output file afterwards to increase the performance.
1805                         DO  j = nysg, nyng
1806                            DO  k = nzb_do, nzt_do
1807                               local_2d_sections_l(is,j,k) = local_2d(j,k)
1808                            ENDDO
1809                         ENDDO
1810#endif
1811                      ENDIF
1812
1813                   ELSE
1814
1815                      IF ( data_output_2d_on_each_pe )  THEN
1816!
1817!--                      Output of partial arrays on each PE. If the cross
1818!--                      section does not reside on the PE, output special
1819!--                      index values.
1820#if defined( __netcdf )
1821                         IF ( myid == 0 )  THEN
1822                            WRITE ( 23 )  time_since_reference_point,          &
1823                                          do2d_yz_time_count(av), av
1824                         ENDIF
1825#endif
1826                         DO  i = 0, io_blocks-1
1827                            IF ( i == io_group )  THEN
1828                               IF ( ( section(is,s) >= nxl  .AND.              &
1829                                      section(is,s) <= nxr )  .OR.             &
1830                                    ( section(is,s) == -1  .AND.               &
1831                                      nxl-1 == -1 ) )                          &
1832                               THEN
1833                                  WRITE (23)  nysg, nyng, nzb_do, nzt_do, nzb, nzt+1
1834                                  WRITE (23)  local_2d
1835                               ELSE
1836                                  WRITE (23)  -1, -1, -1, -1, -1, -1
1837                               ENDIF
1838                            ENDIF
1839#if defined( __parallel )
1840                            CALL MPI_BARRIER( comm2d, ierr )
1841#endif
1842                         ENDDO
1843
1844                      ELSE
1845!
1846!--                      PE0 receives partial arrays from all processors of the
1847!--                      respective cross section and outputs them. Here a
1848!--                      barrier has to be set, because otherwise
1849!--                      "-MPI- FATAL: Remote protocol queue full" may occur.
1850                         CALL MPI_BARRIER( comm2d, ierr )
1851
1852                         ngp = ( nyng-nysg+1 ) * ( nzt_do-nzb_do+1 )
1853                         IF ( myid == 0 )  THEN
1854!
1855!--                         Local array can be relocated directly.
1856                            IF ( ( section(is,s) >= nxl  .AND.                 &
1857                                   section(is,s) <= nxr )   .OR.               &
1858                                 ( section(is,s) == -1  .AND.  nxl-1 == -1 ) ) &
1859                            THEN
1860                               total_2d(nysg:nyng,nzb_do:nzt_do) = local_2d
1861                            ENDIF
1862!
1863!--                         Receive data from all other PEs.
1864                            DO  n = 1, numprocs-1
1865!
1866!--                            Receive index limits first, then array.
1867!--                            Index limits are received in arbitrary order from
1868!--                            the PEs.
1869                               CALL MPI_RECV( ind(1), 4, MPI_INTEGER,          &
1870                                              MPI_ANY_SOURCE, 0, comm2d,       &
1871                                              status, ierr )
1872!
1873!--                            Not all PEs have data for YZ-cross-section.
1874                               IF ( ind(1) /= -9999 )  THEN
1875                                  sender = status(MPI_SOURCE)
1876                                  DEALLOCATE( local_2d )
1877                                  ALLOCATE( local_2d(ind(1):ind(2),            &
1878                                                     ind(3):ind(4)) )
1879                                  CALL MPI_RECV( local_2d(ind(1),ind(3)), ngp, &
1880                                                 MPI_REAL, sender, 1, comm2d,  &
1881                                                 status, ierr )
1882                                  total_2d(ind(1):ind(2),ind(3):ind(4)) =      &
1883                                                                        local_2d
1884                               ENDIF
1885                            ENDDO
1886!
1887!--                         Relocate the local array for the next loop increment
1888                            DEALLOCATE( local_2d )
1889                            ALLOCATE( local_2d(nysg:nyng,nzb_do:nzt_do) )
1890
1891#if defined( __netcdf )
1892                            nc_stat = NF90_PUT_VAR( id_set_yz(av),          &
1893                                                 id_var_do2d(av,if),        &
1894                                                 total_2d(0:ny+1,nzb_do:nzt_do),&
1895                            start = (/ is, 1, 1, do2d_yz_time_count(av) /), &
1896                                             count = (/ 1, ny+2, nzt_do-nzb_do+1, 1 /) )
1897                            CALL netcdf_handle_error( 'data_output_2d', 61 )
1898#endif
1899
1900                         ELSE
1901!
1902!--                         If the cross section resides on the PE, send the
1903!--                         local index limits, otherwise send -9999 to PE0.
1904                            IF ( ( section(is,s) >= nxl  .AND.                 &
1905                                   section(is,s) <= nxr )  .OR.                &
1906                                 ( section(is,s) == -1  .AND.  nxl-1 == -1 ) ) &
1907                            THEN
1908                               ind(1) = nysg; ind(2) = nyng
1909                               ind(3) = nzb_do;   ind(4) = nzt_do
1910                            ELSE
1911                               ind(1) = -9999; ind(2) = -9999
1912                               ind(3) = -9999; ind(4) = -9999
1913                            ENDIF
1914                            CALL MPI_SEND( ind(1), 4, MPI_INTEGER, 0, 0,       &
1915                                           comm2d, ierr )
1916!
1917!--                         If applicable, send data to PE0.
1918                            IF ( ind(1) /= -9999 )  THEN
1919                               CALL MPI_SEND( local_2d(nysg,nzb_do), ngp,         &
1920                                              MPI_REAL, 0, 1, comm2d, ierr )
1921                            ENDIF
1922                         ENDIF
1923!
1924!--                      A barrier has to be set, because otherwise some PEs may
1925!--                      proceed too fast so that PE0 may receive wrong data on
1926!--                      tag 0
1927                         CALL MPI_BARRIER( comm2d, ierr )
1928                      ENDIF
1929
1930                   ENDIF
1931#else
1932#if defined( __netcdf )
1933                   nc_stat = NF90_PUT_VAR( id_set_yz(av),                   &
1934                                           id_var_do2d(av,if),              &
1935                                           local_2d(nys:nyn+1,nzb_do:nzt_do),   &
1936                            start = (/ is, 1, 1, do2d_xz_time_count(av) /), &
1937                                           count = (/ 1, ny+2, nzt_do-nzb_do+1, 1 /) )
1938                   CALL netcdf_handle_error( 'data_output_2d', 452 )
1939#endif
1940#endif
1941                   do2d_yz_n = do2d_yz_n + 1
1942
1943             END SELECT
1944
1945             is = is + 1
1946          ENDDO loop1
1947
1948!
1949!--       For parallel output, all data were collected before on a local array
1950!--       and are written now to the netcdf file. This must be done to increase
1951!--       the performance of the parallel output.
1952#if defined( __netcdf )
1953          IF ( netcdf_data_format > 4 )  THEN
1954
1955                SELECT CASE ( mode )
1956
1957                   CASE ( 'xy' )
1958                      IF ( two_d ) THEN
1959                         nis = 1
1960                         two_d = .FALSE.
1961                      ELSE
1962                         nis = ns
1963                      ENDIF
1964!
1965!--                   Do not output redundant ghost point data except for the
1966!--                   boundaries of the total domain.
1967                      IF ( nxr == nx  .AND.  nyn /= ny )  THEN
1968                         nc_stat = NF90_PUT_VAR( id_set_xy(av),                &
1969                                                 id_var_do2d(av,if),           &
1970                                                 local_2d_sections(nxl:nxr+1,  &
1971                                                    nys:nyn,1:nis),            &
1972                                                 start = (/ nxl+1, nys+1, 1,   &
1973                                                    do2d_xy_time_count(av) /), &
1974                                                 count = (/ nxr-nxl+2,         &
1975                                                            nyn-nys+1, nis, 1  &
1976                                                          /) )
1977                      ELSEIF ( nxr /= nx  .AND.  nyn == ny )  THEN
1978                         nc_stat = NF90_PUT_VAR( id_set_xy(av),                &
1979                                                 id_var_do2d(av,if),           &
1980                                                 local_2d_sections(nxl:nxr,    &
1981                                                    nys:nyn+1,1:nis),          &
1982                                                 start = (/ nxl+1, nys+1, 1,   &
1983                                                    do2d_xy_time_count(av) /), &
1984                                                 count = (/ nxr-nxl+1,         &
1985                                                            nyn-nys+2, nis, 1  &
1986                                                          /) )
1987                      ELSEIF ( nxr == nx  .AND.  nyn == ny )  THEN
1988                         nc_stat = NF90_PUT_VAR( id_set_xy(av),                &
1989                                                 id_var_do2d(av,if),           &
1990                                                 local_2d_sections(nxl:nxr+1,  &
1991                                                    nys:nyn+1,1:nis),          &
1992                                                 start = (/ nxl+1, nys+1, 1,   &
1993                                                    do2d_xy_time_count(av) /), &
1994                                                 count = (/ nxr-nxl+2,         &
1995                                                            nyn-nys+2, nis, 1  &
1996                                                          /) )
1997                      ELSE
1998                         nc_stat = NF90_PUT_VAR( id_set_xy(av),                &
1999                                                 id_var_do2d(av,if),           &
2000                                                 local_2d_sections(nxl:nxr,    &
2001                                                    nys:nyn,1:nis),            &
2002                                                 start = (/ nxl+1, nys+1, 1,   &
2003                                                    do2d_xy_time_count(av) /), &
2004                                                 count = (/ nxr-nxl+1,         &
2005                                                            nyn-nys+1, nis, 1  &
2006                                                          /) )
2007                      ENDIF   
2008
2009                      CALL netcdf_handle_error( 'data_output_2d', 55 )
2010
2011                   CASE ( 'xz' )
2012!
2013!--                   First, all PEs get the information of all cross-sections.
2014!--                   Then the data are written to the output file by all PEs
2015!--                   while NF90_COLLECTIVE is set in subroutine
2016!--                   define_netcdf_header. Although redundant information are
2017!--                   written to the output file in that case, the performance
2018!--                   is significantly better compared to the case where only
2019!--                   the first row of PEs in x-direction (myidx = 0) is given
2020!--                   the output while NF90_INDEPENDENT is set.
2021                      IF ( npey /= 1 )  THEN
2022                         
2023#if defined( __parallel )
2024!
2025!--                      Distribute data over all PEs along y
2026                         ngp = ( nxrg-nxlg+1 ) * ( nzt_do-nzb_do+1 ) * ns
2027                         IF ( collective_wait ) CALL MPI_BARRIER( comm2d, ierr )
2028                         CALL MPI_ALLREDUCE( local_2d_sections_l(nxlg,1,nzb_do),  &
2029                                             local_2d_sections(nxlg,1,nzb_do),    &
2030                                             ngp, MPI_REAL, MPI_SUM, comm1dy,  &
2031                                             ierr )
2032#else
2033                         local_2d_sections = local_2d_sections_l
2034#endif
2035                      ENDIF
2036!
2037!--                   Do not output redundant ghost point data except for the
2038!--                   boundaries of the total domain.
2039                      IF ( nxr == nx )  THEN
2040                         nc_stat = NF90_PUT_VAR( id_set_xz(av),                &
2041                                             id_var_do2d(av,if),               & 
2042                                             local_2d_sections(nxl:nxr+1,1:ns, &
2043                                                nzb_do:nzt_do),                &
2044                                             start = (/ nxl+1, 1, 1,           &
2045                                                do2d_xz_time_count(av) /),     &
2046                                             count = (/ nxr-nxl+2, ns, nzt_do-nzb_do+1,  &
2047                                                        1 /) )
2048                      ELSE
2049                         nc_stat = NF90_PUT_VAR( id_set_xz(av),                &
2050                                             id_var_do2d(av,if),               &
2051                                             local_2d_sections(nxl:nxr,1:ns,   &
2052                                                nzb_do:nzt_do),                &
2053                                             start = (/ nxl+1, 1, 1,           &
2054                                                do2d_xz_time_count(av) /),     &
2055                                             count = (/ nxr-nxl+1, ns, nzt_do-nzb_do+1,  &
2056                                                1 /) )
2057                      ENDIF
2058
2059                      CALL netcdf_handle_error( 'data_output_2d', 57 )
2060
2061                   CASE ( 'yz' )
2062!
2063!--                   First, all PEs get the information of all cross-sections.
2064!--                   Then the data are written to the output file by all PEs
2065!--                   while NF90_COLLECTIVE is set in subroutine
2066!--                   define_netcdf_header. Although redundant information are
2067!--                   written to the output file in that case, the performance
2068!--                   is significantly better compared to the case where only
2069!--                   the first row of PEs in y-direction (myidy = 0) is given
2070!--                   the output while NF90_INDEPENDENT is set.
2071                      IF ( npex /= 1 )  THEN
2072
2073#if defined( __parallel )
2074!
2075!--                      Distribute data over all PEs along x
2076                         ngp = ( nyng-nysg+1 ) * ( nzt-nzb + 2 ) * ns
2077                         IF ( collective_wait ) CALL MPI_BARRIER( comm2d, ierr )
2078                         CALL MPI_ALLREDUCE( local_2d_sections_l(1,nysg,nzb_do),  &
2079                                             local_2d_sections(1,nysg,nzb_do),    &
2080                                             ngp, MPI_REAL, MPI_SUM, comm1dx,  &
2081                                             ierr )
2082#else
2083                         local_2d_sections = local_2d_sections_l
2084#endif
2085                      ENDIF
2086!
2087!--                   Do not output redundant ghost point data except for the
2088!--                   boundaries of the total domain.
2089                      IF ( nyn == ny )  THEN
2090                         nc_stat = NF90_PUT_VAR( id_set_yz(av),                &
2091                                             id_var_do2d(av,if),               &
2092                                             local_2d_sections(1:ns,           &
2093                                                nys:nyn+1,nzb_do:nzt_do),      &
2094                                             start = (/ 1, nys+1, 1,           &
2095                                                do2d_yz_time_count(av) /),     &
2096                                             count = (/ ns, nyn-nys+2,         &
2097                                                        nzt_do-nzb_do+1, 1 /) )
2098                      ELSE
2099                         nc_stat = NF90_PUT_VAR( id_set_yz(av),                &
2100                                             id_var_do2d(av,if),               &
2101                                             local_2d_sections(1:ns,nys:nyn,   &
2102                                                nzb_do:nzt_do),                &
2103                                             start = (/ 1, nys+1, 1,           &
2104                                                do2d_yz_time_count(av) /),     &
2105                                             count = (/ ns, nyn-nys+1,         &
2106                                                        nzt_do-nzb_do+1, 1 /) )
2107                      ENDIF
2108
2109                      CALL netcdf_handle_error( 'data_output_2d', 60 )
2110
2111                   CASE DEFAULT
2112                      message_string = 'unknown cross-section: ' // TRIM( mode )
2113                      CALL message( 'data_output_2d', 'PA0180', 1, 2, 0, 6, 0 )
2114
2115                END SELECT                     
2116
2117          ENDIF
2118#endif
2119       ENDIF
2120
2121       if = if + 1
2122       l = MAX( 2, LEN_TRIM( do2d(av,if) ) )
2123       do2d_mode = do2d(av,if)(l-1:l)
2124
2125    ENDDO
2126
2127!
2128!-- Deallocate temporary arrays.
2129    IF ( ALLOCATED( level_z ) )  DEALLOCATE( level_z )
2130    IF ( netcdf_data_format > 4 )  THEN
2131       DEALLOCATE( local_pf, local_2d, local_2d_sections )
2132       IF( mode == 'xz' .OR. mode == 'yz' ) DEALLOCATE( local_2d_sections_l )
2133    ENDIF
2134#if defined( __parallel )
2135    IF ( .NOT.  data_output_2d_on_each_pe  .AND.  myid == 0 )  THEN
2136       DEALLOCATE( total_2d )
2137    ENDIF
2138#endif
2139
2140!
2141!-- Close plot output file.
2142    file_id = 20 + s
2143
2144    IF ( data_output_2d_on_each_pe )  THEN
2145       DO  i = 0, io_blocks-1
2146          IF ( i == io_group )  THEN
2147             CALL close_file( file_id )
2148          ENDIF
2149#if defined( __parallel )
2150          CALL MPI_BARRIER( comm2d, ierr )
2151#endif
2152       ENDDO
2153    ELSE
2154       IF ( myid == 0 )  CALL close_file( file_id )
2155    ENDIF
2156
2157    CALL cpu_log( log_point(3), 'data_output_2d', 'stop' )
2158
2159 END SUBROUTINE data_output_2d
Note: See TracBrowser for help on using the repository browser.