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

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

last commit documented

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