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

Last change on this file since 2798 was 2798, checked in by suehring, 6 years ago

Bugfix initialization of %pt_surface array; Output of surface temperature also for default-type surfaces

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