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

Last change on this file since 3421 was 3421, checked in by gronemeier, 5 years ago

new surface-data output; renamed output variables (pt to theta, rho_air to rho, rho_ocean to rho_sea_water)

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