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

Last change on this file since 4598 was 4559, checked in by raasch, 4 years ago

files re-formatted to follow the PALM coding standard

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