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

Last change on this file since 3637 was 3637, checked in by knoop, 4 years ago

M Makefile

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