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

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

last commit documented

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